Cantera  2.0
IdealSolnGasVPSS.cpp
Go to the documentation of this file.
1 /**
2  * @file IdealSolnGasVPSS.cpp
3  * Definition file for a derived class of ThermoPhase that assumes either
4  * an ideal gas or ideal solution approximation and handles
5  * variable pressure standard state methods for calculating
6  * thermodynamic properties (see \ref thermoprops and
7  * class \link Cantera::IdealSolnGasVPSS IdealSolnGasVPSS\endlink).
8  */
9 /*
10  * Copyright (2005) Sandia Corporation. Under the terms of
11  * Contract DE-AC04-94AL85000 with Sandia Corporation, the
12  * U.S. Government retains certain rights in this software.
13  */
14 
16 #include "cantera/thermo/VPSSMgr.h"
17 #include "cantera/thermo/PDSS.h"
18 #include "cantera/thermo/mix_defs.h"
21 
22 using namespace std;
23 
24 namespace Cantera
25 {
26 
27 /*
28  * Default constructor
29  */
30 IdealSolnGasVPSS::IdealSolnGasVPSS() :
32  m_idealGas(0),
33  m_formGC(0)
34 {
35 }
36 
37 
38 IdealSolnGasVPSS::IdealSolnGasVPSS(std::string infile, std::string id) :
40  m_idealGas(0),
41  m_formGC(0)
42 {
43  XML_Node* root = get_XML_File(infile);
44  if (id == "-") {
45  id = "";
46  }
47  XML_Node* xphase = get_XML_NameID("phase", std::string("#")+id, root);
48  if (!xphase) {
49  throw CanteraError("newPhase",
50  "Couldn't find phase named \"" + id + "\" in file, " + infile);
51  }
52  importPhase(*xphase, this);
53 }
54 
55 /*
56  * Copy Constructor:
57  *
58  * Note this stuff will not work until the underlying phase
59  * has a working copy constructor.
60  *
61  * The copy constructor just calls the assignment operator
62  * to do the heavy lifting.
63  */
66  m_idealGas(0),
67  m_formGC(0)
68 {
69  *this = b;
70 }
71 
72 /*
73  * operator=()
74  *
75  * Note this stuff will not work until the underlying phase
76  * has a working assignment operator
77  */
80 {
81  if (&b != this) {
82  /*
83  * Mostly, this is a passthrough to the underlying
84  * assignment operator for the ThermoPhae parent object.
85  */
87  /*
88  * However, we have to handle data that we own.
89  */
91  m_formGC = b.m_formGC;
92  }
93  return *this;
94 }
95 
96 /*
97  * ~IdealSolnGasVPSS(): (virtual)
98  *
99  */
101 {
102 }
103 
104 /*
105  * Duplication function.
106  * This calls the copy constructor for this object.
107  */
109 {
110  IdealSolnGasVPSS* vptp = new IdealSolnGasVPSS(*this);
111  return (ThermoPhase*) vptp;
112 }
113 
115 {
116  if (m_idealGas) {
117  return cIdealSolnGasVPSS;
118  }
119  return cIdealSolnGasVPSS_iscv;
120 }
121 
122 
123 /*
124  * ------------Molar Thermodynamic Properties -------------------------
125  */
126 
127 /// Molar enthalpy. Units: J/kmol.
129 {
131  const vector_fp& enth_RT = m_VPSS_ptr->enthalpy_RT();
132  return (GasConstant * temperature() *
133  mean_X(DATA_PTR(enth_RT)));
134 }
135 
136 /// Molar internal energy. Units: J/kmol.
138 {
139  doublereal p0 = pressure();
140  doublereal md = molarDensity();
141  return (enthalpy_mole() - p0 / md);
142 }
143 
144 /// Molar entropy. Units: J/kmol/K.
146 {
148  const vector_fp& entrop_R = m_VPSS_ptr->entropy_R();
149  return GasConstant * (mean_X(DATA_PTR(entrop_R)) - sum_xlogx());
150 
151 }
152 
153 /// Molar Gibbs function. Units: J/kmol.
155 {
156  return enthalpy_mole() - temperature() * entropy_mole();
157 }
158 
159 /// Molar heat capacity at constant pressure. Units: J/kmol/K.
160 doublereal IdealSolnGasVPSS::cp_mole() const
161 {
163  const vector_fp& cp_R = m_VPSS_ptr->cp_R();
164  return GasConstant * (mean_X(DATA_PTR(cp_R)));
165 }
166 
167 /// Molar heat capacity at constant volume. Units: J/kmol/K.
168 doublereal IdealSolnGasVPSS::cv_mole() const
169 {
170  return cp_mole() - GasConstant;
171 
172 }
173 
175 {
176  m_Pcurrent = p;
178  calcDensity();
179 }
180 
182 {
183  /*
184  * Calculate the molarVolume of the solution (m**3 kmol-1)
185  */
186  if (m_idealGas) {
187  double dens = (m_Pcurrent * meanMolecularWeight()
188  /(GasConstant * temperature()));
189  Phase::setDensity(dens);
190  } else {
191  const doublereal* const dtmp = moleFractdivMMW();
192  const vector_fp& vss = m_VPSS_ptr->standardVolumes();
193  double invDens = dot(vss.begin(), vss.end(), dtmp);
194  /*
195  * Set the density in the parent State object directly,
196  * by calling the Phase::setDensity() function.
197  */
198  double dens = 1.0/invDens;
199  Phase::setDensity(dens);
200  }
201 }
202 
204 {
205  if (m_idealGas) {
206  return 1.0 / m_Pcurrent;
207  } else {
208  throw CanteraError("IdealSolnGasVPSS::isothermalCompressibility() ",
209  "not implemented");
210  }
211  return 0.0;
212 }
213 
215 {
216  if (m_idealGas) {
218  } else {
219  const vector_fp& vss = m_VPSS_ptr->standardVolumes();
220  switch (m_formGC) {
221  case 0:
222  for (size_t k = 0; k < m_kk; k++) {
223  c[k] = moleFraction(k);
224  }
225  break;
226  case 1:
227  for (size_t k = 0; k < m_kk; k++) {
228  c[k] = moleFraction(k) / vss[k];
229  }
230  break;
231  case 2:
232  for (size_t k = 0; k < m_kk; k++) {
233  c[k] = moleFraction(k) / vss[0];
234  }
235  break;
236  }
237  }
238 }
239 
240 /*
241  * Returns the standard concentration \f$ C^0_k \f$, which is used to normalize
242  * the generalized concentration.
243  */
244 doublereal IdealSolnGasVPSS::standardConcentration(size_t k) const
245 {
246  if (m_idealGas) {
247  double p = pressure();
248  return p/(GasConstant * temperature());
249  } else {
250  const vector_fp& vss = m_VPSS_ptr->standardVolumes();
251  switch (m_formGC) {
252  case 0:
253  return 1.0;
254  case 1:
255  return 1.0 / vss[k];
256  case 2:
257  return 1.0/ vss[0];
258  }
259  return 0.0;
260 
261  }
262 }
263 
264 /*
265  * Returns the natural logarithm of the standard
266  * concentration of the kth species
267  */
268 doublereal IdealSolnGasVPSS::logStandardConc(size_t k) const
269 {
270  double c = standardConcentration(k);
271  double lc = std::log(c);
272  return lc;
273 }
274 
275 /*
276  *
277  * getUnitsStandardConcentration()
278  *
279  * Returns the units of the standard and general concentrations
280  * Note they have the same units, as their divisor is
281  * defined to be equal to the activity of the kth species
282  * in the solution, which is unitless.
283  *
284  * This routine is used in print out applications where the
285  * units are needed. Usually, MKS units are assumed throughout
286  * the program and in the XML input files.
287  *
288  * uA[0] = kmol units - default = 1
289  * uA[1] = m units - default = -nDim(), the number of spatial
290  * dimensions in the Phase class.
291  * uA[2] = kg units - default = 0;
292  * uA[3] = Pa(pressure) units - default = 0;
293  * uA[4] = Temperature units - default = 0;
294  * uA[5] = time units - default = 0
295  *
296  * For EOS types other than cIdealSolidSolnPhase1, the default
297  * kmol/m3 holds for standard concentration units. For
298  * cIdealSolidSolnPhase0 type, the standard concentration is
299  * unitless.
300  */
301 void IdealSolnGasVPSS::getUnitsStandardConc(double* uA, int, int sizeUA) const
302 {
303  int eos = eosType();
304  if (eos == cIdealSolnGasPhase0) {
305  for (int i = 0; i < sizeUA; i++) {
306  uA[i] = 0.0;
307  }
308  } else {
309  for (int i = 0; i < sizeUA; i++) {
310  if (i == 0) {
311  uA[0] = 1.0;
312  }
313  if (i == 1) {
314  uA[1] = -int(nDim());
315  }
316  if (i == 2) {
317  uA[2] = 0.0;
318  }
319  if (i == 3) {
320  uA[3] = 0.0;
321  }
322  if (i == 4) {
323  uA[4] = 0.0;
324  }
325  if (i == 5) {
326  uA[5] = 0.0;
327  }
328  }
329  }
330 }
331 
332 
333 /*
334  * Get the array of non-dimensional activity coefficients
335  */
337 {
338  for (size_t k = 0; k < m_kk; k++) {
339  ac[k] = 1.0;
340  }
341 }
342 
343 /*
344  * ---- Partial Molar Properties of the Solution -----------------
345  */
346 
347 /*
348  * Get the array of non-dimensional species chemical potentials
349  * These are partial molar Gibbs free energies.
350  * \f$ \mu_k / \hat R T \f$.
351  * Units: unitless
352  *
353  * We close the loop on this function, here, calling
354  * getChemPotentials() and then dividing by RT.
355  */
356 void IdealSolnGasVPSS::getChemPotentials_RT(doublereal* muRT) const
357 {
358  getChemPotentials(muRT);
359  doublereal invRT = 1.0 / _RT();
360  for (size_t k = 0; k < m_kk; k++) {
361  muRT[k] *= invRT;
362  }
363 }
364 
365 void IdealSolnGasVPSS::getChemPotentials(doublereal* mu) const
366 {
368  doublereal xx;
369  doublereal rt = temperature() * GasConstant;
370  for (size_t k = 0; k < m_kk; k++) {
372  mu[k] += rt*(log(xx));
373  }
374 }
375 
376 
378 {
379  getEnthalpy_RT(hbar);
380  doublereal rt = GasConstant * temperature();
381  scale(hbar, hbar+m_kk, hbar, rt);
382 }
383 
384 void IdealSolnGasVPSS::getPartialMolarEntropies(doublereal* sbar) const
385 {
386  getEntropy_R(sbar);
387  doublereal r = GasConstant;
388  scale(sbar, sbar+m_kk, sbar, r);
389  for (size_t k = 0; k < m_kk; k++) {
390  doublereal xx = std::max(SmallNumber, moleFraction(k));
391  sbar[k] += r * (- log(xx));
392  }
393 }
394 
396 {
397  getIntEnergy_RT(ubar);
398  doublereal rt = GasConstant * temperature();
399  scale(ubar, ubar+m_kk, ubar, rt);
400 }
401 
402 void IdealSolnGasVPSS::getPartialMolarCp(doublereal* cpbar) const
403 {
404  getCp_R(cpbar);
405  doublereal r = GasConstant;
406  scale(cpbar, cpbar+m_kk, cpbar, r);
407 }
408 
409 void IdealSolnGasVPSS::getPartialMolarVolumes(doublereal* vbar) const
410 {
411  getStandardVolumes(vbar);
412 }
413 
414 /*
415  * ----- Thermodynamic Values for the Species Reference States ----
416  */
417 
418 
419 
420 
421 /*
422  * Perform initializations after all species have been
423  * added.
424  */
426 {
427  initLengths();
429 }
430 
431 
432 void IdealSolnGasVPSS::setToEquilState(const doublereal* mu_RT)
433 {
434  double tmp, tmp2;
436  const vector_fp& grt = m_VPSS_ptr->Gibbs_RT_ref();
437 
438  /*
439  * Within the method, we protect against inf results if the
440  * exponent is too high.
441  *
442  * If it is too low, we set
443  * the partial pressure to zero. This capability is needed
444  * by the elemental potential method.
445  */
446  doublereal pres = 0.0;
447  double m_p0 = m_VPSS_ptr->refPressure();
448  for (size_t k = 0; k < m_kk; k++) {
449  tmp = -grt[k] + mu_RT[k];
450  if (tmp < -600.) {
451  m_pp[k] = 0.0;
452  } else if (tmp > 500.0) {
453  tmp2 = tmp / 500.;
454  tmp2 *= tmp2;
455  m_pp[k] = m_p0 * exp(500.) * tmp2;
456  } else {
457  m_pp[k] = m_p0 * exp(tmp);
458  }
459  pres += m_pp[k];
460  }
461  // set state
462  setState_PX(pres, &m_pp[0]);
463 }
464 
465 /*
466  * Initialize the internal lengths.
467  * (this is not a virtual function)
468  */
470 {
471  m_kk = nSpecies();
472  m_pp.resize(m_kk, 0.0);
473 }
474 
475 /*
476  * Import and initialize a ThermoPhase object
477  *
478  * param phaseNode This object must be the phase node of a
479  * complete XML tree
480  * description of the phase, including all of the
481  * species data. In other words while "phase" must
482  * point to an XML phase object, it must have
483  * sibling nodes "speciesData" that describe
484  * the species in the phase.
485  * param id ID of the phase. If nonnull, a check is done
486  * to see if phaseNode is pointing to the phase
487  * with the correct id.
488  *
489  * This routine initializes the lengths in the current object and
490  * then calls the parent routine.
491  */
492 void IdealSolnGasVPSS::initThermoXML(XML_Node& phaseNode, std::string id)
493 {
495 
496  if (phaseNode.hasChild("thermo")) {
497  XML_Node& thermoNode = phaseNode.child("thermo");
498  std::string model = thermoNode["model"];
499  if (model == "IdealGasVPSS") {
500  m_idealGas = 1;
501  } else if (model == "IdealSolnVPSS") {
502  m_idealGas = 0;
503  } else {
504  throw CanteraError("IdealSolnGasVPSS::initThermoXML",
505  "Unknown thermo model : " + model);
506  }
507  }
508 
509  /*
510  * Form of the standard concentrations. Must have one of:
511  *
512  * <standardConc model="unity" />
513  * <standardConc model="molar_volume" />
514  * <standardConc model="solvent_volume" />
515  */
516  if (phaseNode.hasChild("standardConc")) {
517  if (m_idealGas) {
518  throw CanteraError("IdealSolnGasVPSS::initThermoXML",
519  "standardConc node for ideal gas");
520  }
521  XML_Node& scNode = phaseNode.child("standardConc");
522  string formStringa = scNode.attrib("model");
523  string formString = lowercase(formStringa);
524  if (formString == "unity") {
525  m_formGC = 0;
526  } else if (formString == "molar_volume") {
527  m_formGC = 1;
528  } else if (formString == "solvent_volume") {
529  m_formGC = 2;
530  } else {
531  throw CanteraError("initThermoXML",
532  "Unknown standardConc model: " + formStringa);
533  }
534  } else {
535  if (!m_idealGas) {
536  throw CanteraError("initThermoXML",
537  "Unspecified standardConc model");
538  }
539  }
540 
541  VPStandardStateTP::initThermoXML(phaseNode, id);
542 }
543 
545 {
547  std::string model = thermoNode["model"];
548  if (model == "IdealGasVPSS") {
549  m_idealGas = 1;
550  } else if (model == "IdealSolnVPSS") {
551  m_idealGas = 0;
552  } else {
553  throw CanteraError("IdealSolnGasVPSS::initThermoXML",
554  "Unknown thermo model : " + model);
555  }
556 }
557 
558 }
559 
560