Loading [MathJax]/extensions/tex2jax.js
Cantera  3.2.0a1
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
MargulesVPSSTP.h
Go to the documentation of this file.
1/**
2 * @file MargulesVPSSTP.h (see @ref thermoprops and class @link
3 * Cantera::MargulesVPSSTP MargulesVPSSTP@endlink).
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
9#ifndef CT_MARGULESVPSSTP_H
10#define CT_MARGULESVPSSTP_H
11
12#include "GibbsExcessVPSSTP.h"
13
14namespace Cantera
15{
16
17//! MargulesVPSSTP is a derived class of GibbsExcessVPSSTP that employs the
18//! Margules approximation for the excess Gibbs free energy
19/*!
20 * MargulesVPSSTP derives from class GibbsExcessVPSSTP which is derived from
21 * VPStandardStateTP, and overloads the virtual methods defined there with ones
22 * that use expressions appropriate for the Margules Excess Gibbs free energy
23 * approximation.
24 *
25 * The independent unknowns are pressure, temperature, and mass fraction.
26 *
27 * ## Specification of Species Standard State Properties
28 *
29 * All species are defined to have standard states that depend upon both the
30 * temperature and the pressure. The Margules approximation assumes symmetric
31 * standard states, where all of the standard state assume that the species are
32 * in pure component states at the temperature and pressure of the solution. I
33 * don't think it prevents, however, some species from being dilute in the
34 * solution.
35 *
36 * ## Specification of Solution Thermodynamic Properties
37 *
38 * The molar excess Gibbs free energy is given by the following formula which is
39 * a sum over interactions *i*. Each of the interactions are binary interactions
40 * involving two of the species in the phase, denoted, *Ai* and *Bi*. This is
41 * the generalization of the Margules formulation for a phase that has more than
42 * 2 species.
43 *
44 * @f[
45 * G^E = \sum_i \left( H_{Ei} - T S_{Ei} \right)
46 * @f]
47 * @f[
48 * H^E_i = n X_{Ai} X_{Bi} \left( h_{o,i} + h_{1,i} X_{Bi} \right)
49 * @f]
50 * @f[
51 * S^E_i = n X_{Ai} X_{Bi} \left( s_{o,i} + s_{1,i} X_{Bi} \right)
52 * @f]
53 *
54 * where n is the total moles in the solution.
55 *
56 * The activity of a species defined in the phase is given by an excess Gibbs
57 * free energy formulation.
58 *
59 * @f[
60 * a_k = \gamma_k X_k
61 * @f]
62 *
63 * where
64 *
65 * @f[
66 * R T \ln( \gamma_k )= \frac{d(n G^E)}{d(n_k)}\Bigg|_{n_i}
67 * @f]
68 *
69 * Taking the derivatives results in the following expression
70 *
71 * @f[
72 * R T \ln( \gamma_k )= \sum_i \left( \left( \delta_{Ai,k} X_{Bi} + \delta_{Bi,k} X_{Ai} - X_{Ai} X_{Bi} \right)
73 * \left( g^E_{o,i} + g^E_{1,i} X_{Bi} \right) +
74 * \left( \delta_{Bi,k} - X_{Bi} \right) X_{Ai} X_{Bi} g^E_{1,i} \right)
75 * @f]
76 * where
77 * @f$ g^E_{o,i} = h_{o,i} - T s_{o,i} @f$ and
78 * @f$ g^E_{1,i} = h_{1,i} - T s_{1,i} @f$ and where
79 * @f$ X_k @f$ is the mole fraction of species *k*.
80 *
81 * This object inherits from the class VPStandardStateTP. Therefore, the
82 * specification and calculation of all standard state and reference state
83 * values are handled at that level. Various functional forms for the standard
84 * state are permissible. The chemical potential for species *k* is equal to
85 *
86 * @f[
87 * \mu_k(T,P) = \mu^o_k(T, P) + R T \ln(\gamma_k X_k)
88 * @f]
89 *
90 * The partial molar entropy for species *k* is given by
91 *
92 * @f[
93 * \tilde{s}_k(T,P) = s^o_k(T,P) - R \ln( \gamma_k X_k )
94 * - R T \frac{d \ln(\gamma_k) }{dT}
95 * @f]
96 *
97 * The partial molar enthalpy for species *k* is given by
98 *
99 * @f[
100 * \tilde{h}_k(T,P) = h^o_k(T,P) - R T^2 \frac{d \ln(\gamma_k)}{dT}
101 * @f]
102 *
103 * The partial molar volume for species *k* is
104 *
105 * @f[
106 * \tilde V_k(T,P) = V^o_k(T,P) + R T \frac{d \ln(\gamma_k) }{dP}
107 * @f]
108 *
109 * The partial molar Heat Capacity for species *k* is
110 *
111 * @f[
112 * \tilde{C}_{p,k}(T,P) = C^o_{p,k}(T,P) - 2 R T \frac{d \ln( \gamma_k )}{dT}
113 * - R T^2 \frac{d^2 \ln(\gamma_k) }{{dT}^2}
114 * @f]
115 *
116 * ## Application within Kinetics Managers
117 *
118 * @f$ C^a_k @f$ are defined such that @f$ a_k = C^a_k / C^s_k, @f$ where
119 * @f$ C^s_k @f$ is a standard concentration defined below and @f$ a_k @f$ are
120 * activities used in the thermodynamic functions. These activity (or
121 * generalized) concentrations are used by kinetics manager classes to compute
122 * the forward and reverse rates of elementary reactions. The activity
123 * concentration,@f$ C^a_k @f$,is given by the following expression.
124 *
125 * @f[
126 * C^a_k = C^s_k X_k = \frac{P}{R T} X_k
127 * @f]
128 *
129 * The standard concentration for species *k* is independent of *k* and equal to
130 *
131 * @f[
132 * C^s_k = C^s = \frac{P}{R T}
133 * @f]
134 *
135 * For example, a bulk-phase binary gas reaction between species j and k,
136 * producing a new gas species l would have the following equation for its rate
137 * of progress variable, @f$ R^1 @f$, which has units of kmol m-3 s-1.
138 *
139 * @f[
140 * R^1 = k^1 C_j^a C_k^a = k^1 (C^s a_j) (C^s a_k)
141 * @f]
142 * where
143 * @f[
144 * C_j^a = C^s a_j \mbox{\quad and \quad} C_k^a = C^s a_k
145 * @f]
146 *
147 * @f$ C_j^a @f$ is the activity concentration of species j, and @f$ C_k^a @f$
148 * is the activity concentration of species k. @f$ C^s @f$ is the standard
149 * concentration. @f$ a_j @f$ is the activity of species j which is equal to the
150 * mole fraction of j.
151 *
152 * The reverse rate constant can then be obtained from the law of microscopic
153 * reversibility and the equilibrium expression for the system.
154 *
155 * @f[
156 * \frac{a_j a_k}{ a_l} = K_a^{o,1} = \exp(\frac{\mu^o_l - \mu^o_j - \mu^o_k}{R T} )
157 * @f]
158 *
159 * @f$ K_a^{o,1} @f$ is the dimensionless form of the equilibrium constant,
160 * associated with the pressure dependent standard states @f$ \mu^o_l(T,P) @f$
161 * and their associated activities, @f$ a_l @f$, repeated here:
162 *
163 * @f[
164 * \mu_l(T,P) = \mu^o_l(T, P) + R T \ln a_l
165 * @f]
166 *
167 * We can switch over to expressing the equilibrium constant in terms of the
168 * reference state chemical potentials
169 *
170 * @f[
171 * K_a^{o,1} = \exp(\frac{\mu^{ref}_l - \mu^{ref}_j - \mu^{ref}_k}{R T} ) * \frac{P_{ref}}{P}
172 * @f]
173 *
174 * The concentration equilibrium constant, @f$ K_c @f$, may be obtained by
175 * changing over to activity concentrations. When this is done:
176 *
177 * @f[
178 * \frac{C^a_j C^a_k}{ C^a_l} = C^o K_a^{o,1} = K_c^1 =
179 * \exp(\frac{\mu^{ref}_l - \mu^{ref}_j - \mu^{ref}_k}{R T} ) * \frac{P_{ref}}{RT}
180 * @f]
181 *
182 * %Kinetics managers will calculate the concentration equilibrium constant, @f$
183 * K_c @f$, using the second and third part of the above expression as a
184 * definition for the concentration equilibrium constant.
185 *
186 * For completeness, the pressure equilibrium constant may be obtained as well
187 *
188 * @f[
189 * \frac{P_j P_k}{ P_l P_{ref}} = K_p^1 = \exp(\frac{\mu^{ref}_l - \mu^{ref}_j - \mu^{ref}_k}{R T} )
190 * @f]
191 *
192 * @f$ K_p @f$ is the simplest form of the equilibrium constant for ideal gases.
193 * However, it isn't necessarily the simplest form of the equilibrium constant
194 * for other types of phases; @f$ K_c @f$ is used instead because it is
195 * completely general.
196 *
197 * The reverse rate of progress may be written down as
198 * @f[
199 * R^{-1} = k^{-1} C_l^a = k^{-1} (C^o a_l)
200 * @f]
201 *
202 * where we can use the concept of microscopic reversibility to write the
203 * reverse rate constant in terms of the forward rate constant and the
204 * concentration equilibrium constant, @f$ K_c @f$.
205 *
206 * @f[
207 * k^{-1} = k^1 K^1_c
208 * @f]
209 *
210 * @f$ k^{-1} @f$ has units of s-1.
211 *
212 * @ingroup thermoprops
213 */
215{
216public:
217 //! Construct a MargulesVPSSTP object from an input file
218 /*!
219 * @param inputFile Name of the input file containing the phase definition.
220 * If blank, an empty phase will be created.
221 * @param id name (ID) of the phase in the input file. If empty, the
222 * first phase definition in the input file will be used.
223 */
224 explicit MargulesVPSSTP(const string& inputFile="", const string& id="");
225
226 string type() const override {
227 return "Margules";
228 }
229
230 //! @name Molar Thermodynamic Properties
231 //! @{
232
233 double cv_mole() const override;
234
235 //! @}
236 //! @name Activities, Standard States, and Activity Concentrations
237 //!
238 //! The activity @f$ a_k @f$ of a species in solution is related to the
239 //! chemical potential by @f[ \mu_k = \mu_k^0(T) + \hat R T \ln a_k. @f] The
240 //! quantity @f$ \mu_k^0(T,P) @f$ is the chemical potential at unit activity,
241 //! which depends only on temperature and pressure.
242 //! @{
243
244 void getLnActivityCoefficients(double* lnac) const override;
245
246 //! @}
247 //! @name Partial Molar Properties of the Solution
248 //! @{
249
250 void getChemPotentials(double* mu) const override;
251
252 //! Returns an array of partial molar enthalpies for the species in the
253 //! mixture.
254 /*!
255 * Units (J/kmol)
256 *
257 * For this phase, the partial molar enthalpies are equal to the standard
258 * state enthalpies modified by the derivative of the molality-based
259 * activity coefficient wrt temperature
260 *
261 * @f[
262 * \bar h_k(T,P) = h^o_k(T,P) - R T^2 \frac{d \ln(\gamma_k)}{dT}
263 * @f]
264 *
265 * @param hbar Vector of returned partial molar enthalpies
266 * (length m_kk, units = J/kmol)
267 */
268 void getPartialMolarEnthalpies(double* hbar) const override;
269
270 //! Returns an array of partial molar entropies for the species in the
271 //! mixture.
272 /*!
273 * Units (J/kmol)
274 *
275 * For this phase, the partial molar enthalpies are equal to the standard
276 * state enthalpies modified by the derivative of the activity coefficient
277 * wrt temperature
278 *
279 * @f[
280 * \bar s_k(T,P) = s^o_k(T,P) - R T^2 \frac{d \ln(\gamma_k)}{dT}
281 * - R \ln( \gamma_k X_k)
282 * - R T \frac{d \ln(\gamma_k) }{dT}
283 * @f]
284 *
285 * @param sbar Vector of returned partial molar entropies
286 * (length m_kk, units = J/kmol/K)
287 */
288 void getPartialMolarEntropies(double* sbar) const override;
289
290 //! Returns an array of partial molar entropies for the species in the
291 //! mixture.
292 /*!
293 * Units (J/kmol)
294 *
295 * For this phase, the partial molar enthalpies are equal to the standard
296 * state enthalpies modified by the derivative of the activity coefficient
297 * wrt temperature
298 *
299 * @f[
300 * ???????????????
301 * \bar s_k(T,P) = s^o_k(T,P) - R T^2 \frac{d \ln(\gamma_k)}{dT}
302 * - R \ln( \gamma_k X_k)
303 * - R T \frac{d \ln(\gamma_k) }{dT}
304 * ???????????????
305 * @f]
306 *
307 * @param cpbar Vector of returned partial molar heat capacities
308 * (length m_kk, units = J/kmol/K)
309 */
310 void getPartialMolarCp(double* cpbar) const override;
311
312 void getPartialMolarVolumes(double* vbar) const override;
313
314 //! Get the array of temperature second derivatives of the log activity
315 //! coefficients
316 /*!
317 * units = 1/Kelvin
318 *
319 * @param d2lnActCoeffdT2 Output vector of temperature 2nd derivatives of
320 * the log Activity Coefficients. length = m_kk
321 */
322 void getd2lnActCoeffdT2(double* d2lnActCoeffdT2) const;
323
324 void getdlnActCoeffdT(double* dlnActCoeffdT) const override;
325
326 //! @}
327 //! @name Initialization
328 //!
329 //! The following methods are used in the process of constructing the phase
330 //! and setting its parameters from a specification in an input file. They
331 //! are not normally used in application programs. To see how they are used,
332 //! see importPhase()
333 //! @{
334
335 void initThermo() override;
336 void getParameters(AnyMap& phaseNode) const override;
337
338 //! Add a binary species interaction with the specified parameters
339 /*!
340 * @param speciesA name of the first species
341 * @param speciesB name of the second species
342 * @param h0 first excess enthalpy coefficient [J/kmol]
343 * @param h1 second excess enthalpy coefficient [J/kmol]
344 * @param s0 first excess entropy coefficient [J/kmol/K]
345 * @param s1 second excess entropy coefficient [J/kmol/K]
346 * @param vh0 first enthalpy coefficient for excess volume [m^3/kmol]
347 * @param vh1 second enthalpy coefficient for excess volume [m^3/kmol]
348 * @param vs0 first entropy coefficient for excess volume [m^3/kmol/K]
349 * @param vs1 second entropy coefficient for excess volume [m^3/kmol/K]
350 */
351 void addBinaryInteraction(const string& speciesA,
352 const string& speciesB, double h0, double h1, double s0, double s1,
353 double vh0, double vh1, double vs0, double vs1);
354
355 //! @}
356 //! @name Derivatives of Thermodynamic Variables needed for Applications
357 //! @{
358
359 void getdlnActCoeffds(const double dTds, const double* const dXds,
360 double* dlnActCoeffds) const override;
361 void getdlnActCoeffdlnX_diag(double* dlnActCoeffdlnX_diag) const override;
362 void getdlnActCoeffdlnN_diag(double* dlnActCoeffdlnN_diag) const override;
363 void getdlnActCoeffdlnN(const size_t ld, double* const dlnActCoeffdlnN) override;
364
365 //! @}
366
367private:
368 //! Initialize lengths of local variables after all species have been
369 //! identified.
370 void initLengths();
371
372 //! Update the activity coefficients
373 /*!
374 * This function will be called to update the internally stored natural
375 * logarithm of the activity coefficients
376 */
377 void s_update_lnActCoeff() const;
378
379 //! Update the derivative of the log of the activity coefficients wrt T
380 /*!
381 * This function will be called to update the internally stored derivative
382 * of the natural logarithm of the activity coefficients wrt temperature.
383 */
384 void s_update_dlnActCoeff_dT() const;
385
386 //! Update the derivative of the log of the activity coefficients wrt
387 //! log(mole fraction)
388 /*!
389 * This function will be called to update the internally stored derivative
390 * of the natural logarithm of the activity coefficients wrt logarithm of
391 * the mole fractions.
392 */
394
395 //! Update the derivative of the log of the activity coefficients wrt
396 //! log(moles) - diagonal only
397 /*!
398 * This function will be called to update the internally stored diagonal
399 * entries for the derivative of the natural logarithm of the activity
400 * coefficients wrt logarithm of the moles.
401 */
403
404 //! Update the derivative of the log of the activity coefficients wrt
405 //! log(moles_m)
406 /*!
407 * This function will be called to update the internally stored derivative
408 * of the natural logarithm of the activity coefficients wrt logarithm of
409 * the mole number of species
410 */
411 void s_update_dlnActCoeff_dlnN() const;
412
413protected:
414 //! number of binary interaction expressions
416
417 //! Enthalpy term for the binary mole fraction interaction of the
418 //! excess Gibbs free energy expression
419 mutable vector<double> m_HE_b_ij;
420
421 //! Enthalpy term for the ternary mole fraction interaction of the
422 //! excess Gibbs free energy expression
423 mutable vector<double> m_HE_c_ij;
424
425 //! Entropy term for the binary mole fraction interaction of the
426 //! excess Gibbs free energy expression
427 mutable vector<double> m_SE_b_ij;
428
429 //! Entropy term for the ternary mole fraction interaction of the
430 //! excess Gibbs free energy expression
431 mutable vector<double> m_SE_c_ij;
432
433 //! Enthalpy term for the binary mole fraction interaction of the
434 //! excess Gibbs free energy expression
435 mutable vector<double> m_VHE_b_ij;
436
437 //! Enthalpy term for the ternary mole fraction interaction of the
438 //! excess Gibbs free energy expression
439 mutable vector<double> m_VHE_c_ij;
440
441 //! Entropy term for the binary mole fraction interaction of the
442 //! excess Gibbs free energy expression
443 mutable vector<double> m_VSE_b_ij;
444
445 //! Entropy term for the ternary mole fraction interaction of the
446 //! excess Gibbs free energy expression
447 mutable vector<double> m_VSE_c_ij;
448
449 //! vector of species indices representing species A in the interaction
450 /*!
451 * Each Margules excess Gibbs free energy term involves two species, A and
452 * B. This vector identifies species A.
453 */
454 vector<size_t> m_pSpecies_A_ij;
455
456 //! vector of species indices representing species B in the interaction
457 /*!
458 * Each Margules excess Gibbs free energy term involves two species, A and
459 * B. This vector identifies species B.
460 */
461 vector<size_t> m_pSpecies_B_ij;
462
463 //! form of the Margules interaction expression
464 /*!
465 * Currently there is only one form.
466 */
468
469 //! form of the temperature dependence of the Margules interaction expression
470 /*!
471 * Currently there is only one form -> constant wrt temperature.
472 */
474};
475
476}
477
478#endif
Header for intermediate ThermoPhase object for phases which employ Gibbs excess free energy based for...
A map of string keys to values whose type can vary at runtime.
Definition AnyMap.h:432
GibbsExcessVPSSTP is a derived class of ThermoPhase that handles variable pressure standard state met...
MargulesVPSSTP is a derived class of GibbsExcessVPSSTP that employs the Margules approximation for th...
void getdlnActCoeffds(const double dTds, const double *const dXds, double *dlnActCoeffds) const override
Get the change in activity coefficients wrt changes in state (temp, mole fraction,...
void getPartialMolarEnthalpies(double *hbar) const override
Returns an array of partial molar enthalpies for the species in the mixture.
void getChemPotentials(double *mu) const override
Get the species chemical potentials. Units: J/kmol.
vector< double > m_VHE_c_ij
Enthalpy term for the ternary mole fraction interaction of the excess Gibbs free energy expression.
vector< double > m_SE_b_ij
Entropy term for the binary mole fraction interaction of the excess Gibbs free energy expression.
int formTempModel_
form of the temperature dependence of the Margules interaction expression
size_t numBinaryInteractions_
number of binary interaction expressions
vector< double > m_SE_c_ij
Entropy term for the ternary mole fraction interaction of the excess Gibbs free energy expression.
void s_update_dlnActCoeff_dlnN_diag() const
Update the derivative of the log of the activity coefficients wrt log(moles) - diagonal only.
void getParameters(AnyMap &phaseNode) const override
Store the parameters of a ThermoPhase object such that an identical one could be reconstructed using ...
string type() const override
String indicating the thermodynamic model implemented.
void initThermo() override
Initialize the ThermoPhase object after all species have been set up.
void getPartialMolarVolumes(double *vbar) const override
Return an array of partial molar volumes for the species in the mixture.
vector< size_t > m_pSpecies_A_ij
vector of species indices representing species A in the interaction
double cv_mole() const override
Molar heat capacity at constant volume. Units: J/kmol/K.
void s_update_dlnActCoeff_dT() const
Update the derivative of the log of the activity coefficients wrt T.
vector< size_t > m_pSpecies_B_ij
vector of species indices representing species B in the interaction
vector< double > m_VSE_c_ij
Entropy term for the ternary mole fraction interaction of the excess Gibbs free energy expression.
void s_update_dlnActCoeff_dlnN() const
Update the derivative of the log of the activity coefficients wrt log(moles_m)
int formMargules_
form of the Margules interaction expression
vector< double > m_HE_b_ij
Enthalpy term for the binary mole fraction interaction of the excess Gibbs free energy expression.
void getdlnActCoeffdT(double *dlnActCoeffdT) const override
Get the array of temperature derivatives of the log activity coefficients.
vector< double > m_VSE_b_ij
Entropy term for the binary mole fraction interaction of the excess Gibbs free energy expression.
void getPartialMolarCp(double *cpbar) const override
Returns an array of partial molar entropies for the species in the mixture.
void initLengths()
Initialize lengths of local variables after all species have been identified.
void getLnActivityCoefficients(double *lnac) const override
Get the array of non-dimensional molar-based ln activity coefficients at the current solution tempera...
void s_update_dlnActCoeff_dlnX_diag() const
Update the derivative of the log of the activity coefficients wrt log(mole fraction)
void s_update_lnActCoeff() const
Update the activity coefficients.
void getdlnActCoeffdlnX_diag(double *dlnActCoeffdlnX_diag) const override
Get the array of ln mole fraction derivatives of the log activity coefficients - diagonal component o...
vector< double > m_VHE_b_ij
Enthalpy term for the binary mole fraction interaction of the excess Gibbs free energy expression.
void addBinaryInteraction(const string &speciesA, const string &speciesB, double h0, double h1, double s0, double s1, double vh0, double vh1, double vs0, double vs1)
Add a binary species interaction with the specified parameters.
vector< double > m_HE_c_ij
Enthalpy term for the ternary mole fraction interaction of the excess Gibbs free energy expression.
void getdlnActCoeffdlnN_diag(double *dlnActCoeffdlnN_diag) const override
Get the array of log species mole number derivatives of the log activity coefficients.
void getd2lnActCoeffdT2(double *d2lnActCoeffdT2) const
Get the array of temperature second derivatives of the log activity coefficients.
void getPartialMolarEntropies(double *sbar) const override
Returns an array of partial molar entropies for the species in the mixture.
void getdlnActCoeffdlnN(const size_t ld, double *const dlnActCoeffdlnN) override
Get the array of derivatives of the log activity coefficients with respect to the log of the species ...
Namespace for the Cantera kernel.
Definition AnyMap.cpp:595