Cantera  2.2.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
GibbsExcessVPSSTP.cpp
Go to the documentation of this file.
1 /**
2  * @file GibbsExcessVPSSTP.cpp
3  * Definitions for intermediate ThermoPhase object for phases which
4  * employ excess Gibbs free energy formulations
5  * (see \ref thermoprops and class \link Cantera::GibbsExcessVPSSTP GibbsExcessVPSSTP\endlink).
6  *
7  * Header file for a derived class of ThermoPhase that handles
8  * variable pressure standard state methods for calculating
9  * thermodynamic properties that are further based upon expressions
10  * for the excess Gibbs free energy expressed as a function of
11  * the mole fractions.
12  */
13 /*
14  * Copyright (2009) Sandia Corporation. Under the terms of
15  * Contract DE-AC04-94AL85000 with Sandia Corporation, the
16  * U.S. Government retains certain rights in this software.
17  */
20 
21 using namespace std;
22 
23 namespace Cantera
24 {
25 
26 GibbsExcessVPSSTP::GibbsExcessVPSSTP(const GibbsExcessVPSSTP& b)
27 {
28  GibbsExcessVPSSTP::operator=(b);
29 }
30 
31 GibbsExcessVPSSTP& GibbsExcessVPSSTP::operator=(const GibbsExcessVPSSTP& b)
32 {
33  if (&b == this) {
34  return *this;
35  }
36 
37  VPStandardStateTP::operator=(b);
38 
39  moleFractions_ = b.moleFractions_;
40  lnActCoeff_Scaled_ = b.lnActCoeff_Scaled_;
41  dlnActCoeffdT_Scaled_ = b.dlnActCoeffdT_Scaled_;
42  d2lnActCoeffdT2_Scaled_ = b.d2lnActCoeffdT2_Scaled_;
43  dlnActCoeffdlnX_diag_ = b.dlnActCoeffdlnX_diag_;
44  dlnActCoeffdlnN_diag_ = b.dlnActCoeffdlnN_diag_;
45  dlnActCoeffdlnN_ = b.dlnActCoeffdlnN_;
46  m_pp = b.m_pp;
47 
48  return *this;
49 }
50 
51 ThermoPhase* GibbsExcessVPSSTP::duplMyselfAsThermoPhase() const
52 {
53  return new GibbsExcessVPSSTP(*this);
54 }
55 
56 void GibbsExcessVPSSTP::setMassFractions(const doublereal* const y)
57 {
58  Phase::setMassFractions(y);
59  getMoleFractions(DATA_PTR(moleFractions_));
60 }
61 
62 void GibbsExcessVPSSTP::setMassFractions_NoNorm(const doublereal* const y)
63 {
64  Phase::setMassFractions_NoNorm(y);
65  getMoleFractions(DATA_PTR(moleFractions_));
66 }
67 
68 void GibbsExcessVPSSTP::setMoleFractions(const doublereal* const x)
69 {
70  Phase::setMoleFractions(x);
71  getMoleFractions(DATA_PTR(moleFractions_));
72 }
73 
74 void GibbsExcessVPSSTP::setMoleFractions_NoNorm(const doublereal* const x)
75 {
76  Phase::setMoleFractions_NoNorm(x);
77  getMoleFractions(DATA_PTR(moleFractions_));
78 }
79 
80 void GibbsExcessVPSSTP::setConcentrations(const doublereal* const c)
81 {
82  Phase::setConcentrations(c);
83  getMoleFractions(DATA_PTR(moleFractions_));
84 }
85 
86 /*
87  * ------------ Mechanical Properties ------------------------------
88  */
89 void GibbsExcessVPSSTP::setPressure(doublereal p)
90 {
91  setState_TP(temperature(), p);
92 }
93 
94 void GibbsExcessVPSSTP::calcDensity()
95 {
96  vector_fp vbar = getPartialMolarVolumesVector();
97  doublereal vtotal = 0.0;
98  for (size_t i = 0; i < m_kk; i++) {
99  vtotal += vbar[i] * moleFractions_[i];
100  }
101  doublereal dd = meanMolecularWeight() / vtotal;
102  Phase::setDensity(dd);
103 }
104 
105 void GibbsExcessVPSSTP::setState_TP(doublereal t, doublereal p)
106 {
107  Phase::setTemperature(t);
108  /*
109  * Store the current pressure
110  */
111  m_Pcurrent = p;
112  /*
113  * update the standard state thermo
114  * -> This involves calling the water function and setting the pressure
115  */
116  updateStandardStateThermo();
117 
118  /*
119  * Calculate the partial molar volumes, and then the density of the fluid
120  */
121  calcDensity();
122 }
123 
124 /*
125  * - Activities, Standard States, Activity Concentrations -----------
126  */
127 void GibbsExcessVPSSTP::getActivityConcentrations(doublereal* c) const
128 {
129  getActivities(c);
130 }
131 
132 doublereal GibbsExcessVPSSTP::standardConcentration(size_t k) const
133 {
134  return 1.0;
135 }
136 
137 doublereal GibbsExcessVPSSTP::logStandardConc(size_t k) const
138 {
139  return 0.0;
140 }
141 
142 void GibbsExcessVPSSTP::getActivities(doublereal* ac) const
143 {
144  getActivityCoefficients(ac);
145  getMoleFractions(DATA_PTR(moleFractions_));
146  for (size_t k = 0; k < m_kk; k++) {
147  ac[k] *= moleFractions_[k];
148  }
149 }
150 
151 void GibbsExcessVPSSTP::getActivityCoefficients(doublereal* const ac) const
152 {
153  getLnActivityCoefficients(ac);
154  //
155  // Protect against or inform about roundoff when taking exponentials
156  //
157  if ((DEBUG_MODE_ENABLED && realNumberRangeBehavior_ == THROWON_OVERFLOW_DEBUGMODEONLY_CTRB) ||
158  (realNumberRangeBehavior_ == THROWON_OVERFLOW_CTRB)) {
159  for (size_t k = 0; k < m_kk; k++) {
160  if (ac[k] > 700.) {
161  throw CanteraError("GibbsExcessVPSSTP::getActivityCoefficients()",
162  "activity coefficient for " + int2str(k) + " is overflowing: ln(ac) = " + fp2str(ac[k]));
163  } else if (ac[k] < -700.) {
164  throw CanteraError("GibbsExcessVPSSTP::getActivityCoefficients()",
165  "activity coefficient for " + int2str(k) + " is underflowing: ln(ac) = " + fp2str(ac[k]));
166  } else {
167  ac[k] = exp(ac[k]);
168  }
169  }
170  } else if (realNumberRangeBehavior_ == CHANGE_OVERFLOW_CTRB) {
171  for (size_t k = 0; k < m_kk; k++) {
172  if (ac[k] > 700.) {
173  ac[k] = exp(700.0);
174  } else if (ac[k] < -700.) {
175  ac[k] = exp(-700.0);
176  } else {
177  ac[k] = exp(ac[k]);
178  }
179  }
180  } else {
181  for (size_t k = 0; k < m_kk; k++) {
182  ac[k] = exp(ac[k]);
183  }
184  if (realNumberRangeBehavior_ == FENV_CHECK_CTRB) {
185 #ifdef HAVE_FENV_H
187  throw CanteraError("GibbsExcessVPSSTP::getActivityCoefficients()",
188  "activity coefficient is over/underflowing");
189  }
190 #else
191  throw CanteraError("GibbsExcessVPSSTP::getActivityCoefficients()",
192  "realNumberRangeBehavior_ == FENV_CHECK_CTRB not supported by compiler");
193 #endif
194  }
195  }
196 }
197 
198 void GibbsExcessVPSSTP::getElectrochemPotentials(doublereal* mu) const
199 {
200  getChemPotentials(mu);
201  double ve = Faraday * electricPotential();
202  for (size_t k = 0; k < m_kk; k++) {
203  mu[k] += ve*charge(k);
204  }
205 }
206 
207 /*
208  * ------------ Partial Molar Properties of the Solution ------------
209  */
210 void GibbsExcessVPSSTP::getPartialMolarVolumes(doublereal* vbar) const
211 {
212  /*
213  * Get the standard state values in m^3 kmol-1
214  */
215  getStandardVolumes(vbar);
216 }
217 
218 const vector_fp& GibbsExcessVPSSTP::getPartialMolarVolumesVector() const
219 {
220  return getStandardVolumes();
221 }
222 
223 double GibbsExcessVPSSTP::checkMFSum(const doublereal* const x) const
224 {
225  doublereal norm = std::accumulate(x, x + m_kk, 0.0);
226  if (fabs(norm - 1.0) > 1.0E-9) {
227  throw CanteraError("GibbsExcessVPSSTP::checkMFSum",
228  "(MF sum - 1) exceeded tolerance of 1.0E-9:" + fp2str(norm));
229  }
230  return norm;
231 }
232 
233 void GibbsExcessVPSSTP::getUnitsStandardConc(double* uA, int k, int sizeUA) const
234 {
235  warn_deprecated("GibbsExcessVPSSTP::getUnitsStandardConc",
236  "To be removed after Cantera 2.2.");
237  //
238  // We assume here that the units of the standard concentration is unitless. In other words activities are
239  // used unchanged in kinetics expressions. This may be changed in implementations of child classes.
240  //
241  for (int i = 0; i < sizeUA; i++) {
242  if (i == 0) {
243  uA[0] = 0.0;
244  }
245  if (i == 1) {
246  uA[1] = 0.0;
247  }
248  if (i == 2) {
249  uA[2] = 0.0;
250  }
251  if (i == 3) {
252  uA[3] = 0.0;
253  }
254  if (i == 4) {
255  uA[4] = 0.0;
256  }
257  if (i == 5) {
258  uA[5] = 0.0;
259  }
260  }
261 }
262 
263 void GibbsExcessVPSSTP::initThermo()
264 {
265  initLengths();
266  VPStandardStateTP::initThermo();
267  getMoleFractions(DATA_PTR(moleFractions_));
268 }
269 
270 void GibbsExcessVPSSTP::initLengths()
271 {
272  moleFractions_.resize(m_kk);
273  lnActCoeff_Scaled_.resize(m_kk);
274  dlnActCoeffdT_Scaled_.resize(m_kk);
275  d2lnActCoeffdT2_Scaled_.resize(m_kk);
276  dlnActCoeffdlnX_diag_.resize(m_kk);
277  dlnActCoeffdlnN_diag_.resize(m_kk);
278  dlnActCoeffdlnN_.resize(m_kk, m_kk);
279  m_pp.resize(m_kk);
280 }
281 
282 } // end of namespace Cantera
std::vector< doublereal > m_pp
Temporary storage space that is fair game.
std::string int2str(const int n, const std::string &fmt)
Convert an int to a string using a format converter.
Definition: stringUtils.cpp:39
For this specification of range behavior, the overflow or underflow calculation is changed...
Definition: ctexceptions.h:74
When an overflow or underflow occurs, Cantera will throw an error.
Definition: ctexceptions.h:77
void warn_deprecated(const std::string &method, const std::string &extra)
Print a warning indicating that method is deprecated.
Definition: global.cpp:78
std::vector< doublereal > d2lnActCoeffdT2_Scaled_
Storage for the current derivative values of the gradients with respect to temperature of the log of ...
std::vector< doublereal > dlnActCoeffdlnN_diag_
Storage for the current derivative values of the gradients with respect to logarithm of the mole frac...
Base class for a phase with thermodynamic properties.
Definition: ThermoPhase.h:97
std::string fp2str(const double x, const std::string &fmt)
Convert a double into a c++ string.
Definition: stringUtils.cpp:28
bool check_FENV_OverUnder_Flow()
Quick check on whether there has been an underflow or overflow condition in the floating point unit...
Header for intermediate ThermoPhase object for phases which employ Gibbs excess free energy based for...
Base class for exceptions thrown by Cantera classes.
Definition: ctexceptions.h:99
std::vector< doublereal > dlnActCoeffdlnX_diag_
Storage for the current derivative values of the gradients with respect to logarithm of the mole frac...
std::vector< doublereal > lnActCoeff_Scaled_
Storage for the current values of the activity coefficients of the species.
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
std::vector< doublereal > moleFractions_
Storage for the current values of the mole fractions of the species.
Contains declarations for string manipulation functions within Cantera.
#define DATA_PTR(vec)
Creates a pointer to the start of the raw data for a vector.
Definition: ct_defs.h:36
Array2D dlnActCoeffdlnN_
Storage for the current derivative values of the gradients with respect to logarithm of the species m...
Cantera will throw an error in debug mode but will not in production mode.
Definition: ctexceptions.h:86
Cantera will use the fenv check capability introduced in C99 to check for overflow and underflow cond...
Definition: ctexceptions.h:82
std::vector< doublereal > dlnActCoeffdT_Scaled_
Storage for the current derivative values of the gradients with respect to temperature of the log of ...