Cantera  3.1.0
Loading...
Searching...
No Matches
SolutionArray.cpp
Go to the documentation of this file.
1/**
2 * @file SolutionArray.cpp
3 * Definition file for class SolutionArray.
4 */
5
6// This file is part of Cantera. See License.txt in the top-level directory or
7// at https://cantera.org/license.txt for license and copyright information.
8
16#include <boost/algorithm/string.hpp>
17#include <boost/range/adaptor/reversed.hpp>
18#include <fstream>
19#include <sstream>
20
21
22namespace ba = boost::algorithm;
23
24const std::map<std::string, std::string> aliasMap = {
25 {"T", "temperature"},
26 {"P", "pressure"},
27 {"D", "density"},
28 {"Y", "mass-fractions"},
29 {"X", "mole-fractions"},
30 {"C", "coverages"},
31 {"U", "specific-internal-energy"},
32 {"V", "specific-volume"},
33 {"H", "specific-enthalpy"},
34 {"S", "specific-entropy"},
35 {"Q", "vapor-fraction"},
36};
37
38namespace Cantera
39{
40
41SolutionArray::SolutionArray(const shared_ptr<Solution>& sol,
42 int size, const AnyMap& meta)
43 : m_sol(sol)
44 , m_size(size)
45 , m_dataSize(size)
46 , m_meta(meta)
47{
48 if (!m_sol || !m_sol->thermo()) {
49 throw CanteraError("SolutionArray::SolutionArray",
50 "Unable to create SolutionArray from invalid Solution object.");
51 }
52 m_stride = m_sol->thermo()->stateSize();
53 m_sol->thermo()->addSpeciesLock();
54 m_data = make_shared<vector<double>>(m_dataSize * m_stride, 0.);
55 m_extra = make_shared<map<string, AnyValue>>();
56 m_order = make_shared<map<int, string>>();
57 for (size_t i = 0; i < m_dataSize; ++i) {
58 m_active.push_back(static_cast<int>(i));
59 }
60 reset();
61 m_apiShape.resize(1);
62 m_apiShape[0] = static_cast<long>(m_dataSize);
63}
64
65SolutionArray::SolutionArray(const SolutionArray& other,
66 const vector<int>& selected)
67 : m_sol(other.m_sol)
68 , m_size(selected.size())
69 , m_dataSize(other.m_data->size())
70 , m_stride(other.m_stride)
71 , m_data(other.m_data)
72 , m_extra(other.m_extra)
73 , m_order(other.m_order)
74 , m_shared(true)
75{
76 m_sol->thermo()->addSpeciesLock();
77 if (!other.m_shared) {
78 // direct slicing is possible
79 m_active = selected;
80 } else {
81 // slicing a previously sliced SolutionArray
82 m_active.clear();
83 m_active.reserve(selected.size());
84 for (auto loc : selected) {
85 m_active.push_back(other.m_active.at(loc));
86 }
87 }
88 for (auto loc : m_active) {
89 if (loc < 0 || loc >= (int)m_dataSize) {
90 IndexError("SolutionArray::SolutionArray", "indices", loc, m_dataSize);
91 }
92 }
93 set<int> unique(selected.begin(), selected.end());
94 if (unique.size() < selected.size()) {
95 throw CanteraError("SolutionArray::SolutionArray", "Indices must be unique.");
96 }
97}
98
99SolutionArray::~SolutionArray()
100{
101 m_sol->thermo()->removeSpeciesLock();
102}
103
104namespace { // restrict scope of helper functions to local translation unit
105
106template<class T>
107void resetSingle(AnyValue& extra, const vector<int>& slice);
108
109template<class T>
110AnyValue getSingle(const AnyValue& extra, const vector<int>& slice);
111
112template<class T>
113void setSingle(AnyValue& extra, const AnyValue& data, const vector<int>& slice);
114
115template<class T>
116void resizeSingle(AnyValue& extra, size_t size, const AnyValue& value);
117
118template<class T>
119void resetMulti(AnyValue& extra, const vector<int>& slice);
120
121template<class T>
122AnyValue getMulti(const AnyValue& extra, const vector<int>& slice);
123
124template<class T>
125void setMulti(AnyValue& extra, const AnyValue& data, const vector<int>& slice);
126
127template<class T>
128void resizeMulti(AnyValue& extra, size_t size, const AnyValue& value);
129
130template<class T>
131void setAuxiliarySingle(size_t loc, AnyValue& extra, const AnyValue& value);
132
133template<class T>
134void setAuxiliaryMulti(size_t loc, AnyValue& extra, const AnyValue& data);
135
136template<class T>
137void setScalar(AnyValue& extra, const AnyValue& data, const vector<int>& slice);
138
139} // end unnamed namespace
140
141void SolutionArray::reset()
142{
143 size_t nState = m_sol->thermo()->stateSize();
144 vector<double> state(nState);
145 m_sol->thermo()->saveState(state); // thermo contains current state
146 for (size_t k = 0; k < m_size; ++k) {
147 std::copy(state.begin(), state.end(), m_data->data() + m_active[k] * m_stride);
148 }
149 for (auto& [key, extra] : *m_extra) {
150 if (extra.is<void>()) {
151 // cannot reset placeholder (uninitialized component)
152 } else if (extra.isVector<double>()) {
153 resetSingle<double>(extra, m_active);
154 } else if (extra.isVector<long int>()) {
155 resetSingle<long int>(extra, m_active);
156 } else if (extra.isVector<string>()) {
157 resetSingle<string>(extra, m_active);
158 } else if (extra.isVector<vector<double>>()) {
159 resetMulti<double>(extra, m_active);
160 } else if (extra.isVector<vector<long int>>()) {
161 resetMulti<long int>(extra, m_active);
162 } else if (extra.isVector<vector<string>>()) {
163 resetMulti<string>(extra, m_active);
164 } else {
165 throw NotImplementedError("SolutionArray::reset",
166 "Unable to reset component '{}' with type '{}'.",
167 key, extra.type_str());
168 }
169 }
170}
171
172void SolutionArray::resize(int size)
173{
174 if (apiNdim() > 1) {
175 throw CanteraError("SolutionArray::resize",
176 "Resize is ambiguous for multi-dimensional arrays; use setApiShape "
177 "instead.");
178 }
179 if (m_data.use_count() > 1) {
180 throw CanteraError("SolutionArray::resize",
181 "Unable to resize as data are shared by multiple objects.");
182 }
183 _resize(static_cast<size_t>(size));
184 m_apiShape[0] = static_cast<long>(size);
185}
186
187void SolutionArray::setApiShape(const vector<long int>& shape)
188{
189 size_t size = 1;
190 for (auto dim : shape) {
191 size *= dim;
192 }
193 if (m_shared && size != m_size) {
194 throw CanteraError("SolutionArray::setApiShape",
195 "Unable to set shape of shared data as sizes are inconsistent:\n"
196 "active size is {} but shape implies {}.", m_size, size);
197 }
198 if (!m_shared && size != m_dataSize) {
199 if (m_data.use_count() > 1) {
200 throw CanteraError("SolutionArray::setApiShape",
201 "Unable to set shape as data are shared by multiple objects.");
202 }
203 _resize(size);
204 }
205 m_apiShape = shape;
206}
207
208void SolutionArray::_resize(size_t size)
209{
210 m_size = size;
211 m_dataSize = size;
212 m_data->resize(m_dataSize * m_stride, 0.);
213 for (auto& [key, data] : *m_extra) {
214 _resizeExtra(key);
215 }
216 m_active.clear();
217 for (size_t i = 0; i < m_dataSize; ++i) {
218 m_active.push_back(static_cast<int>(i));
219 }
220}
221
222namespace { // restrict scope of helper functions to local translation unit
223
224vector<string> doubleColumn(string name, const vector<double>& comp,
225 int rows, int width)
226{
227 // extract data for processing
228 vector<double> data;
229 vector<string> raw;
230 string notation = fmt::format("{{:{}.{}g}}", width, (width - 1) / 2);
231 int csize = static_cast<int>(comp.size());
232 int dots = csize + 1;
233 if (csize <= rows) {
234 for (const auto& val : comp) {
235 data.push_back(val);
236 raw.push_back(boost::trim_copy(fmt::format(fmt::runtime(notation), val)));
237 }
238 } else {
239 dots = (rows + 1) / 2;
240 for (int row = 0; row < dots; row++) {
241 data.push_back(comp[row]);
242 raw.push_back(boost::trim_copy(
243 fmt::format(fmt::runtime(notation), comp[row])));
244 }
245 for (int row = csize - rows / 2; row < csize; row++) {
246 data.push_back(comp[row]);
247 raw.push_back(boost::trim_copy(
248 fmt::format(fmt::runtime(notation), comp[row])));
249 }
250 }
251
252 // determine notation; all entries use identical formatting
253 bool isFloat = false;
254 bool isScientific = false;
255 size_t head = 0;
256 size_t tail = 0;
257 size_t exp = 0;
258 for (const auto& repr : raw) {
259 string name = repr;
260 if (name[0] == '-') {
261 // leading negative sign is not considered
262 name = name.substr(1);
263 }
264 if (name.find('e') != string::npos) {
265 // scientific notation
266 if (!isScientific) {
267 head = 1;
268 tail = name.find('e') - name.find('.');
269 exp = 4; // size of exponent
270 } else {
271 tail = std::max(tail, name.find('e') - name.find('.'));
272 }
273 isFloat = true;
274 isScientific = true;
275 } else if (name.find('.') != string::npos) {
276 // floating point notation
277 isFloat = true;
278 if (!isScientific) {
279 head = std::max(head, name.find('.'));
280 tail = std::max(tail, name.size() - name.find('.'));
281 }
282 } else {
283 head = std::max(head, name.size());
284 }
285 }
286 size_t maxLen = std::max((size_t)4, head + tail + exp + isFloat + 1);
287 size_t over = std::max(0, (int)name.size() - (int)maxLen);
288 if (isScientific) {
289 // at least one entry has scientific notation
290 notation = fmt::format(" {{:>{}.{}e}}", over + maxLen, tail);
291 } else if (isFloat) {
292 // at least one entry is a floating point
293 notation = fmt::format(" {{:>{}.{}f}}", over + maxLen, tail);
294 } else {
295 // all entries are integers
296 notation = fmt::format(" {{:>{}.0f}}", over + maxLen);
297 }
298 maxLen = fmt::format(fmt::runtime(notation), 0.).size();
299
300 // assemble output
301 string section = fmt::format("{{:>{}}}", maxLen);
302 vector<string> col = {fmt::format(fmt::runtime(section), name)};
303 int count = 0;
304 for (const auto& val : data) {
305 col.push_back(fmt::format(fmt::runtime(notation), val));
306 count++;
307 if (count == dots) {
308 col.push_back(fmt::format(fmt::runtime(section), "..."));
309 }
310 }
311 return col;
312}
313
314vector<string> integerColumn(string name, const vector<long int>& comp,
315 int rows, int width)
316{
317 // extract data for processing
318 vector<long int> data;
319 string notation = fmt::format("{{:{}}}", width);
320 size_t maxLen = 2; // minimum column width is 2
321 int csize = static_cast<int>(comp.size());
322 int dots = csize + 1;
323 if (csize <= rows) {
324 for (const auto& val : comp) {
325 data.push_back(val);
326 string formatted = boost::trim_copy(
327 fmt::format(fmt::runtime(notation), val));
328 if (formatted[0] == '-') {
329 formatted = formatted.substr(1);
330 }
331 maxLen = std::max(maxLen, formatted.size());
332 }
333 } else {
334 dots = (rows + 1) / 2;
335 for (int row = 0; row < dots; row++) {
336 data.push_back(comp[row]);
337 string formatted = boost::trim_copy(
338 fmt::format(fmt::runtime(notation), comp[row]));
339 if (formatted[0] == '-') {
340 formatted = formatted.substr(1);
341 }
342 maxLen = std::max(maxLen, formatted.size());
343 }
344 for (int row = csize - rows / 2; row < csize; row++) {
345 data.push_back(comp[row]);
346 string formatted = boost::trim_copy(
347 fmt::format(fmt::runtime(notation), comp[row]));
348 if (formatted[0] == '-') {
349 formatted = formatted.substr(1);
350 }
351 maxLen = std::max(maxLen, formatted.size());
352 }
353 }
354
355 if (name == "") {
356 // index column
357 notation = fmt::format("{{:<{}}}", maxLen);
358 } else {
359 // regular column
360 maxLen = std::max(maxLen, name.size());
361 notation = fmt::format(" {{:>{}}}", maxLen + 1);
362 }
363
364 // assemble output
365 vector<string> col = {fmt::format(fmt::runtime(notation), name)};
366 int count = 0;
367 for (const auto& val : data) {
368 col.push_back(fmt::format(fmt::runtime(notation), val));
369 count++;
370 if (count == dots) {
371 col.push_back(fmt::format(fmt::runtime(notation), ".."));
372 }
373 }
374 return col;
375}
376
377vector<string> stringColumn(string name, const vector<string>& comp,
378 int rows, int width)
379{
380 // extract data for processing
381 vector<string> data;
382 string notation = fmt::format("{{:{}}}", width);
383 size_t maxLen = 3; // minimum column width is 3
384 int csize = static_cast<int>(comp.size());
385 int dots = csize + 1;
386 if (csize <= rows) {
387 for (const auto& val : comp) {
388 data.push_back(val);
389 maxLen = std::max(maxLen,
390 boost::trim_copy(fmt::format(fmt::runtime(notation), val)).size());
391 }
392 } else {
393 dots = (rows + 1) / 2;
394 for (int row = 0; row < dots; row++) {
395 data.push_back(comp[row]);
396 maxLen = std::max(maxLen,
397 boost::trim_copy(
398 fmt::format(fmt::runtime(notation), comp[row])).size());
399 }
400 for (int row = csize - rows / 2; row < csize; row++) {
401 data.push_back(comp[row]);
402 maxLen = std::max(maxLen,
403 boost::trim_copy(
404 fmt::format(fmt::runtime(notation), comp[row])).size());
405 }
406 }
407
408 // assemble output
409 notation = fmt::format(" {{:>{}}}", maxLen);
410 vector<string> col = {fmt::format(fmt::runtime(notation), name)};
411 int count = 0;
412 for (const auto& val : data) {
413 col.push_back(fmt::format(fmt::runtime(notation), val));
414 count++;
415 if (count == dots) {
416 col.push_back(fmt::format(fmt::runtime(notation), "..."));
417 }
418 }
419 return col;
420}
421
422vector<string> formatColumn(string name, const AnyValue& comp, int rows, int width)
423{
424 if (comp.isVector<double>()) {
425 return doubleColumn(name, comp.asVector<double>(), rows, width);
426 }
427 if (comp.isVector<long int>()) {
428 return integerColumn(name, comp.asVector<long int>(), rows, width);
429 }
430 if (comp.isVector<string>()) {
431 return stringColumn(name, comp.asVector<string>(), rows, width);
432 }
433
434 // create alternative representation
435 string repr;
436 int size;
437 if (comp.isVector<vector<double>>()) {
438 repr = "[ <double> ]";
439 size = len(comp.asVector<vector<double>>());
440 } else if (comp.isVector<vector<long int>>()) {
441 repr = "[ <long int> ]";
442 size = len(comp.asVector<vector<long int>>());
443 } else if (comp.isVector<vector<string>>()) {
444 repr = "[ <string> ]";
445 size = len(comp.asVector<vector<string>>());
446 } else {
447 throw CanteraError(
448 "formatColumn", "Encountered invalid data for component '{}'.", name);
449 }
450 size_t maxLen = std::max(repr.size(), name.size());
451
452 // assemble output
453 string notation = fmt::format(" {{:>{}}}", maxLen);
454 repr = fmt::format(fmt::runtime(notation), repr);
455 vector<string> col = {fmt::format(fmt::runtime(notation), name)};
456 if (size <= rows) {
457 for (int row = 0; row < size; row++) {
458 col.push_back(repr);
459 }
460 } else {
461 int dots = (rows + 1) / 2;
462 for (int row = 0; row < dots; row++) {
463 col.push_back(repr);
464 }
465 col.push_back(fmt::format(fmt::runtime(notation), "..."));
466 for (int row = size - rows / 2; row < size; row++) {
467 col.push_back(repr);
468 }
469 }
470 return col;
471}
472
473} // end unnamed namespace
474
475string SolutionArray::info(const vector<string>& keys, int rows, int width)
476{
477 fmt::memory_buffer b;
478 int col_width = 12; // targeted maximum column width
479 vector<string> components;
480 if (keys.size()) {
481 components = keys;
482 } else {
483 components = componentNames();
484 }
485 try {
486 // build columns
487 vector<long int> index;
488 for (const auto ix : m_active) {
489 index.push_back(ix);
490 }
491 vector<vector<string>> cols = {integerColumn("", index, rows, col_width)};
492 vector<vector<string>> tail; // trailing columns in reverse order
493 size_t size = cols.back().size();
494
495 // assemble columns fitting within a maximum width; if this width is exceeded,
496 // a "..." separator is inserted close to the center. Accordingly, the matrix
497 // needs to be assembled from two halves.
498 int front = 0;
499 int back = len(components) - 1;
500 int fLen = len(cols.back()[0]);
501 int bLen = 0;
502 int sep = 5; // separator width
503 bool done = false;
504 while (!done && front <= back) {
505 string key;
506 while (bLen + sep <= fLen && front <= back) {
507 // add trailing columns
508 key = components[back];
509 auto col = formatColumn(key, getComponent(key), rows, col_width);
510 if (len(col[0]) + fLen + bLen + sep > width) {
511 done = true;
512 break;
513 }
514 tail.push_back(col);
515 bLen += len(tail.back()[0]);
516 back--;
517 }
518 if (done || front > back) {
519 break;
520 }
521 while (fLen <= bLen + sep && front <= back) {
522 // add leading columns
523 key = components[front];
524 auto col = formatColumn(key, getComponent(key), rows, col_width);
525 if (len(col[0]) + fLen + bLen + sep > width) {
526 done = true;
527 break;
528 }
529 cols.push_back(col);
530 fLen += len(cols.back()[0]);
531 front++;
532 }
533 }
534 if (cols.size() + tail.size() < components.size() + 1) {
535 // add separator
536 cols.push_back(vector<string>(size + 1, " ..."));
537 }
538 // copy trailing columns
539 cols.insert(cols.end(), tail.rbegin(), tail.rend());
540
541 // assemble formatted output
542 for (size_t row = 0; row < size; row++) {
543 for (const auto& col : cols) {
544 fmt_append(b, col[row]);
545 }
546 fmt_append(b, "\n");
547 }
548
549 // add size information
550 fmt_append(b, "\n[{} rows x {} components; state='{}']",
551 m_size, components.size(), m_sol->thermo()->nativeMode());
552
553 } catch (CanteraError& err) {
554 return to_string(b) + err.what();
555 }
556 return to_string(b);
557}
558
559shared_ptr<ThermoPhase> SolutionArray::thermo()
560{
561 return m_sol->thermo();
562}
563
564vector<string> SolutionArray::componentNames() const
565{
566 vector<string> components;
567 // leading auxiliary components
568 int pos = 0;
569 while (m_order->count(pos)) {
570 components.push_back(m_order->at(pos));
571 pos++;
572 }
573
574 // state information
575 auto phase = m_sol->thermo();
576 for (auto code : phase->nativeMode()) {
577 string name = string(1, code);
578 if (name == "X" || name == "Y") {
579 for (auto& spc : phase->speciesNames()) {
580 components.push_back(spc);
581 }
582 } else {
583 components.push_back(name);
584 }
585 }
586
587 // trailing auxiliary components
588 pos = -1;
589 while (m_order->count(pos)) {
590 components.push_back(m_order->at(pos));
591 pos--;
592 }
593
594 return components;
595}
596
597void SolutionArray::addExtra(const string& name, bool back)
598{
599 if (m_extra->count(name)) {
600 throw CanteraError("SolutionArray::addExtra",
601 "Component '{}' already exists.", name);
602 }
603 (*m_extra)[name] = AnyValue();
604 if (back) {
605 if (m_order->size()) {
606 // add name at end of back components
607 m_order->emplace(m_order->begin()->first - 1, name);
608 } else {
609 // first name after state components
610 m_order->emplace(-1, name);
611 }
612 } else {
613 if (m_order->size()) {
614 // insert name at end of front components
615 m_order->emplace(m_order->rbegin()->first + 1, name);
616 } else {
617 // name in leading position
618 m_order->emplace(0, name);
619 }
620 }
621}
622
623vector<string> SolutionArray::listExtra(bool all) const
624{
625 vector<string> names;
626 int pos = 0;
627 while (m_order->count(pos)) {
628 const auto& name = m_order->at(pos);
629 if (all || !m_extra->at(name).is<void>()) {
630 names.push_back(name);
631 }
632 pos++;
633 }
634
635 // trailing auxiliary components
636 pos = -1;
637 while (m_order->count(pos)) {
638 const auto& name = m_order->at(pos);
639 if (all || !m_extra->at(name).is<void>()) {
640 names.push_back(name);
641 }
642 pos--;
643 }
644 return names;
645}
646
647bool SolutionArray::hasComponent(const string& name) const
648{
649 if (m_extra->count(name)) {
650 // auxiliary data
651 return true;
652 }
653 if (m_sol->thermo()->speciesIndex(name) != npos) {
654 // species
655 return true;
656 }
657 if (name == "X" || name == "Y") {
658 // reserved names
659 return false;
660 }
661 // native state
662 return (m_sol->thermo()->nativeState().count(name));
663}
664
665AnyValue SolutionArray::getComponent(const string& name) const
666{
667 if (!hasComponent(name)) {
668 throw CanteraError("SolutionArray::getComponent",
669 "Unknown component '{}'.", name);
670 }
671
672 AnyValue out;
673 if (m_extra->count(name)) {
674 // extra component
675 const auto& extra = m_extra->at(name);
676 if (extra.is<void>()) {
677 return AnyValue();
678 }
679 if (m_size == m_dataSize) {
680 return extra; // slicing not necessary
681 }
682 if (extra.isVector<long int>()) {
683 return getSingle<long int>(extra, m_active);
684 }
685 if (extra.isVector<double>()) {
686 return getSingle<double>(extra, m_active);
687 }
688 if (extra.isVector<string>()) {
689 return getSingle<string>(extra, m_active);
690 }
691 if (extra.isVector<vector<double>>()) {
692 return getMulti<double>(extra, m_active);
693 }
694 if (extra.isVector<vector<long int>>()) {
695 return getMulti<long int>(extra, m_active);
696 }
697 if (extra.isVector<vector<string>>()) {
698 return getMulti<string>(extra, m_active);
699 }
700 throw NotImplementedError("SolutionArray::getComponent",
701 "Unable to get sliced data for component '{}' with type '{}'.",
702 name, extra.type_str());
703 }
704
705 // component is part of state information
706 vector<double> data(m_size);
707 size_t ix = m_sol->thermo()->speciesIndex(name);
708 if (ix == npos) {
709 // state other than species
710 ix = m_sol->thermo()->nativeState()[name];
711 } else {
712 // species information
713 ix += m_stride - m_sol->thermo()->nSpecies();
714 }
715 for (size_t k = 0; k < m_size; ++k) {
716 data[k] = (*m_data)[m_active[k] * m_stride + ix];
717 }
718 out = data;
719 return out;
720}
721
722bool isSimpleVector(const AnyValue& any) {
723 return any.isVector<double>() || any.isVector<long int>() ||
724 any.isVector<string>() || any.isVector<bool>() ||
725 any.isVector<vector<double>>() || any.isVector<vector<long int>>() ||
726 any.isVector<vector<string>>() || any.isVector<vector<bool>>();
727}
728
729void SolutionArray::setComponent(const string& name, const AnyValue& data)
730{
731 if (!hasComponent(name)) {
732 throw CanteraError("SolutionArray::setComponent",
733 "Unknown component '{}'.", name);
734 }
735 if (m_extra->count(name)) {
736 _setExtra(name, data);
737 return;
738 }
739
740 size_t size = data.vectorSize();
741 if (size == npos) {
742 throw CanteraError("SolutionArray::setComponent",
743 "Invalid type of component '{}': expected simple array type, "
744 "but received '{}'.", name, data.type_str());
745 }
746 if (size != m_size) {
747 throw CanteraError("SolutionArray::setComponent",
748 "Invalid size of component '{}': expected size {} but received {}.",
749 name, m_size, size);
750 }
751
752 auto& vec = data.asVector<double>();
753 size_t ix = m_sol->thermo()->speciesIndex(name);
754 if (ix == npos) {
755 ix = m_sol->thermo()->nativeState()[name];
756 } else {
757 ix += m_stride - m_sol->thermo()->nSpecies();
758 }
759 for (size_t k = 0; k < m_size; ++k) {
760 (*m_data)[m_active[k] * m_stride + ix] = vec[k];
761 }
762}
763
764void SolutionArray::setLoc(int loc, bool restore)
765{
766 size_t loc_ = static_cast<size_t>(loc);
767 if (m_size == 0) {
768 throw CanteraError("SolutionArray::setLoc",
769 "Unable to set location in empty SolutionArray.");
770 } else if (loc < 0) {
771 if (m_loc == npos) {
772 throw CanteraError("SolutionArray::setLoc",
773 "Both current and buffered indices are invalid.");
774 }
775 return;
776 } else if (static_cast<size_t>(m_active[loc_]) == m_loc) {
777 return;
778 } else if (loc_ >= m_size) {
779 throw IndexError("SolutionArray::setLoc", "indices", loc_, m_size - 1);
780 }
781 m_loc = static_cast<size_t>(m_active[loc_]);
782 if (restore) {
783 size_t nState = m_sol->thermo()->stateSize();
784 m_sol->thermo()->restoreState(nState, m_data->data() + m_loc * m_stride);
785 }
786}
787
788void SolutionArray::updateState(int loc)
789{
790 setLoc(loc, false);
791 size_t nState = m_sol->thermo()->stateSize();
792 m_sol->thermo()->saveState(nState, m_data->data() + m_loc * m_stride);
793}
794
795vector<double> SolutionArray::getState(int loc)
796{
797 setLoc(loc);
798 size_t nState = m_sol->thermo()->stateSize();
799 vector<double> out(nState);
800 m_sol->thermo()->saveState(out); // thermo contains current state
801 return out;
802}
803
804void SolutionArray::setState(int loc, const vector<double>& state)
805{
806 size_t nState = m_sol->thermo()->stateSize();
807 if (state.size() != nState) {
808 throw CanteraError("SolutionArray::setState",
809 "Expected array to have length {}, but received an array of length {}.",
810 nState, state.size());
811 }
812 setLoc(loc, false);
813 m_sol->thermo()->restoreState(state);
814 m_sol->thermo()->saveState(nState, m_data->data() + m_loc * m_stride);
815}
816
817void SolutionArray::normalize() {
818 auto phase = m_sol->thermo();
819 auto nativeState = phase->nativeState();
820 if (nativeState.size() < 3) {
821 return;
822 }
823 size_t nState = phase->stateSize();
824 vector<double> out(nState);
825 if (nativeState.count("Y")) {
826 size_t offset = nativeState["Y"];
827 for (int loc = 0; loc < static_cast<int>(m_size); loc++) {
828 setLoc(loc, true); // set location and restore state
829 phase->setMassFractions(m_data->data() + m_loc * m_stride + offset);
830 m_sol->thermo()->saveState(out);
831 setState(loc, out);
832 }
833 } else if (nativeState.count("X")) {
834 size_t offset = nativeState["X"];
835 for (int loc = 0; loc < static_cast<int>(m_size); loc++) {
836 setLoc(loc, true); // set location and restore state
837 phase->setMoleFractions(m_data->data() + m_loc * m_stride + offset);
838 m_sol->thermo()->saveState(out);
839 setState(loc, out);
840 }
841 } else {
842 throw NotImplementedError("SolutionArray::normalize",
843 "Not implemented for mode '{}'.", phase->nativeMode());
844 }
845}
846
847AnyMap SolutionArray::getAuxiliary(int loc)
848{
849 setLoc(loc);
850 AnyMap out;
851 for (const auto& [key, extra] : *m_extra) {
852 if (extra.is<void>()) {
853 out[key] = extra;
854 } else if (extra.isVector<long int>()) {
855 out[key] = extra.asVector<long int>()[m_loc];
856 } else if (extra.isVector<double>()) {
857 out[key] = extra.asVector<double>()[m_loc];
858 } else if (extra.isVector<string>()) {
859 out[key] = extra.asVector<string>()[m_loc];
860 } else if (extra.isVector<vector<long int>>()) {
861 out[key] = extra.asVector<vector<long int>>()[m_loc];
862 } else if (extra.isVector<vector<double>>()) {
863 out[key] = extra.asVector<vector<double>>()[m_loc];
864 } else if (extra.isVector<vector<string>>()) {
865 out[key] = extra.asVector<vector<string>>()[m_loc];
866 } else {
867 throw NotImplementedError("SolutionArray::getAuxiliary",
868 "Unable to retrieve data for component '{}' with type '{}'.",
869 key, extra.type_str());
870 }
871 }
872 return out;
873}
874
875void SolutionArray::setAuxiliary(int loc, const AnyMap& data)
876{
877 setLoc(loc, false);
878 for (const auto& [name, value] : data) {
879 if (!m_extra->count(name)) {
880 throw CanteraError("SolutionArray::setAuxiliary",
881 "Unknown auxiliary component '{}'.", name);
882 }
883 auto& extra = m_extra->at(name);
884 if (extra.is<void>()) {
885 if (m_dataSize > 1) {
886 throw CanteraError("SolutionArray::setAuxiliary",
887 "Unable to set location for type '{}': "
888 "component is not initialized.", name);
889 }
890 _initExtra(name, value);
891 _resizeExtra(name);
892 }
893 try {
894 if (extra.isVector<long int>()) {
895 setAuxiliarySingle<long int>(m_loc, extra, value);
896 } else if (extra.isVector<double>()) {
897 setAuxiliarySingle<double>(m_loc, extra, value);
898 } else if (extra.isVector<string>()) {
899 setAuxiliarySingle<string>(m_loc, extra, value);
900 } else if (extra.isVector<vector<long int>>()) {
901 setAuxiliaryMulti<long int>(m_loc, extra, value);
902 } else if (extra.isVector<vector<double>>()) {
903 setAuxiliaryMulti<double>(m_loc, extra, value);
904 } else if (extra.isVector<vector<string>>()) {
905 setAuxiliaryMulti<string>(m_loc, extra, value);
906 } else {
907 throw CanteraError("SolutionArray::setAuxiliary",
908 "Unable to set entry for type '{}'.", extra.type_str());
909 }
910 } catch (CanteraError& err) {
911 // make failed type conversions traceable
912 throw CanteraError("SolutionArray::setAuxiliary",
913 "Encountered incompatible value for component '{}':\n{}",
914 name, err.getMessage());
915 }
916 }
917}
918
919AnyMap preamble(const string& desc)
920{
921 AnyMap data;
922 if (desc.size()) {
923 data["description"] = desc;
924 }
925 data["generator"] = "Cantera SolutionArray";
926 data["cantera-version"] = CANTERA_VERSION;
927 data["git-commit"] = gitCommit();
928
929 // Add a timestamp indicating the current time
930 time_t aclock;
931 ::time(&aclock); // Get time in seconds
932 struct tm* newtime = localtime(&aclock); // Convert time to struct tm form
933 data["date"] = stripnonprint(asctime(newtime));
934
935 return data;
936}
937
938AnyMap& openField(AnyMap& root, const string& name)
939{
940 if (!name.size()) {
941 return root;
942 }
943
944 // locate field based on 'name'
945 vector<string> tokens;
946 tokenizePath(name, tokens);
947 AnyMap* ptr = &root; // use raw pointer to avoid copying
948 string path = "";
949 for (auto& field : tokens) {
950 path += "/" + field;
951 AnyMap& sub = *ptr;
952 if (sub.hasKey(field) && !sub[field].is<AnyMap>()) {
953 throw CanteraError("openField",
954 "Encountered invalid existing field '{}'.", path);
955 } else if (!sub.hasKey(field)) {
956 sub[field] = AnyMap();
957 }
958 ptr = &sub[field].as<AnyMap>();
959 }
960 return *ptr;
961}
962
963void SolutionArray::writeHeader(const string& fname, const string& name,
964 const string& desc, bool overwrite)
965{
966 Storage file(fname, true);
967 if (file.checkGroup(name, true)) {
968 if (!overwrite) {
969 throw CanteraError("SolutionArray::writeHeader",
970 "Group name '{}' exists; use 'overwrite' argument to overwrite.", name);
971 }
972 file.deleteGroup(name);
973 file.checkGroup(name, true);
974 }
975 file.writeAttributes(name, preamble(desc));
976}
977
978void SolutionArray::writeHeader(AnyMap& root, const string& name,
979 const string& desc, bool overwrite)
980{
981 AnyMap& data = openField(root, name);
982 if (!data.empty()) {
983 if (!overwrite) {
984 throw CanteraError("SolutionArray::writeHeader",
985 "Field name '{}' exists; use 'overwrite' argument to overwrite.", name);
986 }
987 data.clear();
988 }
989 data.update(preamble(desc));
990}
991
992void SolutionArray::writeEntry(const string& fname, bool overwrite, const string& basis)
993{
994 if (apiNdim() != 1) {
995 throw CanteraError("SolutionArray::writeEntry",
996 "Tabular output of CSV data only works for 1D SolutionArray objects.");
997 }
998 set<string> speciesNames;
999 for (const auto& species : m_sol->thermo()->speciesNames()) {
1000 speciesNames.insert(species);
1001 }
1002 bool mole;
1003 if (basis == "") {
1004 const auto& nativeState = m_sol->thermo()->nativeState();
1005 mole = nativeState.find("X") != nativeState.end();
1006 } else if (basis == "X" || basis == "mole") {
1007 mole = true;
1008 } else if (basis == "Y" || basis == "mass") {
1009 mole = false;
1010 } else {
1011 throw CanteraError("SolutionArray::writeEntry",
1012 "Invalid species basis '{}'.", basis);
1013 }
1014
1015 auto names = componentNames();
1016 size_t last = names.size() - 1;
1017 vector<AnyValue> components;
1018 vector<bool> isSpecies;
1019 fmt::memory_buffer header;
1020 for (const auto& key : names) {
1021 string label = key;
1022 size_t col;
1023 if (speciesNames.find(key) == speciesNames.end()) {
1024 // Pre-read component vectors
1025 isSpecies.push_back(false);
1026 components.emplace_back(getComponent(key));
1027 col = components.size() - 1;
1028 if (!components[col].isVector<double>() &&
1029 !components[col].isVector<long int>() &&
1030 !components[col].isVector<string>())
1031 {
1032 throw CanteraError("SolutionArray::writeEntry",
1033 "Multi-dimensional column '{}' is not supported for CSV output.",
1034 key);
1035 }
1036 } else {
1037 // Delay reading species data as base can be either mole or mass
1038 isSpecies.push_back(true);
1039 components.emplace_back(AnyValue());
1040 col = components.size() - 1;
1041 if (mole) {
1042 label = "X_" + label;
1043 } else {
1044 label = "Y_" + label;
1045 }
1046 }
1047 if (label.find("\"") != string::npos || label.find("\n") != string::npos) {
1048 throw NotImplementedError("SolutionArray::writeEntry",
1049 "Detected column name containing double quotes or line feeds: '{}'.",
1050 label);
1051 }
1052 string sep = (col == last) ? "" : ",";
1053 if (label.find(",") != string::npos) {
1054 fmt_append(header, "\"{}\"{}", label, sep);
1055 } else {
1056 fmt_append(header, "{}{}", label, sep);
1057 }
1058 }
1059
1060 // (Most) potential exceptions have been thrown; start writing data to file
1061 if (std::ifstream(fname).good()) {
1062 if (!overwrite) {
1063 throw CanteraError("SolutionArray::writeEntry",
1064 "File '{}' already exists; use option 'overwrite' to replace CSV file.",
1065 fname);
1066 }
1067 std::remove(fname.c_str());
1068 }
1069 std::ofstream output(fname);
1070 output << to_string(header) << std::endl;
1071
1072 size_t maxLen = npos;
1073 vector<double> buf(speciesNames.size(), 0.);
1074 for (int row = 0; row < static_cast<int>(m_size); row++) {
1075 fmt::memory_buffer line;
1076 if (maxLen != npos) {
1077 line.reserve(maxLen);
1078 }
1079 setLoc(row);
1080 if (mole) {
1081 m_sol->thermo()->getMoleFractions(buf.data());
1082 } else {
1083 m_sol->thermo()->getMassFractions(buf.data());
1084 }
1085
1086 size_t idx = 0;
1087 for (size_t col = 0; col < components.size(); col++) {
1088 string sep = (col == last) ? "" : ",";
1089 if (isSpecies[col]) {
1090 fmt_append(line, "{:.9g}{}", buf[idx++], sep);
1091 } else {
1092 auto& data = components[col];
1093 if (data.isVector<double>()) {
1094 fmt_append(line, "{:.9g}{}", data.asVector<double>()[row], sep);
1095 } else if (data.isVector<long int>()) {
1096 fmt_append(line, "{}{}", data.asVector<long int>()[row], sep);
1097 } else if (data.isVector<bool>()) {
1098 fmt_append(line, "{}{}",
1099 static_cast<bool>(data.asVector<bool>()[row]), sep);
1100 } else {
1101 auto value = data.asVector<string>()[row];
1102 if (value.find("\"") != string::npos ||
1103 value.find("\n") != string::npos)
1104 {
1105 throw NotImplementedError("SolutionArray::writeEntry",
1106 "Detected value containing double quotes or line feeds: "
1107 "'{}'", value);
1108 }
1109 if (value.find(",") != string::npos) {
1110 fmt_append(line, "\"{}\"{}", value, sep);
1111 } else {
1112 fmt_append(line, "{}{}", value, sep);
1113 }
1114 }
1115 }
1116 }
1117 output << to_string(line) << std::endl;
1118 maxLen = std::max(maxLen, line.size());
1119 }
1120}
1121
1122void SolutionArray::writeEntry(const string& fname, const string& name,
1123 const string& sub, bool overwrite, int compression)
1124{
1125 if (name == "") {
1126 throw CanteraError("SolutionArray::writeEntry",
1127 "Group name specifying root location must not be empty.");
1128 }
1129 if (m_size < m_dataSize) {
1130 throw NotImplementedError("SolutionArray::writeEntry",
1131 "Unable to save sliced data.");
1132 }
1133 Storage file(fname, true);
1134 if (compression) {
1135 file.setCompressionLevel(compression);
1136 }
1137 string path = name;
1138 if (sub != "") {
1139 path += "/" + sub;
1140 } else {
1141 path += "/data";
1142 }
1143 if (file.checkGroup(path, true)) {
1144 if (!overwrite) {
1145 throw CanteraError("SolutionArray::writeEntry",
1146 "Group name '{}' exists; use 'overwrite' argument to overwrite.", name);
1147 }
1148 file.deleteGroup(path);
1149 file.checkGroup(path, true);
1150 }
1151 file.writeAttributes(path, m_meta);
1152 AnyMap more;
1153 if (apiNdim() == 1) {
1154 more["size"] = int(m_dataSize);
1155 } else {
1156 more["api-shape"] = m_apiShape;
1157 }
1158 more["components"] = componentNames();
1159 file.writeAttributes(path, more);
1160 if (!m_dataSize) {
1161 return;
1162 }
1163
1164 const auto& nativeState = m_sol->thermo()->nativeState();
1165 size_t nSpecies = m_sol->thermo()->nSpecies();
1166 for (auto& [key, offset] : nativeState) {
1167 if (key == "X" || key == "Y") {
1168 vector<vector<double>> prop;
1169 for (size_t i = 0; i < m_size; i++) {
1170 size_t first = offset + i * m_stride;
1171 prop.emplace_back(m_data->begin() + first,
1172 m_data->begin() + first + nSpecies);
1173 }
1174 AnyValue data;
1175 data = prop;
1176 file.writeData(path, key, data);
1177 } else {
1178 auto data = getComponent(key);
1179 file.writeData(path, key, data);
1180 }
1181 }
1182
1183 for (const auto& [key, value] : *m_extra) {
1184 if (isSimpleVector(value)) {
1185 file.writeData(path, key, value);
1186 } else if (value.is<void>()) {
1187 // skip unintialized component
1188 } else {
1189 throw NotImplementedError("SolutionArray::writeEntry",
1190 "Unable to save component '{}' with data type {}.",
1191 key, value.type_str());
1192 }
1193 }
1194}
1195
1196void SolutionArray::writeEntry(AnyMap& root, const string& name, const string& sub,
1197 bool overwrite)
1198{
1199 if (name == "") {
1200 throw CanteraError("SolutionArray::writeEntry",
1201 "Field name specifying root location must not be empty.");
1202 }
1203 if (m_size < m_dataSize) {
1204 throw NotImplementedError("SolutionArray::writeEntry",
1205 "Unable to save sliced data.");
1206 }
1207 string path = name;
1208 if (sub != "") {
1209 path += "/" + sub;
1210 } else {
1211 path += "/data";
1212 }
1213 AnyMap& data = openField(root, path);
1214 bool preexisting = !data.empty();
1215 if (preexisting && !overwrite) {
1216 throw CanteraError("SolutionArray::writeEntry",
1217 "Field name '{}' exists; use 'overwrite' argument to overwrite.", name);
1218 }
1219 if (apiNdim() == 1) {
1220 data["size"] = int(m_dataSize);
1221 } else {
1222 data["api-shape"] = m_apiShape;
1223 }
1224 data.update(m_meta);
1225
1226 if (m_size > 1) {
1227 data["components"] = componentNames();
1228 }
1229
1230 for (auto& [_, key] : *m_order) {
1231 data[key] = m_extra->at(key);
1232 }
1233
1234 auto phase = m_sol->thermo();
1235 if (m_size == 1) {
1236 setLoc(0);
1237 data["temperature"] = phase->temperature();
1238 data["pressure"] = phase->pressure();
1239 auto surf = std::dynamic_pointer_cast<SurfPhase>(phase);
1240 auto nSpecies = phase->nSpecies();
1241 vector<double> values(nSpecies);
1242 if (surf) {
1243 surf->getCoverages(&values[0]);
1244 } else {
1245 phase->getMassFractions(&values[0]);
1246 }
1247 AnyMap items;
1248 for (size_t k = 0; k < nSpecies; k++) {
1249 if (values[k] != 0.0) {
1250 items[phase->speciesName(k)] = values[k];
1251 }
1252 }
1253 if (surf) {
1254 data["coverages"] = std::move(items);
1255 } else {
1256 data["mass-fractions"] = std::move(items);
1257 }
1258 } else if (m_size > 1) {
1259 for (auto& code : phase->nativeMode()) {
1260 string name(1, code);
1261 if (name == "X" || name == "Y") {
1262 for (auto& spc : phase->speciesNames()) {
1263 data[spc] = getComponent(spc);
1264 }
1265 data["basis"] = name == "X" ? "mole" : "mass";
1266 } else {
1267 data[name] = getComponent(name);
1268 }
1269 }
1270 }
1271
1272 static bool reg = AnyMap::addOrderingRules("SolutionArray",
1273 {{"head", "type"}, {"head", "size"}, {"head", "basis"}});
1274 if (reg) {
1275 data["__type__"] = "SolutionArray";
1276 }
1277
1278 // If this is not replacing an existing solution, put it at the end
1279 if (!preexisting) {
1280 data.setLoc(INT_MAX, 0);
1281 }
1282}
1283
1284void SolutionArray::append(const vector<double>& state, const AnyMap& extra)
1285{
1286 if (apiNdim() > 1) {
1287 throw NotImplementedError("SolutionArray::append",
1288 "Unable to append multi-dimensional arrays.");
1289 }
1290
1291 int pos = size();
1292 resize(pos + 1);
1293 try {
1294 setState(pos, state);
1295 setAuxiliary(pos, extra);
1296 } catch (CanteraError& err) {
1297 // undo resize and rethrow error
1298 resize(pos);
1299 throw CanteraError("SolutionArray::append", err.getMessage());
1300 }
1301}
1302
1303void SolutionArray::save(const string& fname, const string& name, const string& sub,
1304 const string& desc, bool overwrite, int compression,
1305 const string& basis)
1306{
1307 if (m_size < m_dataSize) {
1308 throw NotImplementedError("SolutionArray::save",
1309 "Unable to save sliced data.");
1310 }
1311 size_t dot = fname.find_last_of(".");
1312 string extension = (dot != npos) ? toLowerCopy(fname.substr(dot + 1)) : "";
1313 if (extension == "csv") {
1314 if (name != "") {
1315 warn_user("SolutionArray::save",
1316 "Parameter 'name' not used for CSV output.");
1317 }
1318 writeEntry(fname, overwrite, basis);
1319 return;
1320 }
1321 if (basis != "") {
1322 warn_user("SolutionArray::save",
1323 "Argument 'basis' is not used for HDF or YAML output.", basis);
1324 }
1325 if (extension == "h5" || extension == "hdf" || extension == "hdf5") {
1326 writeHeader(fname, name, desc, overwrite);
1327 writeEntry(fname, name, sub, true, compression);
1328 return;
1329 }
1330 if (extension == "yaml" || extension == "yml") {
1331 // Check for an existing file and load it if present
1332 AnyMap data;
1333 if (std::ifstream(fname).good()) {
1334 data = AnyMap::fromYamlFile(fname);
1335 }
1336 writeHeader(data, name, desc, overwrite);
1337 writeEntry(data, name, sub, true);
1338
1339 // Write the output file and remove the now-outdated cached file
1340 std::ofstream out(fname);
1341 out << data.toYamlString();
1342 AnyMap::clearCachedFile(fname);
1343 return;
1344 }
1345 throw CanteraError("SolutionArray::save",
1346 "Unknown file extension '{}'.", extension);
1347}
1348
1349AnyMap SolutionArray::readHeader(const string& fname, const string& name)
1350{
1351 Storage file(fname, false);
1352 file.checkGroup(name);
1353 return file.readAttributes(name, false);
1354}
1355
1356const AnyMap& locateField(const AnyMap& root, const string& name)
1357{
1358 if (!name.size()) {
1359 return root;
1360 }
1361
1362 // locate field based on 'name'
1363 vector<string> tokens;
1364 tokenizePath(name, tokens);
1365 const AnyMap* ptr = &root; // use raw pointer to avoid copying
1366 string path = "";
1367 for (auto& field : tokens) {
1368 path += "/" + field;
1369 const AnyMap& sub = *ptr;
1370 if (!sub.hasKey(field) || !sub[field].is<AnyMap>()) {
1371 throw CanteraError("SolutionArray::locateField",
1372 "No field or solution with name '{}'.", path);
1373 }
1374 ptr = &sub[field].as<AnyMap>();
1375 }
1376 return *ptr;
1377}
1378
1379AnyMap SolutionArray::readHeader(const AnyMap& root, const string& name)
1380{
1381 auto sub = locateField(root, name);
1382 AnyMap header;
1383 for (const auto& [key, value] : sub) {
1384 if (!sub[key].is<AnyMap>()) {
1385 header[key] = value;
1386 }
1387 }
1388 return header;
1389}
1390
1391AnyMap SolutionArray::restore(const string& fname,
1392 const string& name, const string& sub)
1393{
1394 size_t dot = fname.find_last_of(".");
1395 string extension = (dot != npos) ? toLowerCopy(fname.substr(dot + 1)) : "";
1396 AnyMap header;
1397 if (extension == "csv") {
1398 throw NotImplementedError("SolutionArray::restore",
1399 "CSV import not implemented; if using Python, data can be imported via "
1400 "'read_csv' instead.");
1401 }
1402 if (extension == "h5" || extension == "hdf" || extension == "hdf5") {
1403 readEntry(fname, name, sub);
1404 header = readHeader(fname, name);
1405 } else if (extension == "yaml" || extension == "yml") {
1406 const AnyMap& root = AnyMap::fromYamlFile(fname);
1407 readEntry(root, name, sub);
1408 header = readHeader(root, name);
1409 } else {
1410 throw CanteraError("SolutionArray::restore",
1411 "Unknown file extension '{}'; supported extensions include "
1412 "'h5'/'hdf'/'hdf5' and 'yml'/'yaml'.", extension);
1413 }
1414 return header;
1415}
1416
1417void SolutionArray::_initExtra(const string& name, const AnyValue& value)
1418{
1419 if (!m_extra->count(name)) {
1420 throw CanteraError("SolutionArray::_initExtra",
1421 "Component '{}' does not exist.", name);
1422 }
1423 auto& extra = (*m_extra)[name];
1424 if (!extra.is<void>()) {
1425 throw CanteraError("SolutionArray::_initExtra",
1426 "Component '{}' is already initialized.", name);
1427 }
1428 try {
1429 if (value.is<long int>()) {
1430 extra = vector<long int>(m_dataSize, value.as<long int>());
1431 } else if (value.is<double>()) {
1432 extra = vector<double>(m_dataSize, value.as<double>());
1433 } else if (value.is<string>()) {
1434 extra = vector<string>(m_dataSize, value.as<string>());
1435 } else if (value.isVector<long int>()) {
1436 extra = vector<vector<long int>>(m_dataSize, value.asVector<long int>());
1437 } else if (value.isVector<double>()) {
1438 extra = vector<vector<double>>(m_dataSize, value.asVector<double>());
1439 } else if (value.isVector<string>()) {
1440 extra = vector<vector<string>>(m_dataSize, value.asVector<string>());
1441 } else if (value.is<void>()) {
1442 // placeholder for unknown type; settable in setComponent
1443 extra = AnyValue();
1444 } else {
1445 throw NotImplementedError("SolutionArray::_initExtra",
1446 "Unable to initialize component '{}' with type '{}'.",
1447 name, value.type_str());
1448 }
1449 } catch (CanteraError& err) {
1450 // make failed type conversions traceable
1451 throw CanteraError("SolutionArray::_initExtra",
1452 "Encountered incompatible value for initializing component '{}':\n{}",
1453 name, err.getMessage());
1454 }
1455}
1456
1457void SolutionArray::_resizeExtra(const string& name, const AnyValue& value)
1458{
1459 if (!m_extra->count(name)) {
1460 throw CanteraError("SolutionArray::_resizeExtra",
1461 "Component '{}' does not exist.", name);
1462 }
1463 auto& extra = (*m_extra)[name];
1464 if (extra.is<void>()) {
1465 // cannot resize placeholder (uninitialized component)
1466 return;
1467 }
1468
1469 try {
1470 if (extra.isVector<long int>()) {
1471 resizeSingle<long int>(extra, m_dataSize, value);
1472 } else if (extra.isVector<double>()) {
1473 resizeSingle<double>(extra, m_dataSize, value);
1474 } else if (extra.isVector<string>()) {
1475 resizeSingle<string>(extra, m_dataSize, value);
1476 } else if (extra.isVector<vector<double>>()) {
1477 resizeMulti<double>(extra, m_dataSize, value);
1478 } else if (extra.isVector<vector<long int>>()) {
1479 resizeMulti<long int>(extra, m_dataSize, value);
1480 } else if (extra.isVector<vector<string>>()) {
1481 resizeMulti<string>(extra, m_dataSize, value);
1482 } else {
1483 throw NotImplementedError("SolutionArray::_resizeExtra",
1484 "Unable to resize using type '{}'.", extra.type_str());
1485 }
1486 } catch (CanteraError& err) {
1487 // make failed type conversions traceable
1488 throw CanteraError("SolutionArray::_resizeExtra",
1489 "Encountered incompatible value for resizing component '{}':\n{}",
1490 name, err.getMessage());
1491 }
1492}
1493
1494void SolutionArray::_setExtra(const string& name, const AnyValue& data)
1495{
1496 if (!m_extra->count(name)) {
1497 throw CanteraError("SolutionArray::_setExtra",
1498 "Extra component '{}' does not exist.", name);
1499 }
1500
1501 auto& extra = m_extra->at(name);
1502 if (data.is<void>() && m_size == m_dataSize) {
1503 // reset placeholder
1504 extra = AnyValue();
1505 return;
1506 }
1507
1508 // uninitialized component
1509 if (extra.is<void>()) {
1510 if (m_size != m_dataSize) {
1511 throw CanteraError("SolutionArray::_setExtra",
1512 "Unable to replace '{}' for sliced data.", name);
1513 }
1514 if (data.vectorSize() == m_dataSize || data.matrixShape().first == m_dataSize) {
1515 // assign entire component
1516 extra = data;
1517 return;
1518 }
1519 // data contains default value
1520 if (m_dataSize == 0 && (data.isScalar() || data.vectorSize() > 0)) {
1521 throw CanteraError("SolutionArray::_setExtra",
1522 "Unable to initialize '{}' with non-empty values when SolutionArray is "
1523 "empty.", name);
1524 }
1525 if (m_dataSize && data.vectorSize() == 0) {
1526 throw CanteraError("SolutionArray::_setExtra",
1527 "Unable to initialize '{}' with empty array when SolutionArray is not "
1528 "empty." , name);
1529 }
1530 _initExtra(name, data);
1531 _resizeExtra(name, data);
1532 return;
1533 }
1534
1535 if (data.is<long int>()) {
1536 setScalar<long int>(extra, data, m_active);
1537 } else if (data.is<double>()) {
1538 setScalar<double>(extra, data, m_active);
1539 } else if (data.is<string>()) {
1540 setScalar<string>(extra, data, m_active);
1541 } else if (data.isVector<long int>()) {
1542 setSingle<long int>(extra, data, m_active);
1543 } else if (data.isVector<double>()) {
1544 setSingle<double>(extra, data, m_active);
1545 } else if (data.isVector<string>()) {
1546 setSingle<string>(extra, data, m_active);
1547 } else if (data.isVector<vector<long int>>()) {
1548 setMulti<long int>(extra, data, m_active);
1549 } else if (data.isVector<vector<double>>()) {
1550 setMulti<double>(extra, data, m_active);
1551 } else if (data.isVector<vector<string>>()) {
1552 setMulti<string>(extra, data, m_active);
1553 } else {
1554 throw NotImplementedError("SolutionArray::_setExtra",
1555 "Unable to set sliced data for component '{}' with type '{}'.",
1556 name, data.type_str());
1557 }
1558}
1559
1560string SolutionArray::_detectMode(const set<string>& names, bool native)
1561{
1562 // check set of available names against state acronyms defined by Phase::fullStates
1563 string mode = "";
1564 const auto& nativeState = m_sol->thermo()->nativeState();
1565 bool usesNativeState = false;
1566 auto surf = std::dynamic_pointer_cast<SurfPhase>(m_sol->thermo());
1567 bool found;
1568 for (const auto& item : m_sol->thermo()->fullStates()) {
1569 found = true;
1570 string name;
1571 usesNativeState = true;
1572 for (size_t i = 0; i < item.size(); i++) {
1573 // pick i-th letter from "full" state acronym
1574 name = string(1, item[i]);
1575 if (surf && (name == "X" || name == "Y")) {
1576 // override native state to enable detection of surface phases
1577 name = "C";
1578 usesNativeState = false;
1579 break;
1580 }
1581 if (names.count(name)) {
1582 // property is stored using letter acronym
1583 usesNativeState &= nativeState.count(name) > 0;
1584 } else if (aliasMap.count(name) && names.count(aliasMap.at(name))) {
1585 // property is stored using property name
1586 usesNativeState &= nativeState.count(name) > 0;
1587 } else {
1588 found = false;
1589 break;
1590 }
1591 }
1592 if (found) {
1593 mode = (name == "C") ? item.substr(0, 2) + "C" : item;
1594 break;
1595 }
1596 }
1597 if (!found) {
1598 if (surf && names.count("T") && names.count("X") && names.count("density")) {
1599 // Legacy HDF format erroneously uses density (Cantera < 3.0)
1600 return "legacySurf";
1601 }
1602 if (names.count("mass-flux") && names.count("mass-fractions")) {
1603 // Legacy YAML format stores incomplete state information (Cantera < 3.0)
1604 return "legacyInlet";
1605 }
1606 throw CanteraError("SolutionArray::_detectMode",
1607 "Detected incomplete thermodynamic information. Full states for a '{}' "
1608 "phase\nmay be defined by the following modes:\n'{}'\n"
1609 "Available data are: '{}'", m_sol->thermo()->type(),
1610 ba::join(m_sol->thermo()->fullStates(), "', '"), ba::join(names, "', '"));
1611 }
1612 if (usesNativeState && native) {
1613 return "native";
1614 }
1615 return mode;
1616}
1617
1618set<string> SolutionArray::_stateProperties(
1619 const string& mode, bool alias)
1620{
1621 set<string> states;
1622 if (mode == "native") {
1623 for (const auto& [name, offset] : m_sol->thermo()->nativeState()) {
1624 states.insert(alias ? aliasMap.at(name) : name);
1625 }
1626 } else {
1627 for (const auto& m : mode) {
1628 const string name = string(1, m);
1629 states.insert(alias ? aliasMap.at(name) : name);
1630 }
1631 }
1632
1633 return states;
1634}
1635
1636string getName(const set<string>& names, const string& name)
1637{
1638 if (names.count(name)) {
1639 return name;
1640 }
1641 const auto& alias = aliasMap.at(name);
1642 if (names.count(alias)) {
1643 return alias;
1644 }
1645 return name; // let exception be thrown elsewhere
1646}
1647
1648void SolutionArray::readEntry(const string& fname, const string& name,
1649 const string& sub)
1650{
1651 Storage file(fname, false);
1652 if (name == "") {
1653 throw CanteraError("SolutionArray::readEntry",
1654 "Group name specifying root location must not be empty.");
1655 }
1656 string path = name;
1657 if (sub != "" && file.checkGroup(name + "/" + sub, true)) {
1658 path += "/" + sub;
1659 } else if (sub == "" && file.checkGroup(name + "/data", true)) {
1660 // default data location
1661 path += "/data";
1662 }
1663 if (!file.checkGroup(path)) {
1664 throw CanteraError("SolutionArray::readEntry",
1665 "Group name specifying data entry is empty.");
1666 }
1667 m_extra->clear();
1668 auto [size, names] = file.contents(path);
1669 m_meta = file.readAttributes(path, true);
1670 if (m_meta.hasKey("size")) {
1671 // one-dimensional array
1672 resize(m_meta["size"].as<long int>());
1673 m_meta.erase("size");
1674 } else if (m_meta.hasKey("api-shape")) {
1675 // API uses multiple dimensions to interpret C++ SolutionArray
1676 setApiShape(m_meta["api-shape"].asVector<long int>());
1677 m_meta.erase("api-shape");
1678 } else {
1679 // legacy format; size needs to be detected
1680 resize(static_cast<int>(size));
1681 }
1682
1683 if (m_size == 0) {
1684 return;
1685 }
1686
1687 // determine storage mode of state data
1688 string mode = _detectMode(names);
1689 set<string> states = _stateProperties(mode);
1690 if (states.count("C")) {
1691 if (names.count("X")) {
1692 states.erase("C");
1693 states.insert("X");
1694 mode = mode.substr(0, 2) + "X";
1695 } else if (names.count("Y")) {
1696 states.erase("C");
1697 states.insert("Y");
1698 mode = mode.substr(0, 2) + "Y";
1699 }
1700 }
1701
1702 // restore state data
1703 size_t nSpecies = m_sol->thermo()->nSpecies();
1704 size_t nState = m_sol->thermo()->stateSize();
1705 const auto& nativeStates = m_sol->thermo()->nativeState();
1706 if (mode == "native") {
1707 // native state can be written directly into data storage
1708 for (const auto& [name, offset] : nativeStates) {
1709 if (name == "X" || name == "Y") {
1710 AnyValue data;
1711 data = file.readData(path, name, m_size, nSpecies);
1712 auto prop = data.asVector<vector<double>>();
1713 for (size_t i = 0; i < m_dataSize; i++) {
1714 std::copy(prop[i].begin(), prop[i].end(),
1715 m_data->data() + offset + i * m_stride);
1716 }
1717 } else {
1718 AnyValue data;
1719 data = file.readData(path, getName(names, name), m_dataSize, 0);
1720 setComponent(name, data);
1721 }
1722 }
1723 } else if (mode == "TPX") {
1724 AnyValue data;
1725 data = file.readData(path, getName(names, "T"), m_dataSize, 0);
1726 vector<double> T = std::move(data.asVector<double>());
1727 data = file.readData(path, getName(names, "P"), m_dataSize, 0);
1728 vector<double> P = std::move(data.asVector<double>());
1729 data = file.readData(path, "X", m_dataSize, nSpecies);
1730 vector<vector<double>> X = std::move(data.asVector<vector<double>>());
1731 for (size_t i = 0; i < m_dataSize; i++) {
1732 m_sol->thermo()->setMoleFractions_NoNorm(X[i].data());
1733 m_sol->thermo()->setState_TP(T[i], P[i]);
1734 m_sol->thermo()->saveState(nState, m_data->data() + i * m_stride);
1735 }
1736 } else if (mode == "TDX") {
1737 AnyValue data;
1738 data = file.readData(path, getName(names, "T"), m_dataSize, 0);
1739 vector<double> T = std::move(data.asVector<double>());
1740 data = file.readData(path, getName(names, "D"), m_dataSize, 0);
1741 vector<double> D = std::move(data.asVector<double>());
1742 data = file.readData(path, "X", m_dataSize, nSpecies);
1743 vector<vector<double>> X = std::move(data.asVector<vector<double>>());
1744 for (size_t i = 0; i < m_dataSize; i++) {
1745 m_sol->thermo()->setMoleFractions_NoNorm(X[i].data());
1746 m_sol->thermo()->setState_TD(T[i], D[i]);
1747 m_sol->thermo()->saveState(nState, m_data->data() + i * m_stride);
1748 }
1749 } else if (mode == "TPY") {
1750 AnyValue data;
1751 data = file.readData(path, getName(names, "T"), m_dataSize, 0);
1752 vector<double> T = std::move(data.asVector<double>());
1753 data = file.readData(path, getName(names, "P"), m_dataSize, 0);
1754 vector<double> P = std::move(data.asVector<double>());
1755 data = file.readData(path, "Y", m_dataSize, nSpecies);
1756 vector<vector<double>> Y = std::move(data.asVector<vector<double>>());
1757 for (size_t i = 0; i < m_dataSize; i++) {
1758 m_sol->thermo()->setMassFractions_NoNorm(Y[i].data());
1759 m_sol->thermo()->setState_TP(T[i], P[i]);
1760 m_sol->thermo()->saveState(nState, m_data->data() + i * m_stride);
1761 }
1762 } else if (mode == "legacySurf") {
1763 // erroneous TDX mode (should be TPX or TPY) - Sim1D (Cantera 2.5)
1764 AnyValue data;
1765 data = file.readData(path, getName(names, "T"), m_dataSize, 0);
1766 vector<double> T = std::move(data.asVector<double>());
1767 data = file.readData(path, "X", m_dataSize, nSpecies);
1768 vector<vector<double>> X = std::move(data.asVector<vector<double>>());
1769 for (size_t i = 0; i < m_dataSize; i++) {
1770 m_sol->thermo()->setMoleFractions_NoNorm(X[i].data());
1771 m_sol->thermo()->setTemperature(T[i]);
1772 m_sol->thermo()->saveState(nState, m_data->data() + i * m_stride);
1773 }
1774 warn_user("SolutionArray::readEntry",
1775 "Detected legacy HDF format with incomplete state information\nfor name "
1776 "'{}' (pressure missing).", path);
1777 } else if (mode == "") {
1778 throw CanteraError("SolutionArray::readEntry",
1779 "Data are not consistent with full state modes.");
1780 } else {
1781 throw NotImplementedError("SolutionArray::readEntry",
1782 "Import of '{}' data is not supported.", mode);
1783 }
1784
1785 // restore remaining data
1786 if (m_meta.hasKey("components")) {
1787 const auto& components = m_meta["components"].asVector<string>();
1788 bool back = false;
1789 for (const auto& name : components) {
1790 if (hasComponent(name) || name == "X" || name == "Y") {
1791 back = true;
1792 } else {
1793 addExtra(name, back);
1794 AnyValue data;
1795 data = file.readData(path, name, m_dataSize);
1796 setComponent(name, data);
1797 }
1798 }
1799 m_meta.erase("components");
1800 } else {
1801 // data format used by Python h5py export (Cantera 2.5)
1802 warn_user("SolutionArray::readEntry", "Detected legacy HDF format.");
1803 for (const auto& name : names) {
1804 if (!hasComponent(name) && name != "X" && name != "Y") {
1805 addExtra(name);
1806 AnyValue data;
1807 data = file.readData(path, name, m_dataSize);
1808 setComponent(name, data);
1809 }
1810 }
1811 }
1812}
1813
1814void SolutionArray::readEntry(const AnyMap& root, const string& name, const string& sub)
1815{
1816 if (name == "") {
1817 throw CanteraError("SolutionArray::readEntry",
1818 "Field name specifying root location must not be empty.");
1819 }
1820 auto path = locateField(root, name);
1821 if (path.hasKey("generator") && sub != "") {
1822 // read entry from subfolder (since Cantera 3.0)
1823 path = locateField(root, name + "/" + sub);
1824 } else if (sub == "" && path.hasKey("data")) {
1825 // default data location
1826 path = locateField(root, name + "/data");
1827 }
1828
1829 // set size and initialize
1830 long size = 0;
1831 if (path.hasKey("size")) {
1832 // one-dimensional array
1833 resize(path["size"].asInt());
1834 } else if (path.hasKey("api-shape")) {
1835 // API uses multiple dimensions to interpret C++ SolutionArray
1836 auto& shape = path["api-shape"].asVector<long int>();
1837 setApiShape(shape);
1838 } else {
1839 // legacy format (Cantera 2.6)
1840 size = path.getInt("points", 0);
1841 if (!path.hasKey("T") && !path.hasKey("temperature")) {
1842 // overwrite size - Sim1D erroneously assigns '1'
1843 size = (long)0;
1844 }
1845 resize(static_cast<int>(size));
1846 }
1847 m_extra->clear();
1848
1849 // restore data
1850 set<string> exclude = {"size", "api-shape", "points", "X", "Y"};
1851 set<string> names = path.keys();
1852 size_t nState = m_sol->thermo()->stateSize();
1853 if (m_dataSize == 0) {
1854 // no data points
1855 } else if (m_dataSize == 1) {
1856 // single data point
1857 string mode = _detectMode(names, false);
1858 if (mode == "TPY") {
1859 double T = path[getName(names, "T")].asDouble();
1860 double P = path[getName(names, "P")].asDouble();
1861 auto Y = path["mass-fractions"].asMap<double>();
1862 m_sol->thermo()->setState_TPY(T, P, Y);
1863 } else if (mode == "TPC") {
1864 auto surf = std::dynamic_pointer_cast<SurfPhase>(m_sol->thermo());
1865 double T = path[getName(names, "T")].asDouble();
1866 double P = path["pressure"].asDouble();
1867 m_sol->thermo()->setState_TP(T, P);
1868 auto cov = path["coverages"].asMap<double>();
1869 surf->setCoveragesByName(cov);
1870 } else if (mode == "legacyInlet") {
1871 // missing property - Sim1D (Cantera 2.6)
1872 mode = "TPY";
1873 double T = path[getName(names, "T")].asDouble();
1874 auto Y = path["mass-fractions"].asMap<double>();
1875 m_sol->thermo()->setState_TPY(T, m_sol->thermo()->pressure(), Y);
1876 warn_user("SolutionArray::readEntry",
1877 "Detected legacy YAML format with incomplete state information\n"
1878 "for name '{}' (pressure missing).", name + "/" + sub);
1879 } else if (mode == "") {
1880 throw CanteraError("SolutionArray::readEntry",
1881 "Data are not consistent with full state modes.");
1882 } else {
1883 throw NotImplementedError("SolutionArray::readEntry",
1884 "Import of '{}' data is not supported.", mode);
1885 }
1886 m_sol->thermo()->saveState(nState, m_data->data());
1887 auto props = _stateProperties(mode, true);
1888 exclude.insert(props.begin(), props.end());
1889 } else {
1890 // multiple data points
1891 if (path.hasKey("components")) {
1892 const auto& components = path["components"].asVector<string>();
1893 bool back = false;
1894 for (const auto& name : components) {
1895 if (hasComponent(name)) {
1896 back = true;
1897 } else {
1898 addExtra(name, back);
1899 }
1900 setComponent(name, path[name]);
1901 exclude.insert(name);
1902 }
1903 } else {
1904 // legacy YAML format does not provide for list of components
1905 for (const auto& [name, value] : path) {
1906 if (value.isVector<double>()) {
1907 const vector<double>& data = value.asVector<double>();
1908 if (data.size() == m_dataSize) {
1909 if (!hasComponent(name)) {
1910 addExtra(name);
1911 }
1912 setComponent(name, value);
1913 exclude.insert(name);
1914 }
1915 }
1916 }
1917 }
1918
1919 // check that state data are complete
1920 const auto& nativeState = m_sol->thermo()->nativeState();
1921 set<string> props;
1922 set<string> missingProps;
1923 for (const auto& [name, offset] : nativeState) {
1924 if (exclude.count(name)) {
1925 props.insert(name);
1926 } else {
1927 missingProps.insert(name);
1928 }
1929 }
1930
1931 set<string> TY = {"T", "Y"};
1932 if (props == TY && missingProps.count("D") && path.hasKey("pressure")) {
1933 // missing "D" - Sim1D (Cantera 2.6)
1934 double P = path["pressure"].asDouble();
1935 const size_t offset_T = nativeState.find("T")->second;
1936 const size_t offset_D = nativeState.find("D")->second;
1937 const size_t offset_Y = nativeState.find("Y")->second;
1938 for (size_t i = 0; i < m_dataSize; i++) {
1939 double T = (*m_data)[offset_T + i * m_stride];
1940 m_sol->thermo()->setState_TPY(
1941 T, P, m_data->data() + offset_Y + i * m_stride);
1942 (*m_data)[offset_D + i * m_stride] = m_sol->thermo()->density();
1943 }
1944 } else if (missingProps.size()) {
1945 throw CanteraError("SolutionArray::readEntry",
1946 "Incomplete state information: missing '{}'.",
1947 ba::join(missingProps, "', '"));
1948 }
1949 }
1950
1951 // add meta data
1952 for (const auto& [name, value] : path) {
1953 if (!exclude.count(name)) {
1954 m_meta[name] = value;
1955 }
1956 }
1957}
1958
1959namespace { // restrict scope of helper functions to local translation unit
1960
1961template<class T>
1962AnyValue getSingle(const AnyValue& extra, const vector<int>& slice)
1963{
1964 vector<T> data(slice.size());
1965 const auto& vec = extra.asVector<T>();
1966 for (size_t k = 0; k < slice.size(); ++k) {
1967 data[k] = vec[slice[k]];
1968 }
1969 AnyValue out;
1970 out = data;
1971 return out;
1972}
1973
1974template<class T>
1975AnyValue getMulti(const AnyValue& extra, const vector<int>& slice)
1976{
1977 vector<vector<T>> data(slice.size());
1978 const auto& vec = extra.asVector<vector<T>>();
1979 for (size_t k = 0; k < slice.size(); ++k) {
1980 data[k] = vec[slice[k]];
1981 }
1982 AnyValue out;
1983 out = data;
1984 return out;
1985}
1986
1987template<class T>
1988void setScalar(AnyValue& extra, const AnyValue& data, const vector<int>& slice)
1989{
1990 T value = data.as<T>();
1991 if (extra.isVector<T>()) {
1992 auto& vec = extra.asVector<T>();
1993 for (size_t k = 0; k < slice.size(); ++k) {
1994 vec[slice[k]] = value;
1995 }
1996 } else {
1997 throw CanteraError("SolutionArray::setScalar",
1998 "Incompatible input data: unable to assign '{}' data to '{}'.",
1999 data.type_str(), extra.type_str());
2000 }
2001}
2002
2003template<class T>
2004void setSingle(AnyValue& extra, const AnyValue& data, const vector<int>& slice)
2005{
2006 size_t size = slice.size();
2007 if (extra.vectorSize() == size && data.vectorSize() == size) {
2008 extra = data; // no slicing necessary; type can change
2009 return;
2010 }
2011 if (extra.matrixShape().first == size && data.vectorSize() == size) {
2012 extra = data; // no slicing necessary; shape and type can change
2013 return;
2014 }
2015 if (extra.type_str() != data.type_str()) {
2016 // do not allow changing of data type when slicing
2017 throw CanteraError("SolutionArray::setSingle",
2018 "Incompatible types: expected '{}' but received '{}'.",
2019 extra.type_str(), data.type_str());
2020 }
2021 const auto& vData = data.asVector<T>();
2022 if (vData.size() != size) {
2023 throw CanteraError("SolutionArray::setSingle",
2024 "Invalid input data size: expected {} entries but received {}.",
2025 size, vData.size());
2026 }
2027 auto& vec = extra.asVector<T>();
2028 for (size_t k = 0; k < size; ++k) {
2029 vec[slice[k]] = vData[k];
2030 }
2031}
2032
2033template<class T>
2034void setMulti(AnyValue& extra, const AnyValue& data, const vector<int>& slice)
2035{
2036 if (!data.isMatrix<T>()) {
2037 throw CanteraError("SolutionArray::setMulti",
2038 "Invalid input data shape: inconsistent number of columns.");
2039 }
2040 size_t size = slice.size();
2041 auto [rows, cols] = data.matrixShape();
2042 if (extra.matrixShape().first == size && rows == size) {
2043 extra = data; // no slicing necessary; type can change
2044 return;
2045 }
2046 if (extra.vectorSize() == size && rows == size) {
2047 extra = data; // no slicing necessary; shape and type can change
2048 return;
2049 }
2050 if (extra.type_str() != data.type_str()) {
2051 // do not allow changing of data type when slicing
2052 throw CanteraError("SolutionArray::setMulti",
2053 "Incompatible types: expected '{}' but received '{}'.",
2054 extra.type_str(), data.type_str());
2055 }
2056 if (rows != size) {
2057 throw CanteraError("SolutionArray::setMulti",
2058 "Invalid input data shape: expected {} rows but received {}.",
2059 size, rows);
2060 }
2061 if (extra.matrixShape().second != cols) {
2062 throw CanteraError("SolutionArray::setMulti",
2063 "Invalid input data shape: expected {} columns but received {}.",
2064 extra.matrixShape().second, cols);
2065 }
2066 const auto& vData = data.asVector<vector<T>>();
2067 auto& vec = extra.asVector<vector<T>>();
2068 for (size_t k = 0; k < slice.size(); ++k) {
2069 vec[slice[k]] = vData[k];
2070 }
2071}
2072
2073template<class T>
2074void resizeSingle(AnyValue& extra, size_t size, const AnyValue& value)
2075{
2076 T defaultValue;
2077 if (value.is<void>()) {
2078 defaultValue = vector<T>(1)[0];
2079 } else {
2080 defaultValue = value.as<T>();
2081 }
2082 extra.asVector<T>().resize(size, defaultValue);
2083}
2084
2085template<class T>
2086void resizeMulti(AnyValue& extra, size_t size, const AnyValue& value)
2087{
2088 vector<T> defaultValue;
2089 if (value.is<void>()) {
2090 defaultValue = vector<T>(extra.matrixShape().second);
2091 } else {
2092 defaultValue = value.as<vector<T>>();
2093 }
2094 extra.asVector<vector<T>>().resize(size, defaultValue);
2095}
2096
2097template<class T>
2098void resetSingle(AnyValue& extra, const vector<int>& slice)
2099{
2100 T defaultValue = vector<T>(1)[0];
2101 vector<T>& data = extra.asVector<T>();
2102 for (size_t k = 0; k < slice.size(); ++k) {
2103 data[slice[k]] = defaultValue;
2104 }
2105}
2106
2107template<class T>
2108void resetMulti(AnyValue& extra, const vector<int>& slice)
2109{
2110 vector<T> defaultValue = vector<T>(extra.matrixShape().second);
2111 vector<vector<T>>& data = extra.asVector<vector<T>>();
2112 for (size_t k = 0; k < slice.size(); ++k) {
2113 data[slice[k]] = defaultValue;
2114 }
2115}
2116
2117template<class T>
2118void setAuxiliarySingle(size_t loc, AnyValue& extra, const AnyValue& value)
2119{
2120 extra.asVector<T>()[loc] = value.as<T>();
2121}
2122
2123template<class T>
2124void setAuxiliaryMulti(size_t loc, AnyValue& extra, const AnyValue& data)
2125{
2126 const auto& value = data.asVector<T>();
2127 auto& vec = extra.asVector<vector<T>>();
2128 if (value.size() != vec[loc].size()) {
2129 throw CanteraError("SolutionArray::setAuxiliaryMulti",
2130 "New element size {} does not match existing column size {}.",
2131 value.size(), vec[loc].size());
2132 }
2133 vec[loc] = value;
2134}
2135
2136} // end unnamed namespace
2137
2138}
Header for a simple thermodynamics model of a surface phase derived from ThermoPhase,...
Header file for class ThermoPhase, the base class for phases with thermodynamic properties,...
void setLoc(int line, int column)
For values which are derived from an input file, set the line and column of this value in that file.
Definition AnyMap.cpp:617
A map of string keys to values whose type can vary at runtime.
Definition AnyMap.h:431
bool empty() const
Return boolean indicating whether AnyMap is empty.
Definition AnyMap.cpp:1468
void clear()
Erase all items in the mapping.
Definition AnyMap.cpp:1488
const AnyValue & at(const string &key) const
Get the value of the item stored in key.
Definition AnyMap.cpp:1458
void update(const AnyMap &other, bool keepExisting=true)
Add items from other to this AnyMap.
Definition AnyMap.cpp:1493
A wrapper for a variable whose type is determined at runtime.
Definition AnyMap.h:87
bool isVector() const
Returns true if the held value is a vector of the specified type, such as vector<double>.
Definition AnyMap.inl.h:75
pair< size_t, size_t > matrixShape() const
Returns rows and columns of a matrix.
Definition AnyMap.cpp:714
size_t vectorSize() const
Returns size of the held vector.
Definition AnyMap.cpp:698
bool isScalar() const
Returns true if the held value is a scalar type (such as double, long int, string,...
Definition AnyMap.cpp:694
bool is() const
Returns true if the held value is of the specified type.
Definition AnyMap.inl.h:68
const vector< T > & asVector(size_t nMin=npos, size_t nMax=npos) const
Return the held value, if it is a vector of type T.
Definition AnyMap.inl.h:109
const T & as() const
Get the value of this key as the specified type.
Definition AnyMap.inl.h:16
string type_str() const
Returns a string specifying the type of the held value.
Definition AnyMap.cpp:686
Base class for exceptions thrown by Cantera classes.
const char * what() const override
Get a description of the error.
virtual string getMessage() const
Method overridden by derived classes to format the error message.
An array index is out of range.
An error indicating that an unimplemented function has been called.
A wrapper class handling storage to HDF.
Definition Storage.h:39
pair< size_t, set< string > > contents(const string &id) const
Retrieve contents of file from a specified location.
Definition Storage.cpp:609
void writeAttributes(const string &id, const AnyMap &meta)
Write attributes to a specified location.
Definition Storage.cpp:627
void deleteGroup(const string &id)
Delete group.
Definition Storage.cpp:603
AnyMap readAttributes(const string &id, bool recursive) const
Read attributes from a specified location.
Definition Storage.cpp:621
AnyValue readData(const string &id, const string &name, size_t rows, size_t cols=npos) const
Read dataset from a specified location.
Definition Storage.cpp:633
bool checkGroup(const string &id, bool permissive=false)
Check whether path location exists.
Definition Storage.cpp:597
void setCompressionLevel(int level)
Set compression level (0..9)
Definition Storage.cpp:585
void writeData(const string &id, const string &name, const AnyValue &data)
Write dataset to a specified location.
Definition Storage.cpp:640
void fmt_append(fmt::memory_buffer &b, const std::string &tmpl, Args... args)
Versions 6.2.0 and 6.2.1 of fmtlib do not include this define before they include windows....
Definition fmt.h:29
void tokenizePath(const string &in_val, vector< string > &v)
This function separates a string up into tokens according to the location of path separators.
string toLowerCopy(const string &input)
Convert to lower case.
U len(const T &container)
Get the size of a container, cast to a signed integer type.
Definition utilities.h:198
double dot(InputIter x_begin, InputIter x_end, InputIter2 y_begin)
Function that calculates a templated inner product.
Definition utilities.h:82
void warn_user(const string &method, const string &msg, const Args &... args)
Print a user warning raised from method as CanteraWarning.
Definition global.h:263
Namespace for the Cantera kernel.
Definition AnyMap.cpp:595
const size_t npos
index returned by functions to indicate "no position"
Definition ct_defs.h:180
offset
Offsets of solution components in the 1D solution array.
Definition Flow1D.h:24
Contains declarations for string manipulation functions within Cantera.
Various templated functions that carry out common vector and polynomial operations (see Templated Arr...