Cantera  2.0
IdealSolidSolnPhase.h
Go to the documentation of this file.
1 /**
2  * @file IdealSolidSolnPhase.h
3  * Header file for an ideal solid solution model
4  * with incompressible thermodynamics (see \ref thermoprops and
5  * \link Cantera::IdealSolidSolnPhase IdealSolidSolnPhase\endlink).
6  *
7  * This class inherits from the Cantera class ThermoPhase
8  * and implements an ideal solid solution model with incompressible
9  * thermodynamics.
10  */
11 /*
12  * Copyright 2006 Sandia Corporation. Under the terms of Contract
13  * DE-AC04-94AL85000, with Sandia Corporation, the U.S. Government
14  * retains certain rights in this software.
15  */
16 
17 #ifndef CT_IDEALSOLIDSOLNPHASE_H
18 #define CT_IDEALSOLIDSOLNPHASE_H
19 
20 #include "mix_defs.h"
21 #include "ThermoPhase.h"
22 #include "ThermoFactory.h"
23 #include "SpeciesThermo.h"
24 
25 
26 namespace Cantera
27 {
28 
29 /*!
30  * @name CONSTANTS - Models for the Standard State of IdealSolidSolnPhase's
31  */
32 //@{
33 const int cIdealSolidSolnPhase0 = 5010;
34 const int cIdealSolidSolnPhase1 = 5011;
35 const int cIdealSolidSolnPhase2 = 5012;
36 //@}
37 
38 /**
39  * Class IdealSolidSolnPhase represents a condensed phase ideal
40  * solution compound. The phase and the pure species phases which
41  * comprise the standard states of the species are assumed to have
42  * zero volume expansivity and zero isothermal compressibility.
43  * Each species does, however, have constant but distinct partial
44  * molar volumes equal to their pure species molar volumes.
45  * The class derives from class ThermoPhase,
46  * and overloads the virtual methods defined there with ones that
47  * use expressions appropriate for ideal solution mixtures.
48  * File name for the XML datafile containing information
49  * for this phase
50  * The generalized concentrations can have three different forms
51  * depending on the value of the member attribute m_formGC, which
52  * is supplied in the constructor and in the XML file.
53  * <TABLE>
54  * <TR><TD> m_formGC </TD><TD> GeneralizedConc </TD><TD> StandardConc </TD></TR>
55  * <TR><TD> 0 </TD><TD> X_k </TD><TD> 1.0 </TD></TR>
56  * <TR><TD> 1 </TD><TD> X_k / V_k </TD><TD> 1.0 / V_k </TD></TR>
57  * <TR><TD> 2 </TD><TD> X_k / V_N </TD><TD> 1.0 / V_N </TD></TR>
58  * </TABLE>
59  * The value and form of the generalized concentration will affect
60  * reaction rate constants involving species in this phase.
61  *
62  * @ingroup thermoprops
63  */
65 {
66 
67 public:
68 
69  /**
70  * Constructor for IdealSolidSolnPhase.
71  * The generalized concentrations can have three different forms
72  * depending on the value of the member attribute m_formGC, which
73  * is supplied in the constructor or read from the xml data file.
74  * <TABLE>
75  * <TR><TD> m_formGC </TD><TD> GeneralizedConc </TD><TD> StandardConc </TD></TR>
76  * <TR><TD> 0 </TD><TD> X_k </TD><TD> 1.0 </TD></TR>
77  * <TR><TD> 1 </TD><TD> X_k / V_k </TD><TD> 1.0 / V_k </TD></TR>
78  * <TR><TD> 2 </TD><TD> X_k / V_N </TD><TD> 1.0 / V_N </TD></TR>
79  * </TABLE>
80  *
81  * @param formCG This parameter initializes the m_formGC variable. The default
82  * is a value of 0.
83  */
84  IdealSolidSolnPhase(int formCG=0);
85 
86 
87  //! Construct and initialize an IdealSolidSolnPhase ThermoPhase object
88  //! directly from an ASCII input file
89  /*!
90  *
91  * This constructor will also fully initialize the object.
92  * The generalized concentrations can have three different forms
93  * depending on the value of the member attribute m_formGC, which
94  * is supplied in the constructor or read from the xml data file.
95  *
96  * <TABLE>
97  * <TR><TD> m_formGC </TD><TD> GeneralizedConc </TD><TD> StandardConc </TD></TR>
98  * <TR><TD> 0 </TD><TD> X_k </TD><TD> 1.0 </TD></TR>
99  * <TR><TD> 1 </TD><TD> X_k / V_k </TD><TD> 1.0 / V_k </TD></TR>
100  * <TR><TD> 2 </TD><TD> X_k / V_N </TD><TD> 1.0 / V_N </TD></TR>
101  * </TABLE>
102  *
103  * @param infile File name for the XML datafile containing information
104  * for this phase
105  * @param id The name of this phase. This is used to look up
106  * the phase in the XML datafile.
107  * @param formCG This parameter initializes the m_formGC variable. The default
108  * is a value of 0.
109  */
110  IdealSolidSolnPhase(std::string infile, std::string id="", int formCG=0);
111 
112  //! Construct and initialize an IdealSolidSolnPhase ThermoPhase object
113  //! directly from an XML database
114  /*!
115  * The generalized concentrations can have three different forms
116  * depending on the value of the member attribute m_formGC, which
117  * is supplied in the constructor and/or read from the data file.
118  *
119  * <TABLE>
120  * <TR><TD> m_formGC </TD><TD> GeneralizedConc </TD><TD> StandardConc </TD></TR>
121  * <TR><TD> 0 </TD><TD> X_k </TD><TD> 1.0 </TD></TR>
122  * <TR><TD> 1 </TD><TD> X_k / V_k </TD><TD> 1.0 / V_k </TD></TR>
123  * <TR><TD> 2 </TD><TD> X_k / V_N </TD><TD> 1.0 / V_N </TD></TR>
124  * </TABLE>
125  *
126  * @param root XML tree containing a description of the phase.
127  * The tree must be positioned at the XML element
128  * named phase with id, "id", on input to this routine.
129  * @param id The name of this phase. This is used to look up
130  * the phase in the XML datafile.
131  * @param formCG This parameter initializes the m_formGC variable. The default
132  * is a value of 0.
133  */
134  IdealSolidSolnPhase(XML_Node& root, std::string id="", int formCG=0);
135 
136  /*!
137  * Copy Constructor
138  */
140 
141  /*!
142  * Assignment operator
143  */
145 
146  /*!
147  * Base Class Duplication Function
148  * -> given a pointer to ThermoPhase, this function can
149  * duplicate the object. (note has to be a separate function
150  * not the copy constructor, because it has to be
151  * a virtual function)
152  */
153  virtual ThermoPhase* duplMyselfAsThermoPhase() const;
154 
155  //! Destructor
156  virtual ~IdealSolidSolnPhase() {}
157 
158  /**
159  * Equation of state flag. Returns a value depending upon the value of
160  * m_formGC, which is defined at instantiation.
161  */
162  virtual int eosType() const;
163 
164  /**
165  * @name Molar Thermodynamic Properties of the Solution ------------------------
166  * @{
167  */
168 
169  /**
170  * Molar enthalpy of the solution. Units: J/kmol.
171  * For an ideal, constant partial molar volume solution mixture with
172  * pure species phases which exhibit zero volume expansivity and
173  * zero isothermal compressibility:
174  * \f[
175  * \hat h(T,P) = \sum_k X_k \hat h^0_k(T) + (P - P_{ref}) (\sum_k X_k \hat V^0_k)
176  * \f]
177  * The reference-state pure-species enthalpies at the reference pressure Pref
178  * \f$ \hat h^0_k(T) \f$, are computed by the species thermodynamic
179  * property manager. They are polynomial functions of temperature.
180  * @see SpeciesThermo
181  */
182  virtual doublereal enthalpy_mole() const;
183 
184  /**
185  * Molar internal energy of the solution. Units: J/kmol.
186  * For an ideal, constant partial molar volume solution mixture with
187  * pure species phases which exhibit zero volume expansivity and
188  * zero isothermal compressibility:
189  * \f[
190  * \hat u(T,X) = \hat h(T,P,X) - p \hat V
191  * = \sum_k X_k \hat h^0_k(T) - P_{ref} (\sum_k{X_k \hat V^0_k})
192  * \f]
193  * and is a function only of temperature.
194  * The reference-state pure-species enthalpies
195  * \f$ \hat h^0_k(T) \f$ are computed by the species thermodynamic
196  * property manager.
197  * @see SpeciesThermo
198  */
199  virtual doublereal intEnergy_mole() const;
200 
201  /**
202  * Molar entropy of the solution. Units: J/kmol/K.
203  * For an ideal, constant partial molar volume solution mixture with
204  * pure species phases which exhibit zero volume expansivity:
205  * \f[
206  * \hat s(T, P, X_k) = \sum_k X_k \hat s^0_k(T) - \hat R \sum_k X_k log(X_k)
207  * \f]
208  * The reference-state pure-species entropies
209  * \f$ \hat s^0_k(T,p_{ref}) \f$ are computed by the species thermodynamic
210  * property manager. The pure species entropies are independent of
211  * pressure since the volume expansivities are equal to zero.
212  * @see SpeciesThermo
213  */
214  virtual doublereal entropy_mole() const;
215 
216  /**
217  * Molar gibbs free energy of the solution. Units: J/kmol.
218  * For an ideal, constant partial molar volume solution mixture with
219  * pure species phases which exhibit zero volume expansivity:
220  * \f[
221  * \hat g(T, P) = \sum_k X_k \hat g^0_k(T,P) + \hat R T \sum_k X_k log(X_k)
222  * \f]
223  * The reference-state pure-species gibbs free energies
224  * \f$ \hat g^0_k(T) \f$ are computed by the species thermodynamic
225  * property manager, while the standard state gibbs free energies
226  * \f$ \hat g^0_k(T,P) \f$ are computed by the member function, gibbs_RT().
227  * @see SpeciesThermo
228  */
229  virtual doublereal gibbs_mole() const;
230 
231  /**
232  * Molar heat capacity at constant pressure of the solution.
233  * Units: J/kmol/K.
234  * For an ideal, constant partial molar volume solution mixture with
235  * pure species phases which exhibit zero volume expansivity:
236  * \f[
237  * \hat c_p(T,P) = \sum_k X_k \hat c^0_{p,k}(T) .
238  * \f]
239  * The heat capacity is independent of pressure.
240  * The reference-state pure-species heat capacities
241  * \f$ \hat c^0_{p,k}(T) \f$ are computed by the species thermodynamic
242  * property manager.
243  * @see SpeciesThermo
244  */
245  virtual doublereal cp_mole() const;
246 
247  /**
248  * Molar heat capacity at constant volume of the solution.
249  * Units: J/kmol/K.
250  * For an ideal, constant partial molar volume solution mixture with
251  * pure species phases which exhibit zero volume expansivity:
252  * \f[ \hat c_v(T,P) = \hat c_p(T,P) \f]
253  * The two heat capacities are equal.
254  */
255  virtual doublereal cv_mole() const {
256  return cp_mole();
257  }
258 
259  //@}
260  /** @name Mechanical Equation of State Properties ------------------------------------
261  *
262  * In this equation of state implementation, the density is a
263  * function only of the mole fractions. Therefore, it can't be
264  * an independent variable. Instead, the pressure is used as the
265  * independent variable. Functions which try to set the thermodynamic
266  * state by calling setDensity() may cause an exception to be
267  * thrown.
268  */
269  //@{
270 
271  /**
272  * Pressure. Units: Pa.
273  * For this incompressible system, we return the internally stored
274  * independent value of the pressure.
275  */
276  virtual doublereal pressure() const {
277  return m_Pcurrent;
278  }
279 
280  /**
281  * Set the pressure at constant temperature. Units: Pa.
282  * This method sets a constant within the object.
283  * The mass density is not a function of pressure.
284  *
285  * @param p Input Pressure (Pa)
286  */
287  virtual void setPressure(doublereal p);
288 
289  /**
290  * Calculate the density of the mixture using the partial
291  * molar volumes and mole fractions as input
292  *
293  * The formula for this is
294  *
295  * \f[
296  * \rho = \frac{\sum_k{X_k W_k}}{\sum_k{X_k V_k}}
297  * \f]
298  *
299  * where \f$X_k\f$ are the mole fractions, \f$W_k\f$ are
300  * the molecular weights, and \f$V_k\f$ are the pure species
301  * molar volumes.
302  *
303  * Note, the basis behind this formula is that in an ideal
304  * solution the partial molar volumes are equal to the pure
305  * species molar volumes. We have additionally specified
306  * in this class that the pure species molar volumes are
307  * independent of temperature and pressure.
308  *
309  * NOTE: This is a non-virtual function, which is not a
310  * member of the ThermoPhase base class.
311  */
312  void calcDensity();
313 
314  /**
315  * Overwritten setDensity() function is necessary because the
316  * density is not an independent variable.
317  *
318  * This function will now throw an error condition
319  *
320  * @internal May have to adjust the strategy here to make
321  * the eos for these materials slightly compressible, in order
322  * to create a condition where the density is a function of
323  * the pressure.
324  *
325  * This function will now throw an error condition.
326  *
327  * NOTE: This is a virtual function that overwrites the State.h
328  * class
329  *
330  * @param rho Input density
331  */
332  virtual void setDensity(const doublereal rho);
333 
334  /**
335  * Overwritten setMolarDensity() function is necessary because the
336  * density is not an independent variable.
337  *
338  * This function will now throw an error condition.
339  *
340  * NOTE: This is virtual function that overwrites the State.h
341  * class
342  *
343  * @param rho Input Density
344  */
345  virtual void setMolarDensity(const doublereal rho);
346 
347  //! Set the mole fractions
348  /*!
349  * @param x Input vector of mole fractions.
350  * Length: m_kk.
351  */
352  virtual void setMoleFractions(const doublereal* const x);
353 
354  //! Set the mole fractions, but don't normalize them to one.
355  /*!
356  * @param x Input vector of mole fractions.
357  * Length: m_kk.
358  */
359  virtual void setMoleFractions_NoNorm(const doublereal* const x);
360 
361  //! Set the mass fractions, and normalize them to one.
362  /*!
363  * @param y Input vector of mass fractions.
364  * Length: m_kk.
365  */
366  virtual void setMassFractions(const doublereal* const y);
367 
368  //! Set the mass fractions, but don't normalize them to one
369  /*!
370  * @param y Input vector of mass fractions.
371  * Length: m_kk.
372  */
373  virtual void setMassFractions_NoNorm(const doublereal* const y);
374 
375  //! Set the concentration,
376  /*!
377  * @param c Input vector of concentrations.
378  * Length: m_kk.
379  */
380  virtual void setConcentrations(const doublereal* const c);
381 
382 
383  //@}
384 
385  /**
386  * @name Chemical Potentials and Activities -----------------------------------------
387  *
388  * The activity \f$a_k\f$ of a species in solution is
389  * related to the chemical potential by
390  * \f[
391  * \mu_k(T,P,X_k) = \mu_k^0(T,P)
392  * + \hat R T \log a_k.
393  * \f]
394  * The quantity \f$\mu_k^0(T,P)\f$ is
395  * the standard state chemical potential at unit activity.
396  * It may depend on the pressure and the temperature. However,
397  * it may not depend on the mole fractions of the species
398  * in the solid solution.
399  *
400  * The activities are related to the generalized
401  * concentrations, \f$\tilde C_k\f$, and standard
402  * concentrations, \f$C^0_k\f$, by the following formula:
403  *
404  * \f[
405  * a_k = \frac{\tilde C_k}{C^0_k}
406  * \f]
407  * The generalized concentrations are used in the kinetics classes
408  * to describe the rates of progress of reactions involving the
409  * species. Their formulation depends upon the specification
410  * of the rate constants for reaction, especially the units used
411  * in specifying the rate constants. The bridge between the
412  * thermodynamic equilibrium expressions that use a_k and the
413  * kinetics expressions which use the generalized concentrations
414  * is provided by the multiplicative factor of the
415  * standard concentrations.
416  * @{
417  */
418 
419  /**
420  * This method returns the array of generalized
421  * concentrations. The generalized concentrations are used
422  * in the evaluation of the rates of progress for reactions
423  * involving species in this phase. The generalized
424  * concentration divided by the standard concentration is also
425  * equal to the activity of species.
426  *
427  * For this implentation the activity is defined to be the
428  * mole fraction of the species. The generalized concentration
429  * is defined to be equal to the mole fraction divided by
430  * the partial molar volume. The generalized concentrations
431  * for species in this phase therefore have units of
432  * kmol m<SUP>-3</SUP>. Rate constants must reflect this fact.
433  *
434  * On a general note, the following must be true.
435  * For an ideal solution, the generalized concentration must consist
436  * of the mole fraction multiplied by a constant. The constant may be
437  * fairly arbitrarily chosen, with differences adsorbed into the
438  * reaction rate expression. 1/V_N, 1/V_k, or 1 are equally good,
439  * as long as the standard concentration is adjusted accordingly.
440  * However, it must be a constant (and not the concentration, btw,
441  * which is a function of the mole fractions) in order for the
442  * ideal solution properties to hold at the same time having the
443  * standard concentration to be independent of the mole fractions.
444  *
445  * In this implementation the form of the generalized concentrations
446  * depend upon the member attribute, m_formGC:
447  *
448  * <TABLE>
449  * <TR><TD> m_formGC </TD><TD> GeneralizedConc </TD><TD> StandardConc </TD></TR>
450  * <TR><TD> 0 </TD><TD> X_k </TD><TD> 1.0 </TD></TR>
451  * <TR><TD> 1 </TD><TD> X_k / V_k </TD><TD> 1.0 / V_k </TD></TR>
452  * <TR><TD> 2 </TD><TD> X_k / V_N </TD><TD> 1.0 / V_N </TD></TR>
453  * </TABLE>
454  *
455  * HKM Note: We have absorbed the pressure dependence of the pure species
456  * state into the thermodynamics functions. Therefore the
457  * standard state on which the activities are based depend
458  * on both temperature and pressure. If we hadn't, it would have
459  * appeared in this function in a very awkward exp[] format.
460  *
461  * @param c Pointer to array of doubles of length m_kk, which on exit
462  * will contain the generalized concentrations.
463  */
464  virtual void getActivityConcentrations(doublereal* c) const;
465 
466  /**
467  * The standard concentration \f$ C^0_k \f$ used to normalize
468  * the generalized concentration.
469  * In many cases, this quantity
470  * will be the same for all species in a phase.
471  * However, for this case, we will return a distinct concentration
472  * for each species. This is the inverse of the species molar
473  * volume. Units for the standard concentration are
474  * kmol m<SUP>-3</SUP>.
475  *
476  * @param k Species number: this is a require parameter,
477  * a change from the ThermoPhase base class, where it was
478  * an optional parameter.
479  */
480  virtual doublereal standardConcentration(size_t k) const;
481 
482  /**
483  * The reference (ie standard) concentration \f$ C^0_k \f$ used to normalize
484  * the generalized concentration. In many cases, this quantity
485  * will be the same for all species in a phase.
486  * However, for this case, we will return a distinct concentration
487  * for each species. (clone of the standard concentration ->
488  * suggest changing the name). This is the inverse of the species molar
489  * volume.
490  *
491  * @param k Species index.
492  */
493  virtual doublereal referenceConcentration(int k) const;
494 
495  /**
496  * Returns the log of the standard concentration of the kth species
497  *
498  * @param k Species number: this is a require parameter,
499  * a change from the ThermoPhase base class, where it was
500  * an optional parameter.
501  */
502  virtual doublereal logStandardConc(size_t k) const;
503 
504  /**
505  * Returns the units of the standard and general concentrations
506  * Note they have the same units, as their divisor is
507  * defined to be equal to the activity of the kth species
508  * in the solution, which is unitless.
509  *
510  * This routine is used in print out applications where the
511  * units are needed. Usually, MKS units are assumed throughout
512  * the program and in the XML input files.
513  *
514  * @param uA Output vector containing the units
515  * uA[0] = kmol units - default = 1
516  * uA[1] = m units - default = -nDim(), the number of spatial
517  * dimensions in the Phase class.
518  * uA[2] = kg units - default = 0;
519  * uA[3] = Pa(pressure) units - default = 0;
520  * uA[4] = Temperature units - default = 0;
521  * uA[5] = time units - default = 0
522  * @param k species index. Defaults to 0.
523  * @param sizeUA output int containing the size of the vector.
524  * Currently, this is equal to 6.
525  *
526  * For EOS types other than cIdealSolidSolnPhase0, the default
527  * kmol/m3 holds for standard concentration units. For
528  * cIdealSolidSolnPhase0 type, the standard concentration is
529  * unitless.
530  */
531  virtual void getUnitsStandardConc(double* uA, int k = 0,
532  int sizeUA = 6) const;
533 
534 
535  //! Get the array of species activity coefficients
536  /*!
537  * @param ac output vector of activity coefficients. Length: m_kk
538  */
539  virtual void getActivityCoefficients(doublereal* ac) const;
540 
541  /**
542  * Get the species chemical potentials. Units: J/kmol.
543  *
544  * This function returns a vector of chemical potentials of the
545  * species in solution.
546  * \f[
547  * \mu_k = \mu^{ref}_k(T) + V_k * (p - p_o) + R T ln(X_k)
548  * \f]
549  * or another way to phrase this is
550  * \f[
551  * \mu_k = \mu^o_k(T,p) + R T ln(X_k)
552  * \f]
553  * where \f$ \mu^o_k(T,p) = \mu^{ref}_k(T) + V_k * (p - p_o)\f$
554  *
555  * @param mu Output vector of chemical potentials.
556  */
557  virtual void getChemPotentials(doublereal* mu) const;
558 
559  /**
560  * Get the array of non-dimensional species solution
561  * chemical potentials at the current T and P
562  * \f$\mu_k / \hat R T \f$.
563  * \f[
564  * \mu^0_k(T,P) = \mu^{ref}_k(T) + (P - P_{ref}) * V_k + RT ln(X_k)
565  * \f]
566  * where \f$V_k\f$ is the molar volume of pure species <I>k</I>.
567  * \f$ \mu^{ref}_k(T)\f$ is the chemical potential of pure
568  * species <I>k</I> at the reference pressure, \f$P_{ref}\f$.
569  *
570  * @param mu Output vector of dimensionless chemical potentials. Length = m_kk.
571  */
572  virtual void getChemPotentials_RT(doublereal* mu) const;
573 
574  //@}
575  /// @name Partial Molar Properties of the Solution -----------------------------
576  //@{
577 
578 
579  //! Returns an array of partial molar enthalpies for the species in the mixture.
580  /*!
581  * Units (J/kmol)
582  * For this phase, the partial molar enthalpies are equal to the
583  * pure species enthalpies
584  * \f[
585  * \bar h_k(T,P) = \hat h^{ref}_k(T) + (P - P_{ref}) \hat V^0_k
586  * \f]
587  * The reference-state pure-species enthalpies, \f$ \hat h^{ref}_k(T) \f$,
588  * at the reference pressure,\f$ P_{ref} \f$,
589  * are computed by the species thermodynamic
590  * property manager. They are polynomial functions of temperature.
591  * @see SpeciesThermo
592  *
593  * @param hbar Output vector containing partial molar enthalpies.
594  * Length: m_kk.
595  */
596  virtual void getPartialMolarEnthalpies(doublereal* hbar) const;
597 
598  /**
599  * Returns an array of partial molar entropies of the species in the
600  * solution. Units: J/kmol/K.
601  * For this phase, the partial molar entropies are equal to the
602  * pure species entropies plus the ideal solution contribution.
603  * \f[
604  * \bar s_k(T,P) = \hat s^0_k(T) - R log(X_k)
605  * \f]
606  * The reference-state pure-species entropies,\f$ \hat s^{ref}_k(T) \f$,
607  * at the reference pressure, \f$ P_{ref} \f$, are computed by the
608  * species thermodynamic
609  * property manager. They are polynomial functions of temperature.
610  * @see SpeciesThermo
611  *
612  * @param sbar Output vector containing partial molar entropies.
613  * Length: m_kk.
614  */
615  virtual void getPartialMolarEntropies(doublereal* sbar) const;
616 
617  /**
618  * Returns an array of partial molar Heat Capacities at constant
619  * pressure of the species in the
620  * solution. Units: J/kmol/K.
621  * For this phase, the partial molar heat capacities are equal
622  * to the standard state heat capacities.
623  *
624  * @param cpbar Output vector of partial heat capacities. Length: m_kk.
625  */
626  virtual void getPartialMolarCp(doublereal* cpbar) const;
627 
628  /**
629  * returns an array of partial molar volumes of the species
630  * in the solution. Units: m^3 kmol-1.
631  *
632  * For this solution, thepartial molar volumes are equal to the
633  * constant species molar volumes.
634  *
635  * @param vbar Output vector of partial molar volumes. Length: m_kk.
636  */
637  virtual void getPartialMolarVolumes(doublereal* vbar) const;
638 
639  //@}
640  /// @name Properties of the Standard State of the Species in the Solution -------------------------------------
641  //@{
642 
643 
644  /**
645  * Get the standard state chemical potentials of the species.
646  * This is the array of chemical potentials at unit activity
647  * \f$ \mu^0_k(T,P) \f$.
648  * We define these here as the chemical potentials of the pure
649  * species at the temperature and pressure of the solution.
650  * This function is used in the evaluation of the
651  * equilibrium constant Kc. Therefore, Kc will also depend
652  * on T and P. This is the norm for liquid and solid systems.
653  *
654  * units = J / kmol
655  *
656  * @param mu0 Output vector of standard state chemical potentials.
657  * Length: m_kk.
658  */
659  virtual void getStandardChemPotentials(doublereal* mu0) const {
660  getPureGibbs(mu0);
661  }
662 
663 
664  //! Get the array of nondimensional Enthalpy functions for the standard state species
665  //! at the current <I>T</I> and <I>P</I> of the solution.
666  /*!
667  * We assume an incompressible constant partial molar
668  * volume here:
669  * \f[
670  * h^0_k(T,P) = h^{ref}_k(T) + (P - P_{ref}) * V_k
671  * \f]
672  * where \f$V_k\f$ is the molar volume of pure species <I>k</I>.
673  * \f$ h^{ref}_k(T)\f$ is the enthalpy of the pure
674  * species <I>k</I> at the reference pressure, \f$P_{ref}\f$.
675  *
676  * @param hrt Vector of length m_kk, which on return hrt[k]
677  * will contain the nondimensional
678  * standard state enthalpy of species k.
679  */
680  void getEnthalpy_RT(doublereal* hrt) const;
681 
682 
683  //! Get the nondimensional Entropies for the species
684  //! standard states at the current T and P of the solution.
685  /*!
686  * Note, this is equal to the reference state entropies
687  * due to the zero volume expansivity:
688  * i.e., (dS/dP)_T = (dV/dT)_P = 0.0
689  *
690  * @param sr Vector of length m_kk, which on return sr[k]
691  * will contain the nondimensional
692  * standard state entropy for species k.
693  */
694  void getEntropy_R(doublereal* sr) const;
695 
696  /**
697  * Get the nondimensional gibbs function for the species
698  * standard states at the current T and P of the solution.
699  *
700  * \f[
701  * \mu^0_k(T,P) = \mu^{ref}_k(T) + (P - P_{ref}) * V_k
702  * \f]
703  * where \f$V_k\f$ is the molar volume of pure species <I>k</I>.
704  * \f$ \mu^{ref}_k(T)\f$ is the chemical potential of pure
705  * species <I>k</I> at the reference pressure, \f$P_{ref}\f$.
706  *
707  * @param grt Vector of length m_kk, which on return sr[k]
708  * will contain the nondimensional
709  * standard state gibbs function for species k.
710  */
711  virtual void getGibbs_RT(doublereal* grt) const;
712 
713  /**
714  * Get the Gibbs functions for the pure species
715  * at the current <I>T</I> and <I>P</I> of the solution.
716  * We assume an incompressible constant partial molar
717  * volume here:
718  * \f[
719  * \mu^0_k(T,P) = \mu^{ref}_k(T) + (P - P_{ref}) * V_k
720  * \f]
721  * where \f$V_k\f$ is the molar volume of pure species <I>k</I>.
722  * \f$ \mu^{ref}_k(T)\f$ is the chemical potential of pure
723  * species <I>k</I> at the reference pressure, \f$P_{ref}\f$.
724  *
725  * @param gpure Output vector of Gibbs functions for species
726  * Length: m_kk.
727  */
728  virtual void getPureGibbs(doublereal* gpure) const;
729 
730 
731  //! Returns the vector of nondimensional
732  //! internal Energies of the standard state at the current
733  //! temperature and pressure of the solution for each species.
734  /*!
735  *
736  * @param urt Output vector of standard state nondimensional internal energies.
737  * Length: m_kk.
738  */
739  virtual void getIntEnergy_RT(doublereal* urt) const;
740 
741  /**
742  * Get the nondimensional heat capacity at constant pressure
743  * function for the species
744  * standard states at the current T and P of the solution.
745  * \f[
746  * Cp^0_k(T,P) = Cp^{ref}_k(T)
747  * \f]
748  * where \f$V_k\f$ is the molar volume of pure species <I>k</I>.
749  * \f$ Cp^{ref}_k(T)\f$ is the constant pressure heat capacity
750  * of species <I>k</I> at the reference pressure, \f$p_{ref}\f$.
751  *
752  * @param cpr Vector of length m_kk, which on return cpr[k]
753  * will contain the nondimensional
754  * constant pressure heat capacity for species k.
755  */
756  void getCp_R(doublereal* cpr) const;
757 
758  /**
759  * Get the molar volumes of each species in their standard
760  * states at the current
761  * <I>T</I> and <I>P</I> of the solution.
762  * units = m^3 / kmol
763  *
764  * @param vol Output vector of standard state volumes.
765  * Length: m_kk.
766  */
767  virtual void getStandardVolumes(doublereal* vol) const;
768 
769 
770  //@}
771  /// @name Thermodynamic Values for the Species Reference States ------
772  //@{
773 
774 
775  /**
776  * Returns the vector of nondimensional
777  * enthalpies of the reference state at the current temperature
778  * of the solution and the reference pressure for the species.
779  *
780  * @param hrt Output vector containing reference nondimensional enthalpies.
781  * Length: m_kk.
782  */
783  virtual void getEnthalpy_RT_ref(doublereal* hrt) const;
784 
785  /**
786  * Returns the vector of nondimensional
787  * enthalpies of the reference state at the current temperature
788  * of the solution and the reference pressure for the species.
789  *
790  * @param grt Output vector containing reference nondimensional Gibbs free energies.
791  * Length: m_kk.
792  */
793  virtual void getGibbs_RT_ref(doublereal* grt) const;
794 
795  /**
796  * Returns the vector of the
797  * gibbs function of the reference state at the current temperature
798  * of the solution and the reference pressure for the species.
799  * units = J/kmol
800  *
801  * @param g Output vector containing reference Gibbs free energies.
802  * Length: m_kk.
803  */
804  virtual void getGibbs_ref(doublereal* g) const;
805 
806  /**
807  * Returns the vector of nondimensional
808  * entropies of the reference state at the current temperature
809  * of the solution and the reference pressure for the species.
810  *
811  * @param er Output vector containing reference nondimensional entropies.
812  * Length: m_kk.
813  */
814  virtual void getEntropy_R_ref(doublereal* er) const;
815 
816  /**
817  * Returns the vector of nondimensional
818  * internal Energies of the reference state at the current temperature
819  * of the solution and the reference pressure for each species.
820  *
821  * @param urt Output vector containing reference nondimensional internal energies.
822  * Length: m_kk.
823  */
824  virtual void getIntEnergy_RT_ref(doublereal* urt) const;
825 
826  /**
827  * Returns the vector of nondimensional
828  * constant pressure heat capacities of the reference state
829  * at the current temperature of the solution
830  * and reference pressure for the species.
831  *
832  * @param cprt Output vector containing reference nondimensional heat capacities.
833  * Length: m_kk.
834  */
835  virtual void getCp_R_ref(doublereal* cprt) const;
836 
837  /**
838  * Returns a reference to the vector of nondimensional
839  * enthalpies of the reference state at the current temperature.
840  * Real reason for its existence is that it also checks
841  * to see if a recalculation of the reference thermodynamics
842  * functions needs to be done.
843  */
844  const vector_fp& enthalpy_RT_ref() const;
845 
846  /**
847  * Returns a reference to the vector of nondimensional
848  * enthalpies of the reference state at the current temperature.
849  * Real reason for its existence is that it also checks
850  * to see if a recalculation of the reference thermodynamics
851  * functions needs to be done.
852  */
853  const vector_fp& gibbs_RT_ref() const {
854  _updateThermo();
855  return m_g0_RT;
856  }
857 
858  /**
859  * Returns a reference to the vector of nondimensional
860  * enthalpies of the reference state at the current temperature.
861  * Real reason for its existence is that it also checks
862  * to see if a recalculation of the reference thermodynamics
863  * functions needs to be done.
864  */
865  const vector_fp& expGibbs_RT_ref() const;
866 
867  /**
868  * Returns a reference to the vector of nondimensional
869  * enthalpies of the reference state at the current temperature.
870  * Real reason for its existence is that it also checks
871  * to see if a recalculation of the reference thermodynamics
872  * functions needs to be done.
873  */
874  const vector_fp& entropy_R_ref() const;
875 
876  /**
877  * Returns a reference to the vector of nondimensional
878  * enthalpies of the reference state at the current temperature.
879  * Real reason for its existence is that it also checks
880  * to see if a recalculation of the reference thermodynamics
881  * functions needs to be done.
882  */
883  const vector_fp& cp_R_ref() const {
884  _updateThermo();
885  return m_cp0_R;
886  }
887 
888  virtual void setPotentialEnergy(int k, doublereal pe) {
889  m_pe[k] = pe;
890  _updateThermo();
891  }
892 
893  virtual doublereal potentialEnergy(int k) const {
894  return m_pe[k];
895  }
896  //@}
897  /// @name Utility Functions -----------------------------------------------
898  //@{
899 
900 
901 
902  /**
903  * Initialization of an IdealSolidSolnPhase phase using an
904  * xml file
905  *
906  * This routine is a precursor to constructPhaseXML(XML_Node*)
907  * routine, which does most of the work.
908  *
909  * @param infile XML file containing the description of the
910  * phase
911  *
912  * @param id Optional parameter identifying the name of the
913  * phase. If none is given, the first XML
914  * phase element will be used.
915  */
916  void constructPhaseFile(std::string infile, std::string id="");
917 
918  /**
919  * Import and initialize an IdealSolidSolnPhase phase
920  * specification in an XML tree into the current object.
921  * Here we read an XML description of the phase.
922  * We import descriptions of the elements that make up the
923  * species in a phase.
924  * We import information about the species, including their
925  * reference state thermodynamic polynomials. We then freeze
926  * the state of the species.
927  * This routine calls importPhase() to do most of its work.
928  * Then, importPhase() calls initThermoXML() to finish
929  * off the work.
930  *
931 
932  * @param phaseNode This object must be the phase node of a
933  * complete XML tree
934  * description of the phase, including all of the
935  * species data. In other words while "phase" must
936  * point to an XML phase object, it must have
937  * sibling nodes "speciesData" that describe
938  * the species in the phase.
939  * @param id ID of the phase. If nonnull, a check is done
940  * to see if phaseNode is pointing to the phase
941  * with the correct id.
942  */
943  void constructPhaseXML(XML_Node& phaseNode, std::string id="");
944 
945  /**
946  * Initialization of an IdealSolidSolnPhase phase:
947  * Note this function is pretty much useless because it doesn't
948  * get the xml tree passed to it. Suggest a change.
949  */
950  virtual void initThermo();
951 
952  /**
953  * @internal
954  * Import and initialize a ThermoPhase object
955  * using an XML tree.
956  * Here we read extra information about the XML description
957  * of a phase. Regular information about elements and species
958  * and their reference state thermodynamic information
959  * have already been read at this point.
960  * For example, we do not need to call this function for
961  * ideal gas equations of state.
962  * This function is called from importPhase()
963  * after the elements and the
964  * species are initialized with default ideal solution
965  * level data.
966  *
967  * @param phaseNode This object must be the phase node of a
968  * complete XML tree
969  * description of the phase, including all of the
970  * species data. In other words while "phase" must
971  * point to an XML phase object, it must have
972  * sibling nodes "speciesData" that describe
973  * the species in the phase.
974  * @param id ID of the phase. If nonnull, a check is done
975  * to see if phaseNode is pointing to the phase
976  * with the correct id.
977  */
978  virtual void initThermoXML(XML_Node& phaseNode, std::string id);
979 
980 
981  /**
982  * Set mixture to an equilibrium state consistent with specified
983  * element potentials and the temperature.
984  *
985  * @param lambda_RT vector of non-dimensional element potentials
986  * \f$ \lambda_m/RT \f$.
987  *
988  */
989  virtual void setToEquilState(const doublereal* lambda_RT);
990 
991 
992  /**
993  * Report the molar volume of species k
994  *
995  * units - \f$ m^3 kmol^-1 \f$
996  *
997  * @param k species index
998  */
999  double speciesMolarVolume(int k) const;
1000 
1001  /**
1002  * Fill in a return vector containing the species molar volumes.
1003  *
1004  * units - \f$ m^3 kmol^-1 \f$
1005  *
1006  * @param smv output vector containing species molar volumes.
1007  * Length: m_kk.
1008  */
1009  void getSpeciesMolarVolumes(doublereal* smv) const;
1010 
1011  //@}
1012 
1013 protected:
1014 
1015  /**
1016  * Format for the generalized concentrations
1017  * 0 = C_k = X_k. (default)
1018  * 1 = C_k = X_k / V_k
1019  * 2 = C_k = X_k / V_N
1020  */
1022  /**
1023  * m_mm = Number of distinct elements defined in species in this
1024  * phase
1025  */
1026  size_t m_mm;
1027 
1028  /**
1029  * Maximum temperature that this phase can accurately describe
1030  * the thermodynamics.
1031  */
1032  doublereal m_tmin;
1033 
1034  /**
1035  * Minimum temperature that this phase can accurately describe
1036  * the thermodynamics.
1037  */
1038  doublereal m_tmax;
1039  /**
1040  * Value of the reference pressure for all species in this phase.
1041  * The T dependent polynomials are evaluated at the reference
1042  * pressure. Note, because this is a single value, all species
1043  * are required to have the same reference pressure.
1044  */
1045  doublereal m_Pref;
1046 
1047  /**
1048  * m_Pcurrent = The current pressure
1049  * Since the density isn't a function of pressure, but only of the
1050  * mole fractions, we need to independently specify the pressure.
1051  * The density variable which is inherited as part of the State class,
1052  * m_dens, is always kept current whenever T, P, or X[] change.
1053  */
1054  doublereal m_Pcurrent;
1055 
1056  //! Vector of molar volumes for each species in the solution
1057  /**
1058  * Species molar volumes \f$ m^3 kmol^-1 \f$
1059  */
1061 
1062  /**
1063  * Value of the temperature at which the thermodynamics functions
1064  * for the reference state of the species were last evaluated.
1065  */
1066  mutable doublereal m_tlast;
1067 
1068  /**
1069  * Vector containing the species reference enthalpies at T = m_tlast
1070  */
1072 
1073  /**
1074  * Vector containing the species reference constant pressure
1075  * heat capacities at T = m_tlast
1076  */
1078 
1079  /**
1080  * Vector containing the species reference Gibbs functions
1081  * at T = m_tlast
1082  */
1084 
1085  /**
1086  * Vector containing the species reference entropies
1087  * at T = m_tlast
1088  */
1090 
1091  /**
1092  * Vector containing the species reference exp(-G/RT) functions
1093  * at T = m_tlast
1094  */
1096 
1097  /**
1098  * Vector of potential energies for the species.
1099  */
1100  mutable vector_fp m_pe;
1101 
1102  /**
1103  * Temporary array used in equilibrium calculations
1104  */
1105  mutable vector_fp m_pp;
1106 
1107 private:
1108  /// @name Utility Functions ------------------------------------------
1109  //@{
1110  /**
1111  * This function gets called for every call to functions in this
1112  * class. It checks to see whether the temperature has changed and
1113  * thus the reference thermodynamics functions for all of the species
1114  * must be recalculated.
1115  * If the temperature has changed, the species thermo manager is called
1116  * to recalculate G, Cp, H, and S at the current temperature.
1117  */
1118  void _updateThermo() const;
1119 
1120  /**
1121  * This internal function adjusts the lengths of arrays
1122  */
1123  void initLengths();
1124 
1125  //@}
1126 };
1127 }
1128 
1129 #endif
1130 
1131 
1132 
1133 
1134