Cantera  2.4.0
ThermoPhase.cpp
Go to the documentation of this file.
1 /**
2  * @file ThermoPhase.cpp
3  * Definition file for class ThermoPhase, the base class for phases with
4  * thermodynamic properties
5  * (see class \link Cantera::ThermoPhase ThermoPhase\endlink).
6  */
7 
8 // This file is part of Cantera. See License.txt in the top-level directory or
9 // at http://www.cantera.org/license.txt for license and copyright information.
10 
17 #include "cantera/base/ctml.h"
18 
19 #include <iomanip>
20 #include <fstream>
21 
22 using namespace std;
23 
24 namespace Cantera
25 {
26 
27 ThermoPhase::ThermoPhase() :
28  m_speciesData(0),
29  m_phi(0.0),
30  m_hasElementPotentials(false),
31  m_chargeNeutralityNecessary(false),
32  m_ssConvention(cSS_CONVENTION_TEMPERATURE),
33  m_tlast(0.0)
34 {
35 }
36 
37 ThermoPhase::~ThermoPhase()
38 {
39  for (size_t k = 0; k < m_speciesData.size(); k++) {
40  delete m_speciesData[k];
41  }
42 }
43 
44 void ThermoPhase::resetHf298(size_t k) {
45  if (k != npos) {
47  } else {
48  for (size_t k = 0; k < nSpecies(); k++) {
50  }
51  }
53 }
54 
56 {
57  return cAC_CONVENTION_MOLAR;
58 }
59 
61 {
62  return m_ssConvention;
63 }
64 
65 doublereal ThermoPhase::logStandardConc(size_t k) const
66 {
67  return log(standardConcentration(k));
68 }
69 
70 void ThermoPhase::getActivities(doublereal* a) const
71 {
73  for (size_t k = 0; k < nSpecies(); k++) {
74  a[k] /= standardConcentration(k);
75  }
76 }
77 
78 void ThermoPhase::getLnActivityCoefficients(doublereal* lnac) const
79 {
81  for (size_t k = 0; k < m_kk; k++) {
82  lnac[k] = std::log(lnac[k]);
83  }
84 }
85 
86 void ThermoPhase::getElectrochemPotentials(doublereal* mu) const
87 {
89  double ve = Faraday * electricPotential();
90  for (size_t k = 0; k < m_kk; k++) {
91  mu[k] += ve*charge(k);
92  }
93 }
94 
95 void ThermoPhase::setState_TPX(doublereal t, doublereal p, const doublereal* x)
96 {
98  setState_TP(t,p);
99 }
100 
101 void ThermoPhase::setState_TPX(doublereal t, doublereal p, const compositionMap& x)
102 {
104  setState_TP(t,p);
105 }
106 
107 void ThermoPhase::setState_TPX(doublereal t, doublereal p, const std::string& x)
108 {
110  setState_TP(t,p);
111 }
112 
113 void ThermoPhase::setState_TPY(doublereal t, doublereal p, const doublereal* y)
114 {
115  setMassFractions(y);
116  setState_TP(t,p);
117 }
118 
119 void ThermoPhase::setState_TPY(doublereal t, doublereal p, const compositionMap& y)
120 {
122  setState_TP(t,p);
123 }
124 
125 void ThermoPhase::setState_TPY(doublereal t, doublereal p, const std::string& y)
126 {
128  setState_TP(t,p);
129 }
130 
131 void ThermoPhase::setState_TP(doublereal t, doublereal p)
132 {
133  setTemperature(t);
134  setPressure(p);
135 }
136 
137 void ThermoPhase::setState_RPX(doublereal rho, doublereal p, const doublereal* x)
138 {
139  setMoleFractions(x);
140  setState_RP(rho, p);
141 }
142 
143 void ThermoPhase::setState_RPX(doublereal rho, doublereal p, const compositionMap& x)
144 {
146  setState_RP(rho,p);
147 }
148 
149 void ThermoPhase::setState_RPX(doublereal rho, doublereal p, const std::string& x)
150 {
152  setState_RP(rho,p);
153 }
154 
155 void ThermoPhase::setState_RPY(doublereal rho, doublereal p, const doublereal* y)
156 {
157  setMassFractions(y);
158  setState_RP(rho,p);
159 }
160 
161 void ThermoPhase::setState_RPY(doublereal rho, doublereal p, const compositionMap& y)
162 {
164  setState_RP(rho,p);
165 }
166 
167 void ThermoPhase::setState_RPY(doublereal rho, doublereal p, const std::string& y)
168 {
170  setState_RP(rho,p);
171 }
172 
173 void ThermoPhase::setState_PX(doublereal p, doublereal* x)
174 {
175  setMoleFractions(x);
176  setPressure(p);
177 }
178 
179 void ThermoPhase::setState_PY(doublereal p, doublereal* y)
180 {
181  setMassFractions(y);
182  setPressure(p);
183 }
184 
185 void ThermoPhase::setState_HP(double Htarget, double p, double rtol)
186 {
187  setState_HPorUV(Htarget, p, rtol, false);
188 }
189 
190 void ThermoPhase::setState_UV(double u, double v, double rtol)
191 {
192  setState_HPorUV(u, v, rtol, true);
193 }
194 
195 void ThermoPhase::setState_conditional_TP(doublereal t, doublereal p, bool set_p)
196 {
197  setTemperature(t);
198  if (set_p) {
199  setPressure(p);
200  }
201 }
202 
203 void ThermoPhase::setState_HPorUV(double Htarget, double p,
204  double rtol, bool doUV)
205 {
206  doublereal dt;
207  doublereal v = 0.0;
208 
209  // Assign the specific volume or pressure and make sure it's positive
210  if (doUV) {
211  doublereal v = p;
212  if (v < 1.0E-300) {
213  throw CanteraError("setState_HPorUV (UV)",
214  "Input specific volume is too small or negative. v = {}", v);
215  }
216  setDensity(1.0/v);
217  } else {
218  if (p < 1.0E-300) {
219  throw CanteraError("setState_HPorUV (HP)",
220  "Input pressure is too small or negative. p = {}", p);
221  }
222  setPressure(p);
223  }
224  double Tmax = maxTemp() + 0.1;
225  double Tmin = minTemp() - 0.1;
226 
227  // Make sure we are within the temperature bounds at the start
228  // of the iteration
229  double Tnew = temperature();
230  double Tinit = Tnew;
231  if (Tnew > Tmax) {
232  Tnew = Tmax - 1.0;
233  } else if (Tnew < Tmin) {
234  Tnew = Tmin + 1.0;
235  }
236  if (Tnew != Tinit) {
237  setState_conditional_TP(Tnew, p, !doUV);
238  }
239 
240  double Hnew = (doUV) ? intEnergy_mass() : enthalpy_mass();
241  double Cpnew = (doUV) ? cv_mass() : cp_mass();
242  double Htop = Hnew;
243  double Ttop = Tnew;
244  double Hbot = Hnew;
245  double Tbot = Tnew;
246 
247  bool ignoreBounds = false;
248  // Unstable phases are those for which cp < 0.0. These are possible for
249  // cases where we have passed the spinodal curve.
250  bool unstablePhase = false;
251  // Counter indicating the last temperature point where the
252  // phase was unstable
253  double Tunstable = -1.0;
254  bool unstablePhaseNew = false;
255 
256  // Newton iteration
257  for (int n = 0; n < 500; n++) {
258  double Told = Tnew;
259  double Hold = Hnew;
260  double cpd = Cpnew;
261  if (cpd < 0.0) {
262  unstablePhase = true;
263  Tunstable = Tnew;
264  }
265  // limit step size to 100 K
266  dt = clip((Htarget - Hold)/cpd, -100.0, 100.0);
267 
268  // Calculate the new T
269  Tnew = Told + dt;
270 
271  // Limit the step size so that we are convergent This is the step that
272  // makes it different from a Newton's algorithm
273  if ((dt > 0.0 && unstablePhase) || (dt <= 0.0 && !unstablePhase)) {
274  if (Hbot < Htarget && Tnew < (0.75 * Tbot + 0.25 * Told)) {
275  dt = 0.75 * (Tbot - Told);
276  Tnew = Told + dt;
277  }
278  } else if (Htop > Htarget && Tnew > (0.75 * Ttop + 0.25 * Told)) {
279  dt = 0.75 * (Ttop - Told);
280  Tnew = Told + dt;
281  }
282 
283  // Check Max and Min values
284  if (Tnew > Tmax && !ignoreBounds) {
285  setState_conditional_TP(Tmax, p, !doUV);
286  double Hmax = (doUV) ? intEnergy_mass() : enthalpy_mass();
287  if (Hmax >= Htarget) {
288  if (Htop < Htarget) {
289  Ttop = Tmax;
290  Htop = Hmax;
291  }
292  } else {
293  Tnew = Tmax + 1.0;
294  ignoreBounds = true;
295  }
296  }
297  if (Tnew < Tmin && !ignoreBounds) {
298  setState_conditional_TP(Tmin, p, !doUV);
299  double Hmin = (doUV) ? intEnergy_mass() : enthalpy_mass();
300  if (Hmin <= Htarget) {
301  if (Hbot > Htarget) {
302  Tbot = Tmin;
303  Hbot = Hmin;
304  }
305  } else {
306  Tnew = Tmin - 1.0;
307  ignoreBounds = true;
308  }
309  }
310 
311  // Try to keep phase within its region of stability
312  // -> Could do a lot better if I calculate the
313  // spinodal value of H.
314  for (int its = 0; its < 10; its++) {
315  Tnew = Told + dt;
316  if (Tnew < Told / 3.0) {
317  Tnew = Told / 3.0;
318  dt = -2.0 * Told / 3.0;
319  }
320  setState_conditional_TP(Tnew, p, !doUV);
321  if (doUV) {
322  Hnew = intEnergy_mass();
323  Cpnew = cv_mass();
324  } else {
325  Hnew = enthalpy_mass();
326  Cpnew = cp_mass();
327  }
328  if (Cpnew < 0.0) {
329  unstablePhaseNew = true;
330  Tunstable = Tnew;
331  } else {
332  unstablePhaseNew = false;
333  break;
334  }
335  if (unstablePhase == false && unstablePhaseNew == true) {
336  dt *= 0.25;
337  }
338  }
339 
340  if (Hnew == Htarget) {
341  return;
342  } else if (Hnew > Htarget && (Htop < Htarget || Hnew < Htop)) {
343  Htop = Hnew;
344  Ttop = Tnew;
345  } else if (Hnew < Htarget && (Hbot > Htarget || Hnew > Hbot)) {
346  Hbot = Hnew;
347  Tbot = Tnew;
348  }
349  // Convergence in H
350  double Herr = Htarget - Hnew;
351  double acpd = std::max(fabs(cpd), 1.0E-5);
352  double denom = std::max(fabs(Htarget), acpd * Tnew);
353  double HConvErr = fabs((Herr)/denom);
354  if (HConvErr < rtol || fabs(dt/Tnew) < rtol) {
355  return;
356  }
357  }
358  // We are here when there hasn't been convergence
359 
360  // Formulate a detailed error message, since questions seem to arise often
361  // about the lack of convergence.
362  string ErrString = "No convergence in 500 iterations\n";
363  if (doUV) {
364  ErrString += fmt::format(
365  "\tTarget Internal Energy = {}\n"
366  "\tCurrent Specific Volume = {}\n"
367  "\tStarting Temperature = {}\n"
368  "\tCurrent Temperature = {}\n"
369  "\tCurrent Internal Energy = {}\n"
370  "\tCurrent Delta T = {}\n",
371  Htarget, v, Tinit, Tnew, Hnew, dt);
372  } else {
373  ErrString += fmt::format(
374  "\tTarget Enthalpy = {}\n"
375  "\tCurrent Pressure = {}\n"
376  "\tStarting Temperature = {}\n"
377  "\tCurrent Temperature = {}\n"
378  "\tCurrent Enthalpy = {}\n"
379  "\tCurrent Delta T = {}\n",
380  Htarget, p, Tinit, Tnew, Hnew, dt);
381  }
382  if (unstablePhase) {
383  ErrString += fmt::format(
384  "\t - The phase became unstable (Cp < 0) T_unstable_last = {}\n",
385  Tunstable);
386  }
387  if (doUV) {
388  throw CanteraError("setState_HPorUV (UV)", ErrString);
389  } else {
390  throw CanteraError("setState_HPorUV (HP)", ErrString);
391  }
392 }
393 
394 void ThermoPhase::setState_SP(double Starget, double p, double rtol)
395 {
396  setState_SPorSV(Starget, p, rtol, false);
397 }
398 
399 void ThermoPhase::setState_SV(double Starget, double v, double rtol)
400 {
401  setState_SPorSV(Starget, v, rtol, true);
402 }
403 
404 void ThermoPhase::setState_SPorSV(double Starget, double p,
405  double rtol, bool doSV)
406 {
407  doublereal v = 0.0;
408  doublereal dt;
409  if (doSV) {
410  v = p;
411  if (v < 1.0E-300) {
412  throw CanteraError("setState_SPorSV (SV)",
413  "Input specific volume is too small or negative. v = {}", v);
414  }
415  setDensity(1.0/v);
416  } else {
417  if (p < 1.0E-300) {
418  throw CanteraError("setState_SPorSV (SP)",
419  "Input pressure is too small or negative. p = {}", p);
420  }
421  setPressure(p);
422  }
423  double Tmax = maxTemp() + 0.1;
424  double Tmin = minTemp() - 0.1;
425 
426  // Make sure we are within the temperature bounds at the start
427  // of the iteration
428  double Tnew = temperature();
429  double Tinit = Tnew;
430  if (Tnew > Tmax) {
431  Tnew = Tmax - 1.0;
432  } else if (Tnew < Tmin) {
433  Tnew = Tmin + 1.0;
434  }
435  if (Tnew != Tinit) {
436  setState_conditional_TP(Tnew, p, !doSV);
437  }
438 
439  double Snew = entropy_mass();
440  double Cpnew = (doSV) ? cv_mass() : cp_mass();
441  double Stop = Snew;
442  double Ttop = Tnew;
443  double Sbot = Snew;
444  double Tbot = Tnew;
445 
446  bool ignoreBounds = false;
447  // Unstable phases are those for which Cp < 0.0. These are possible for
448  // cases where we have passed the spinodal curve.
449  bool unstablePhase = false;
450  double Tunstable = -1.0;
451  bool unstablePhaseNew = false;
452 
453  // Newton iteration
454  for (int n = 0; n < 500; n++) {
455  double Told = Tnew;
456  double Sold = Snew;
457  double cpd = Cpnew;
458  if (cpd < 0.0) {
459  unstablePhase = true;
460  Tunstable = Tnew;
461  }
462  // limit step size to 100 K
463  dt = clip((Starget - Sold)*Told/cpd, -100.0, 100.0);
464  Tnew = Told + dt;
465 
466  // Limit the step size so that we are convergent
467  if ((dt > 0.0 && unstablePhase) || (dt <= 0.0 && !unstablePhase)) {
468  if (Sbot < Starget && Tnew < Tbot) {
469  dt = 0.75 * (Tbot - Told);
470  Tnew = Told + dt;
471  }
472  } else if (Stop > Starget && Tnew > Ttop) {
473  dt = 0.75 * (Ttop - Told);
474  Tnew = Told + dt;
475  }
476 
477  // Check Max and Min values
478  if (Tnew > Tmax && !ignoreBounds) {
479  setState_conditional_TP(Tmax, p, !doSV);
480  double Smax = entropy_mass();
481  if (Smax >= Starget) {
482  if (Stop < Starget) {
483  Ttop = Tmax;
484  Stop = Smax;
485  }
486  } else {
487  Tnew = Tmax + 1.0;
488  ignoreBounds = true;
489  }
490  } else if (Tnew < Tmin && !ignoreBounds) {
491  setState_conditional_TP(Tmin, p, !doSV);
492  double Smin = entropy_mass();
493  if (Smin <= Starget) {
494  if (Sbot > Starget) {
495  Tbot = Tmin;
496  Sbot = Smin;
497  }
498  } else {
499  Tnew = Tmin - 1.0;
500  ignoreBounds = true;
501  }
502  }
503 
504  // Try to keep phase within its region of stability
505  // -> Could do a lot better if I calculate the
506  // spinodal value of H.
507  for (int its = 0; its < 10; its++) {
508  Tnew = Told + dt;
509  setState_conditional_TP(Tnew, p, !doSV);
510  Cpnew = (doSV) ? cv_mass() : cp_mass();
511  Snew = entropy_mass();
512  if (Cpnew < 0.0) {
513  unstablePhaseNew = true;
514  Tunstable = Tnew;
515  } else {
516  unstablePhaseNew = false;
517  break;
518  }
519  if (unstablePhase == false && unstablePhaseNew == true) {
520  dt *= 0.25;
521  }
522  }
523 
524  if (Snew == Starget) {
525  return;
526  } else if (Snew > Starget && (Stop < Starget || Snew < Stop)) {
527  Stop = Snew;
528  Ttop = Tnew;
529  } else if (Snew < Starget && (Sbot > Starget || Snew > Sbot)) {
530  Sbot = Snew;
531  Tbot = Tnew;
532  }
533  // Convergence in S
534  double Serr = Starget - Snew;
535  double acpd = std::max(fabs(cpd), 1.0E-5);
536  double denom = std::max(fabs(Starget), acpd * Tnew);
537  double SConvErr = fabs((Serr * Tnew)/denom);
538  if (SConvErr < rtol || fabs(dt/Tnew) < rtol) {
539  return;
540  }
541  }
542  // We are here when there hasn't been convergence
543 
544  // Formulate a detailed error message, since questions seem to arise often
545  // about the lack of convergence.
546  string ErrString = "No convergence in 500 iterations\n";
547  if (doSV) {
548  ErrString += fmt::format(
549  "\tTarget Entropy = {}\n"
550  "\tCurrent Specific Volume = {}\n"
551  "\tStarting Temperature = {}\n"
552  "\tCurrent Temperature = {}\n"
553  "\tCurrent Entropy = {}\n"
554  "\tCurrent Delta T = {}\n",
555  Starget, v, Tinit, Tnew, Snew, dt);
556  } else {
557  ErrString += fmt::format(
558  "\tTarget Entropy = {}\n"
559  "\tCurrent Pressure = {}\n"
560  "\tStarting Temperature = {}\n"
561  "\tCurrent Temperature = {}\n"
562  "\tCurrent Entropy = {}\n"
563  "\tCurrent Delta T = {}\n",
564  Starget, p, Tinit, Tnew, Snew, dt);
565  }
566  if (unstablePhase) {
567  ErrString += fmt::format("\t - The phase became unstable (Cp < 0) T_unstable_last = {}\n",
568  Tunstable);
569  }
570  if (doSV) {
571  throw CanteraError("setState_SPorSV (SV)", ErrString);
572  } else {
573  throw CanteraError("setState_SPorSV (SP)", ErrString);
574  }
575 }
576 
578 {
579  return m_spthermo;
580 }
581 
582 void ThermoPhase::initThermoFile(const std::string& inputFile,
583  const std::string& id)
584 {
585  XML_Node* fxml = get_XML_File(inputFile);
586  XML_Node* fxml_phase = findXMLPhase(fxml, id);
587  if (!fxml_phase) {
588  throw CanteraError("ThermoPhase::initThermoFile",
589  "ERROR: Can not find phase named {} in file"
590  " named {}", id, inputFile);
591  }
592  importPhase(*fxml_phase, this);
593 }
594 
595 void ThermoPhase::initThermoXML(XML_Node& phaseNode, const std::string& id)
596 {
597  if (phaseNode.hasChild("state")) {
598  setStateFromXML(phaseNode.child("state"));
599  }
600 }
601 
603 {
604  // Check to see that all of the species thermo objects have been initialized
605  if (!m_spthermo.ready(m_kk)) {
606  throw CanteraError("ThermoPhase::initThermo()",
607  "Missing species thermo data");
608  }
609 }
610 
611 bool ThermoPhase::addSpecies(shared_ptr<Species> spec)
612 {
613  if (!spec->thermo) {
614  throw CanteraError("ThermoPhase::addSpecies",
615  "Species {} has no thermo data", spec->name);
616  }
617  bool added = Phase::addSpecies(spec);
618  if (added) {
619  spec->thermo->validate(spec->name);
620  m_spthermo.install_STIT(m_kk-1, spec->thermo);
621  }
622  return added;
623 }
624 
625 void ThermoPhase::modifySpecies(size_t k, shared_ptr<Species> spec)
626 {
627  if (!spec->thermo) {
628  throw CanteraError("ThermoPhase::modifySpecies",
629  "Species {} has no thermo data", spec->name);
630  }
631  Phase::modifySpecies(k, spec);
632  if (speciesName(k) != spec->name) {
633  throw CanteraError("ThermoPhase::modifySpecies",
634  "New species '{}' does not match existing species '{}' at index {}",
635  spec->name, speciesName(k), k);
636  }
637  spec->thermo->validate(spec->name);
638  m_spthermo.modifySpecies(k, spec->thermo);
639 }
640 
641 void ThermoPhase::saveSpeciesData(const size_t k, const XML_Node* const data)
642 {
643  if (m_speciesData.size() < (k + 1)) {
644  m_speciesData.resize(k+1, 0);
645  }
646  m_speciesData[k] = new XML_Node(*data);
647 }
648 
649 const std::vector<const XML_Node*> & ThermoPhase::speciesData() const
650 {
651  if (m_speciesData.size() != m_kk) {
652  throw CanteraError("ThermoPhase::speciesData",
653  "m_speciesData is the wrong size");
654  }
655  return m_speciesData;
656 }
657 
659 {
660  string comp = getChildValue(state,"moleFractions");
661  if (comp != "") {
663  } else {
664  comp = getChildValue(state,"massFractions");
665  if (comp != "") {
667  }
668  }
669  if (state.hasChild("temperature")) {
670  double t = getFloat(state, "temperature", "temperature");
671  setTemperature(t);
672  }
673  if (state.hasChild("pressure")) {
674  double p = getFloat(state, "pressure", "pressure");
675  setPressure(p);
676  }
677  if (state.hasChild("density")) {
678  double rho = getFloat(state, "density", "density");
679  setDensity(rho);
680  }
681 }
682 
685  m_tlast += 0.1234;
686 }
687 
688 void ThermoPhase::equilibrate(const std::string& XY, const std::string& solver,
689  double rtol, int max_steps, int max_iter,
690  int estimate_equil, int log_level)
691 {
692  if (solver == "auto" || solver == "element_potential") {
693  vector_fp initial_state;
694  saveState(initial_state);
695  debuglog("Trying ChemEquil solver\n", log_level);
696  try {
697  ChemEquil E;
698  E.options.maxIterations = max_steps;
699  E.options.relTolerance = rtol;
700  int ret = E.equilibrate(*this, XY.c_str(), log_level-1);
701  if (ret < 0) {
702  throw CanteraError("ThermoPhase::equilibrate",
703  "ChemEquil solver failed. Return code: {}", ret);
704  }
705  m_lambdaRRT.resize(nElements());
706  for (size_t m = 0; m < nElements(); m++) {
707  m_lambdaRRT[m] = E.elementPotentials()[m] / RT();
708  }
709  m_hasElementPotentials = true;
710  debuglog("ChemEquil solver succeeded\n", log_level);
711  return;
712  } catch (std::exception& err) {
713  debuglog("ChemEquil solver failed.\n", log_level);
714  debuglog(err.what(), log_level);
715  restoreState(initial_state);
716  if (solver == "auto") {
717  } else {
718  throw;
719  }
720  }
721  }
722 
723  if (solver == "auto" || solver == "vcs" || solver == "gibbs") {
724  MultiPhase M;
725  M.addPhase(this, 1.0);
726  M.init();
727  M.equilibrate(XY, solver, rtol, max_steps, max_iter,
728  estimate_equil, log_level);
729  return;
730  }
731 
732  if (solver != "auto") {
733  throw CanteraError("ThermoPhase::equilibrate",
734  "Invalid solver specified: '{}'", solver);
735  }
736 }
737 
739 {
740  warn_deprecated("ThermoPhase::setElementPotentials",
741  "To be removed after Cantera 2.4");
742  size_t mm = nElements();
743  if (lambda.size() < mm) {
744  throw CanteraError("setElementPotentials", "lambda too small");
745  }
746  if (!m_hasElementPotentials) {
747  m_lambdaRRT.resize(mm);
748  }
749  scale(lambda.begin(), lambda.end(), m_lambdaRRT.begin(), 1.0/RT());
750  m_hasElementPotentials = true;
751 }
752 
753 bool ThermoPhase::getElementPotentials(doublereal* lambda) const
754 {
755  warn_deprecated("ThermoPhase::getElementPotentials",
756  "To be removed after Cantera 2.4");
758  scale(m_lambdaRRT.begin(), m_lambdaRRT.end(), lambda, RT());
759  }
760  return m_hasElementPotentials;
761 }
762 
763 void ThermoPhase::getdlnActCoeffdlnN(const size_t ld, doublereal* const dlnActCoeffdlnN)
764 {
765  for (size_t m = 0; m < m_kk; m++) {
766  for (size_t k = 0; k < m_kk; k++) {
767  dlnActCoeffdlnN[ld * k + m] = 0.0;
768  }
769  }
770  return;
771 }
772 
773 void ThermoPhase::getdlnActCoeffdlnN_numderiv(const size_t ld, doublereal* const dlnActCoeffdlnN)
774 {
775  double deltaMoles_j = 0.0;
776  double pres = pressure();
777 
778  // Evaluate the current base activity coefficients if necessary
779  vector_fp ActCoeff_Base(m_kk);
780  getActivityCoefficients(ActCoeff_Base.data());
781  vector_fp Xmol_Base(m_kk);
782  getMoleFractions(Xmol_Base.data());
783 
784  // Make copies of ActCoeff and Xmol_ for use in taking differences
785  vector_fp ActCoeff(m_kk);
786  vector_fp Xmol(m_kk);
787  double v_totalMoles = 1.0;
788  double TMoles_base = v_totalMoles;
789 
790  // Loop over the columns species to be deltad
791  for (size_t j = 0; j < m_kk; j++) {
792  // Calculate a value for the delta moles of species j
793  // -> Note Xmol_[] and Tmoles are always positive or zero quantities.
794  // -> experience has shown that you always need to make the deltas
795  // greater than needed to change the other mole fractions in order
796  // to capture some effects.
797  double moles_j_base = v_totalMoles * Xmol_Base[j];
798  deltaMoles_j = 1.0E-7 * moles_j_base + v_totalMoles * 1.0E-13 + 1.0E-150;
799 
800  // Now, update the total moles in the phase and all of the mole
801  // fractions based on this.
802  v_totalMoles = TMoles_base + deltaMoles_j;
803  for (size_t k = 0; k < m_kk; k++) {
804  Xmol[k] = Xmol_Base[k] * TMoles_base / v_totalMoles;
805  }
806  Xmol[j] = (moles_j_base + deltaMoles_j) / v_totalMoles;
807 
808  // Go get new values for the activity coefficients.
809  // -> Note this calls setState_PX();
810  setState_PX(pres, Xmol.data());
811  getActivityCoefficients(ActCoeff.data());
812 
813  // Calculate the column of the matrix
814  double* const lnActCoeffCol = dlnActCoeffdlnN + ld * j;
815  for (size_t k = 0; k < m_kk; k++) {
816  lnActCoeffCol[k] = (2*moles_j_base + deltaMoles_j) *(ActCoeff[k] - ActCoeff_Base[k]) /
817  ((ActCoeff[k] + ActCoeff_Base[k]) * deltaMoles_j);
818  }
819  // Revert to the base case Xmol_, v_totalMoles
820  v_totalMoles = TMoles_base;
821  Xmol = Xmol_Base;
822  }
823 
824  setState_PX(pres, Xmol_Base.data());
825 }
826 
827 std::string ThermoPhase::report(bool show_thermo, doublereal threshold) const
828 {
829  fmt::memory_buffer b;
830  try {
831  if (name() != "") {
832  format_to(b, "\n {}:\n", name());
833  }
834  format_to(b, "\n");
835  format_to(b, " temperature {:12.6g} K\n", temperature());
836  format_to(b, " pressure {:12.6g} Pa\n", pressure());
837  format_to(b, " density {:12.6g} kg/m^3\n", density());
838  format_to(b, " mean mol. weight {:12.6g} amu\n", meanMolecularWeight());
839 
840  doublereal phi = electricPotential();
841  if (phi != 0.0) {
842  format_to(b, " potential {:12.6g} V\n", phi);
843  }
844  if (show_thermo) {
845  format_to(b, "\n");
846  format_to(b, " 1 kg 1 kmol\n");
847  format_to(b, " ----------- ------------\n");
848  format_to(b, " enthalpy {:12.5g} {:12.4g} J\n",
850  format_to(b, " internal energy {:12.5g} {:12.4g} J\n",
852  format_to(b, " entropy {:12.5g} {:12.4g} J/K\n",
854  format_to(b, " Gibbs function {:12.5g} {:12.4g} J\n",
855  gibbs_mass(), gibbs_mole());
856  format_to(b, " heat capacity c_p {:12.5g} {:12.4g} J/K\n",
857  cp_mass(), cp_mole());
858  try {
859  format_to(b, " heat capacity c_v {:12.5g} {:12.4g} J/K\n",
860  cv_mass(), cv_mole());
861  } catch (NotImplementedError&) {
862  format_to(b, " heat capacity c_v <not implemented> \n");
863  }
864  }
865 
866  vector_fp x(m_kk);
867  vector_fp y(m_kk);
868  vector_fp mu(m_kk);
869  getMoleFractions(&x[0]);
870  getMassFractions(&y[0]);
871  getChemPotentials(&mu[0]);
872  int nMinor = 0;
873  doublereal xMinor = 0.0;
874  doublereal yMinor = 0.0;
875  format_to(b, "\n");
876  if (show_thermo) {
877  format_to(b, " X "
878  " Y Chem. Pot. / RT\n");
879  format_to(b, " ------------- "
880  "------------ ------------\n");
881  for (size_t k = 0; k < m_kk; k++) {
882  if (abs(x[k]) >= threshold) {
883  if (abs(x[k]) > SmallNumber) {
884  format_to(b, "{:>18s} {:12.6g} {:12.6g} {:12.6g}\n",
885  speciesName(k), x[k], y[k], mu[k]/RT());
886  } else {
887  format_to(b, "{:>18s} {:12.6g} {:12.6g}\n",
888  speciesName(k), x[k], y[k]);
889  }
890  } else {
891  nMinor++;
892  xMinor += x[k];
893  yMinor += y[k];
894  }
895  }
896  } else {
897  format_to(b, " X Y\n");
898  format_to(b, " ------------- ------------\n");
899  for (size_t k = 0; k < m_kk; k++) {
900  if (abs(x[k]) >= threshold) {
901  format_to(b, "{:>18s} {:12.6g} {:12.6g}\n",
902  speciesName(k), x[k], y[k]);
903  } else {
904  nMinor++;
905  xMinor += x[k];
906  yMinor += y[k];
907  }
908  }
909  }
910  if (nMinor) {
911  format_to(b, " [{:+5d} minor] {:12.6g} {:12.6g}\n",
912  nMinor, xMinor, yMinor);
913  }
914  } catch (CanteraError& err) {
915  return to_string(b) + err.what();
916  }
917  return to_string(b);
918 }
919 
920 void ThermoPhase::reportCSV(std::ofstream& csvFile) const
921 {
922  int tabS = 15;
923  int tabM = 30;
924  csvFile.precision(8);
925  vector_fp X(nSpecies());
926  getMoleFractions(&X[0]);
927  std::vector<std::string> pNames;
928  std::vector<vector_fp> data;
929  getCsvReportData(pNames, data);
930 
931  csvFile << setw(tabS) << "Species,";
932  for (size_t i = 0; i < pNames.size(); i++) {
933  csvFile << setw(tabM) << pNames[i] << ",";
934  }
935  csvFile << endl;
936  for (size_t k = 0; k < nSpecies(); k++) {
937  csvFile << setw(tabS) << speciesName(k) + ",";
938  if (X[k] > SmallNumber) {
939  for (size_t i = 0; i < pNames.size(); i++) {
940  csvFile << setw(tabM) << data[i][k] << ",";
941  }
942  csvFile << endl;
943  } else {
944  for (size_t i = 0; i < pNames.size(); i++) {
945  csvFile << setw(tabM) << 0 << ",";
946  }
947  csvFile << endl;
948  }
949  }
950 }
951 
952 void ThermoPhase::getCsvReportData(std::vector<std::string>& names,
953  std::vector<vector_fp>& data) const
954 {
955  names.clear();
956  data.assign(10, vector_fp(nSpecies()));
957 
958  names.push_back("X");
959  getMoleFractions(&data[0][0]);
960 
961  names.push_back("Y");
962  getMassFractions(&data[1][0]);
963 
964  names.push_back("Chem. Pot (J/kmol)");
965  getChemPotentials(&data[2][0]);
966 
967  names.push_back("Activity");
968  getActivities(&data[3][0]);
969 
970  names.push_back("Act. Coeff.");
971  getActivityCoefficients(&data[4][0]);
972 
973  names.push_back("Part. Mol Enthalpy (J/kmol)");
974  getPartialMolarEnthalpies(&data[5][0]);
975 
976  names.push_back("Part. Mol. Entropy (J/K/kmol)");
977  getPartialMolarEntropies(&data[6][0]);
978 
979  names.push_back("Part. Mol. Energy (J/kmol)");
980  getPartialMolarIntEnergies(&data[7][0]);
981 
982  names.push_back("Part. Mol. Cp (J/K/kmol");
983  getPartialMolarCp(&data[8][0]);
984 
985  names.push_back("Part. Mol. Cv (J/K/kmol)");
986  getPartialMolarVolumes(&data[9][0]);
987 }
988 
989 }
std::map< std::string, doublereal > compositionMap
Map connecting a string name with a double.
Definition: ct_defs.h:149
Class ChemEquil implements a chemical equilibrium solver for single-phase solutions.
Definition: ChemEquil.h:82
virtual void getLnActivityCoefficients(doublereal *lnac) const
Get the array of non-dimensional molar-based ln activity coefficients at the current solution tempera...
Definition: ThermoPhase.cpp:78
EquilOpt options
Options controlling how the calculation is carried out.
Definition: ChemEquil.h:136
virtual void setState_SV(double s, double v, double tol=1e-9)
Set the specific entropy (J/kg/K) and specific volume (m^3/kg).
size_t nElements() const
Number of elements.
Definition: Phase.cpp:88
virtual bool addSpecies(shared_ptr< Species > spec)
const vector_fp & elementPotentials() const
Definition: ChemEquil.h:128
void restoreState(const vector_fp &state)
Restore a state saved on a previous call to saveState.
Definition: Phase.cpp:233
virtual void setState_RPX(doublereal rho, doublereal p, const doublereal *x)
Set the density (kg/m**3), pressure (Pa) and mole fractions.
XML_Node * get_XML_File(const std::string &file, int debug)
Return a pointer to the XML tree for a Cantera input file.
Definition: global.cpp:105
CTML ("Cantera Markup Language") is the variant of XML that Cantera uses to store data...
MultiSpeciesThermo m_spthermo
Pointer to the calculation manager for species reference-state thermodynamic properties.
Definition: ThermoPhase.h:1607
virtual void modifySpecies(size_t k, shared_ptr< Species > spec)
Modify the thermodynamic data associated with a species.
XML_Node * findXMLPhase(XML_Node *root, const std::string &idtarget)
Search an XML_Node tree for a named phase XML_Node.
Definition: xml.cpp:1038
void getMassFractions(doublereal *const y) const
Get the species mass fractions.
Definition: Phase.cpp:508
std::string getChildValue(const XML_Node &parent, const std::string &nameString)
This function reads a child node with the name, nameString, and returns its XML value as the return s...
Definition: ctml.cpp:131
doublereal relTolerance
Relative tolerance.
Definition: ChemEquil.h:34
bool ready(size_t nSpecies)
Check if data for all species (0 through nSpecies-1) has been installed.
doublereal temperature() const
Temperature (K).
Definition: Phase.h:601
const std::vector< const XML_Node * > & speciesData() const
Return a pointer to the vector of XML nodes containing the species data for this phase.
std::vector< const XML_Node * > m_speciesData
Vector of pointers to the species databases.
Definition: ThermoPhase.h:1615
virtual doublereal cv_mole() const
Molar heat capacity at constant volume. Units: J/kmol/K.
Definition: ThermoPhase.h:230
An error indicating that an unimplemented function has been called.
Definition: ctexceptions.h:187
virtual void resetHf298(const size_t k)
Restore the original heat of formation of one or more species.
virtual void setState_TPX(doublereal t, doublereal p, const doublereal *x)
Set the temperature (K), pressure (Pa), and mole fractions.
Definition: ThermoPhase.cpp:95
virtual void modifySpecies(size_t index, shared_ptr< SpeciesThermoInterpType > spec)
Modify the species thermodynamic property parameterization for a species.
virtual void getPartialMolarEntropies(doublereal *sbar) const
Returns an array of partial molar entropies of the species in the solution.
Definition: ThermoPhase.h:497
doublereal cp_mass() const
Specific heat at constant pressure. Units: J/kg/K.
Definition: ThermoPhase.h:734
void saveState(vector_fp &state) const
Save the current internal state of the phase.
Definition: Phase.cpp:221
virtual void initThermo()
Initialize the ThermoPhase object after all species have been set up.
const size_t npos
index returned by functions to indicate "no position"
Definition: ct_defs.h:165
virtual void getChemPotentials(doublereal *mu) const
Get the species chemical potentials. Units: J/kmol.
Definition: ThermoPhase.h:461
const int cAC_CONVENTION_MOLAR
Standard state uses the molar convention.
Definition: ThermoPhase.h:26
bool getElementPotentials(doublereal *lambda) const
Returns the element potentials stored in the ThermoPhase object.
Headers for the factory class that can create known ThermoPhase objects (see Thermodynamic Properties...
void addPhase(ThermoPhase *p, doublereal moles)
Add a phase to the mixture.
Definition: MultiPhase.cpp:49
const char * what() const
Get a description of the error.
const int cSS_CONVENTION_TEMPERATURE
Standard state uses the molar convention.
Definition: ThermoPhase.h:36
Class XML_Node is a tree-based representation of the contents of an XML file.
Definition: xml.h:97
void setState_conditional_TP(doublereal t, doublereal p, bool set_p)
Helper function used by setState_HPorUV and setState_SPorSV.
virtual void getActivityConcentrations(doublereal *c) const
This method returns an array of generalized concentrations.
Definition: ThermoPhase.h:370
T clip(const T &value, const T &lower, const T &upper)
Clip value such that lower <= value <= upper.
Definition: global.h:269
void warn_deprecated(const std::string &method, const std::string &extra)
Print a warning indicating that method is deprecated.
Definition: global.cpp:54
size_t nSpecies() const
Returns the number of species in the phase.
Definition: Phase.h:266
STL namespace.
virtual doublereal density() const
Density (kg/m^3).
Definition: Phase.h:607
vector_fp m_lambdaRRT
Vector of element potentials.
Definition: ThermoPhase.h:1622
virtual doublereal minTemp(size_t k=npos) const
Minimum temperature for which the thermodynamic data for the species or phase are valid...
Definition: ThermoPhase.h:131
doublereal enthalpy_mass() const
Specific enthalpy. Units: J/kg.
Definition: ThermoPhase.h:714
void setState_HPorUV(doublereal h, doublereal p, doublereal tol=1e-9, bool doUV=false)
Carry out work in HP and UV calculations.
virtual void setState_HP(double h, double p, double tol=1e-9)
Set the internally stored specific enthalpy (J/kg) and pressure (Pa) of the phase.
virtual doublereal enthalpy_mole() const
Molar enthalpy. Units: J/kmol.
Definition: ThermoPhase.h:205
virtual void getPartialMolarVolumes(doublereal *vbar) const
Return an array of partial molar volumes for the species in the mixture.
Definition: ThermoPhase.h:528
void equilibrate(const std::string &XY, const std::string &solver="auto", double rtol=1e-9, int max_steps=50000, int max_iter=100, int estimate_equil=0, int log_level=0)
Equilibrate a ThermoPhase object.
void saveSpeciesData(const size_t k, const XML_Node *const data)
Store a reference pointer to the XML tree containing the species data for this phase.
virtual MultiSpeciesThermo & speciesThermo(int k=-1)
Return a changeable reference to the calculation manager for species reference-state thermodynamic pr...
doublereal m_tlast
last value of the temperature processed by reference state
Definition: ThermoPhase.h:1643
virtual doublereal standardConcentration(size_t k=0) const
Return the standard concentration for the kth species.
Definition: ThermoPhase.h:391
Headers for the MultiPhase object that is used to set up multiphase equilibrium problems (see Equilfu...
virtual void setState_PY(doublereal p, doublereal *y)
Set the internally stored pressure (Pa) and mass fractions.
void setState_SPorSV(double s, double p, double tol=1e-9, bool doSV=false)
Carry out work in SP and SV calculations.
doublereal RT() const
Return the Gas Constant multiplied by the current temperature.
Definition: ThermoPhase.h:748
virtual void getActivities(doublereal *a) const
Get the array of non-dimensional activities at the current solution temperature, pressure, and solution concentration.
Definition: ThermoPhase.cpp:70
Pure Virtual Base class for individual species reference state thermodynamic managers and text for th...
virtual void invalidateCache()
Invalidate any cached values which are normally updated only when a change in state is detected...
virtual void setStateFromXML(const XML_Node &state)
Set the initial state of the phase to the conditions specified in the state XML element.
virtual doublereal intEnergy_mole() const
Molar internal energy. Units: J/kmol.
Definition: ThermoPhase.h:210
A class for multiphase mixtures.
Definition: MultiPhase.h:57
void getElectrochemPotentials(doublereal *mu) const
Get the species electrochemical potentials.
Definition: ThermoPhase.cpp:86
virtual doublereal gibbs_mole() const
Molar Gibbs function. Units: J/kmol.
Definition: ThermoPhase.h:220
virtual void getCsvReportData(std::vector< std::string > &names, std::vector< vector_fp > &data) const
Fills names and data with the column names and species thermo properties to be included in the output...
std::string speciesName(size_t k) const
Name of the species with index k.
Definition: Phase.cpp:191
virtual void setState_TPY(doublereal t, doublereal p, const doublereal *y)
Set the internally stored temperature (K), pressure (Pa), and mass fractions of the phase...
doublereal electricPotential() const
Returns the electric potential of this phase (V).
Definition: ThermoPhase.h:302
void setMoleFractionsByName(const compositionMap &xMap)
Set the species mole fractions by name.
Definition: Phase.cpp:292
int m_ssConvention
Contains the standard state convention.
Definition: ThermoPhase.h:1640
doublereal entropy_mass() const
Specific entropy. Units: J/kg/K.
Definition: ThermoPhase.h:724
virtual int activityConvention() const
This method returns the convention used in specification of the activities, of which there are curren...
Definition: ThermoPhase.cpp:55
void setMassFractionsByName(const compositionMap &yMap)
Set the species mass fractions by name.
Definition: Phase.cpp:336
virtual void modifySpecies(size_t k, shared_ptr< Species > spec)
Modify the thermodynamic data associated with a species.
Definition: Phase.cpp:786
Base class for exceptions thrown by Cantera classes.
Definition: ctexceptions.h:65
virtual void getPartialMolarEnthalpies(doublereal *hbar) const
Returns an array of partial molar enthalpies for the species in the mixture.
Definition: ThermoPhase.h:487
virtual int standardStateConvention() const
This method returns the convention used in specification of the standard state, of which there are cu...
Definition: ThermoPhase.cpp:60
virtual void invalidateCache()
Invalidate any cached values which are normally updated only when a change in state is detected...
Definition: Phase.cpp:830
virtual void getPartialMolarIntEnergies(doublereal *ubar) const
Return an array of partial molar internal energies for the species in the mixture.
Definition: ThermoPhase.h:507
virtual doublereal cp_mole() const
Molar heat capacity at constant pressure. Units: J/kmol/K.
Definition: ThermoPhase.h:225
virtual void setMoleFractions(const doublereal *const x)
Set the mole fractions to the specified values.
Definition: Phase.cpp:251
virtual doublereal logStandardConc(size_t k=0) const
Natural logarithm of the standard concentration of the kth species.
Definition: ThermoPhase.cpp:65
virtual void setState_UV(double u, double v, double tol=1e-9)
Set the specific internal energy (J/kg) and specific volume (m^3/kg).
virtual void resetHf298(const size_t k=npos)
Restore the original heat of formation of one or more species.
Definition: ThermoPhase.cpp:44
virtual void setState_RP(doublereal rho, doublereal p)
Set the density (kg/m**3) and pressure (Pa) at constant composition.
Definition: ThermoPhase.h:1059
virtual bool addSpecies(shared_ptr< Species > spec)
Add a Species to this Phase.
Definition: Phase.cpp:697
void getMoleFractions(doublereal *const x) const
Get the species mole fraction vector.
Definition: Phase.cpp:466
void importPhase(XML_Node &phase, ThermoPhase *th)
Import a phase information into an empty ThermoPhase object.
virtual void setState_TP(doublereal t, doublereal p)
Set the temperature (K) and pressure (Pa)
virtual void getActivityCoefficients(doublereal *ac) const
Get the array of non-dimensional molar-based activity coefficients at the current solution temperatur...
Definition: ThermoPhase.h:420
virtual doublereal pressure() const
Return the thermodynamic pressure (Pa).
Definition: ThermoPhase.h:245
int equilibrate(thermo_t &s, const char *XY, int loglevel=0)
Definition: ChemEquil.cpp:308
virtual void setState_SP(double s, double p, double tol=1e-9)
Set the specific entropy (J/kg/K) and pressure (Pa).
bool hasChild(const std::string &ch) const
Tests whether the current node has a child node with a particular name.
Definition: xml.cpp:536
void debuglog(const std::string &msg, int loglevel)
Write a message to the log only if loglevel > 0.
Definition: global.h:135
XML_Node & child(const size_t n) const
Return a changeable reference to the n&#39;th child of the current node.
Definition: xml.cpp:546
doublereal cv_mass() const
Specific heat at constant volume. Units: J/kg/K.
Definition: ThermoPhase.h:739
virtual std::string report(bool show_thermo=true, doublereal threshold=-1e-14) const
returns a summary of the state of the phase as a string
virtual void setPressure(doublereal p)
Set the internally stored pressure (Pa) at constant temperature and composition.
Definition: ThermoPhase.h:770
virtual void initThermoXML(XML_Node &phaseNode, const std::string &id)
Import and initialize a ThermoPhase object using an XML tree.
const doublereal SmallNumber
smallest number to compare to zero.
Definition: ct_defs.h:126
std::vector< double > vector_fp
Turn on the use of stl vectors for the basic array type within cantera Vector of doubles.
Definition: ct_defs.h:157
doublereal meanMolecularWeight() const
The mean molecular weight. Units: (kg/kmol)
Definition: Phase.h:661
virtual void setTemperature(const doublereal temp)
Set the internally stored temperature of the phase (K).
Definition: Phase.h:637
void scale(InputIter begin, InputIter end, OutputIter out, S scale_factor)
Multiply elements of an array by a scale factor.
Definition: utilities.h:130
doublereal gibbs_mass() const
Specific Gibbs function. Units: J/kg.
Definition: ThermoPhase.h:729
void init()
Process phases and build atomic composition array.
Definition: MultiPhase.cpp:116
std::string name() const
Return the name of the phase.
Definition: Phase.cpp:78
virtual void setState_PX(doublereal p, doublereal *x)
Set the pressure (Pa) and mole fractions.
Contains declarations for string manipulation functions within Cantera.
doublereal getFloat(const XML_Node &parent, const std::string &name, const std::string &type)
Get a floating-point value from a child element.
Definition: ctml.cpp:164
virtual void setMassFractions(const doublereal *const y)
Set the mass fractions to the specified values and normalize them.
Definition: Phase.cpp:311
void setElementPotentials(const vector_fp &lambda)
Stores the element potentials in the ThermoPhase object.
int maxIterations
Maximum number of iterations.
Definition: ChemEquil.h:36
size_t m_kk
Number of species in the phase.
Definition: Phase.h:788
virtual void getdlnActCoeffdlnN(const size_t ld, doublereal *const dlnActCoeffdlnN)
Get the array of derivatives of the log activity coefficients with respect to the log of the species ...
A species thermodynamic property manager for a phase.
virtual void initThermoFile(const std::string &inputFile, const std::string &id)
virtual void reportCSV(std::ofstream &csvFile) const
returns a summary of the state of the phase to a comma separated file.
virtual void getPartialMolarCp(doublereal *cpbar) const
Return an array of partial molar heat capacities for the species in the mixture.
Definition: ThermoPhase.h:518
virtual void install_STIT(size_t index, shared_ptr< SpeciesThermoInterpType > stit)
Install a new species thermodynamic property parameterization for one species.
Namespace for the Cantera kernel.
Definition: AnyMap.cpp:8
doublereal intEnergy_mass() const
Specific internal energy. Units: J/kg.
Definition: ThermoPhase.h:719
bool m_hasElementPotentials
Boolean indicating whether there is a valid set of saved element potentials for this phase...
Definition: ThermoPhase.h:1627
virtual void setState_RPY(doublereal rho, doublereal p, const doublereal *y)
Set the density (kg/m**3), pressure (Pa) and mass fractions.
virtual doublereal maxTemp(size_t k=npos) const
Maximum temperature for which the thermodynamic data for the species are valid.
Definition: ThermoPhase.h:184
Header file for class ThermoPhase, the base class for phases with thermodynamic properties, and the text for the Module thermoprops (see Thermodynamic Properties and class ThermoPhase).
virtual doublereal entropy_mole() const
Molar entropy. Units: J/kmol/K.
Definition: ThermoPhase.h:215
Chemical equilibrium.
doublereal charge(size_t k) const
Dimensionless electrical charge of a single molecule of species k The charge is normalized by the the...
Definition: Phase.h:577
virtual void setDensity(const doublereal density_)
Set the internally stored density (kg/m^3) of the phase.
Definition: Phase.h:622