diff --git a/Cantera/src/ThermoPhase.h b/Cantera/src/ThermoPhase.h
index f7e1d41d4..08abc04b5 100755
--- a/Cantera/src/ThermoPhase.h
+++ b/Cantera/src/ThermoPhase.h
@@ -115,7 +115,7 @@ namespace Cantera {
* - IdealGasPDSS in thermo/IdealGasPDSS.h
* - MolalityVPSSTP in thermo/MolalityVPSSTP.h
* - HMWSoln in thermo/HMWSoln.h
- * - WaterTP in thermo/WaterTP.h
+ * - WaterSSTP in thermo/WaterSSTP.h
* .
*
* @see newPhase(std::string file, std::string id) Description for how to
diff --git a/Cantera/src/thermo/IdealSolidSolnPhase.cpp b/Cantera/src/thermo/IdealSolidSolnPhase.cpp
index ffc340c16..88c2e2bd8 100644
--- a/Cantera/src/thermo/IdealSolidSolnPhase.cpp
+++ b/Cantera/src/thermo/IdealSolidSolnPhase.cpp
@@ -850,7 +850,7 @@ namespace Cantera {
copy(_s.begin(), _s.end(), sr);
}
- /**
+ /*
* Returns the vector of nondimensional
* internal Energies of the standard state at the current temperature
* of the solution and current pressure for each species.
@@ -870,7 +870,7 @@ namespace Cantera {
}
}
- /**
+ /*
* Get the nondimensional heat capacity at constant pressure
* function for the species
* standard states at the current T and P of the solution.
@@ -891,7 +891,7 @@ namespace Cantera {
copy(_cpr.begin(), _cpr.end(), cpr);
}
- /**
+ /*
* Get the molar volumes of each species in their standard
* states at the current
* T and P of the solution.
@@ -907,7 +907,7 @@ namespace Cantera {
* Thermodynamic Values for the Species Reference States
*********************************************************************/
- /**
+ /*
* Returns the vector of non-dimensional Enthalpy function
* of the reference state at the current temperature
* of the solution and the reference pressure for the species.
diff --git a/Cantera/src/thermo/IdealSolidSolnPhase.h b/Cantera/src/thermo/IdealSolidSolnPhase.h
index 3f5eeff19..762c3a96e 100644
--- a/Cantera/src/thermo/IdealSolidSolnPhase.h
+++ b/Cantera/src/thermo/IdealSolidSolnPhase.h
@@ -667,10 +667,10 @@ namespace Cantera {
getPureGibbs(mu0);
}
- /**
- * Get the array of nondimensional Enthalpy functions for the
- * standard state species
- * at the current T and P of the solution.
+
+ //! Get the array of nondimensional Enthalpy functions for the standard state species
+ //! at the current T and P of the solution.
+ /*!
* We assume an incompressible constant partial molar
* volume here:
* \f[
@@ -686,10 +686,10 @@ namespace Cantera {
*/
void getEnthalpy_RT(doublereal* hrt) const;
- /**
- * Get the nondimensional Entropies for the species
- * standard states at the current T and P of the solution.
- *
+
+ //! Get the nondimensional Entropies for the species
+ //! standard states at the current T and P of the solution.
+ /*!
* Note, this is equal to the reference state entropies
* due to the zero volume expansivity:
* i.e., (dS/dP)_T = (dV/dT)_P = 0.0
@@ -699,6 +699,7 @@ namespace Cantera {
* standard state entropy for species k.
*/
void getEntropy_R(doublereal* sr) const;
+
/**
* Get the nondimensional gibbs function for the species
* standard states at the current T and P of the solution.
@@ -733,10 +734,11 @@ namespace Cantera {
*/
virtual void getPureGibbs(doublereal* gpure) const;
- /**
- * Returns the vector of nondimensional
- * internal Energies of the standard state at the current
- * temperature and pressure of the solution for each species.
+
+ //! Returns the vector of nondimensional
+ //! internal Energies of the standard state at the current
+ //! temperature and pressure of the solution for each species.
+ /*!
*
* @param urt Output vector of standard state nondimensional internal energies.
* Length: m_kk.
diff --git a/Cantera/src/thermo/Makefile.in b/Cantera/src/thermo/Makefile.in
index 570a8a164..d03a926cb 100644
--- a/Cantera/src/thermo/Makefile.in
+++ b/Cantera/src/thermo/Makefile.in
@@ -34,16 +34,16 @@ ELECTRO_OBJ = SingleSpeciesTP.o StoichSubstanceSSTP.o \
MolalityVPSSTP.o VPStandardStateTP.o \
IdealSolidSolnPhase.o IdealMolalSoln.o \
WaterPropsIAPWSphi.o WaterPropsIAPWS.o WaterProps.o \
- PDSS.o WaterPDSS.o WaterTP.o \
+ PDSS.o WaterPDSS.o \
HMWSoln.o HMWSoln_input.o DebyeHuckel.o \
- IdealGasPDSS.o
+ IdealGasPDSS.o WaterSSTP.o
ELECTRO_H = SingleSpeciesTP.h StoichSubstanceSSTP.h \
MolalityVPSSTP.h VPStandardStateTP.h \
IdealSolidSolnPhase.h IdealMolalSoln.h \
WaterPropsIAPWSphi.h WaterPropsIAPWS.h WaterProps.h \
- PDSS.h WaterPDSS.h WaterTP.h HMWSoln.h electrolytes.h \
- DebyeHuckel.h IdealGasPDSS.o
+ PDSS.h WaterPDSS.h HMWSoln.h electrolytes.h \
+ DebyeHuckel.h IdealGasPDSS.h WaterSSTP.h
endif
ifeq ($(do_issp),1)
ISSP_OBJ = IdealSolidSolnPhase.o
diff --git a/Cantera/src/thermo/SingleSpeciesTP.cpp b/Cantera/src/thermo/SingleSpeciesTP.cpp
index f9223dbe9..e3f31a92f 100644
--- a/Cantera/src/thermo/SingleSpeciesTP.cpp
+++ b/Cantera/src/thermo/SingleSpeciesTP.cpp
@@ -144,7 +144,7 @@ namespace Cantera {
return cpbar;
}
- /**
+ /*
* cv_mole():
*
* Molar heat capacity at constant volume of the mixture.
@@ -246,6 +246,14 @@ namespace Cantera {
sbar[0] *= GasConstant;
}
+ /**
+ * Get the species partial molar Heat Capacities. Units: J/kmol K.
+ */
+ void SingleSpeciesTP::getPartialMolarCp(doublereal* cpbar) const {
+ getCp_R(cpbar);
+ cpbar[0] *= GasConstant;
+ }
+
/**
* Get the species partial molar volumes. Units: m^3/kmol.
*/
diff --git a/Cantera/src/thermo/SingleSpeciesTP.h b/Cantera/src/thermo/SingleSpeciesTP.h
index ba583222f..b8fe82641 100644
--- a/Cantera/src/thermo/SingleSpeciesTP.h
+++ b/Cantera/src/thermo/SingleSpeciesTP.h
@@ -304,6 +304,17 @@ namespace Cantera {
*/
void getPartialMolarEntropies(doublereal* sbar) const;
+ //! Get the species partial molar heat capacties. Units: J/kmol/K.
+ /*!
+ * This function is resolved here by calling the standard state
+ * thermo function.
+ *
+ * @param cpbar Output vector of species partial molar heat capacities
+ * Length = 1. units are J/kmol/K.
+ */
+ void getPartialMolarCp(doublereal* cpbar) const;
+
+
//! Get the species partial molar volumes. Units: m^3/kmol.
/*!
* This function is resolved here by calling the density function.
diff --git a/Cantera/src/thermo/WaterPropsIAPWS.cpp b/Cantera/src/thermo/WaterPropsIAPWS.cpp
index 92b34721e..af06a784a 100644
--- a/Cantera/src/thermo/WaterPropsIAPWS.cpp
+++ b/Cantera/src/thermo/WaterPropsIAPWS.cpp
@@ -103,22 +103,22 @@ double WaterPropsIAPWS::helmholtzFE() const{
*/
double WaterPropsIAPWS::pressure(double temperature, double rho) {
calcDim(temperature, rho);
- double retn = pressure_rhoRT();
- return (retn * rho * Rgas * temperature);
+ double retn = pressureM_rhoRT();
+ return (retn * rho * Rgas * temperature/M_water);
}
double WaterPropsIAPWS::pressure() const{
- double retn = pressure_rhoRT();
+ double retn = pressureM_rhoRT();
double rho = delta * Rho_c;
double temperature = T_c / tau;
- return (retn * rho * Rgas * temperature);
+ return (retn * rho * Rgas * temperature/M_water);
}
/*
* Calculates the pressure in dimensionless form
- * p/(rhoRT) at the currently stored tau and delta values
+ * pM/(rhoRT) at the currently stored tau and delta values
*/
-double WaterPropsIAPWS::pressure_rhoRT() const {
- double retn = m_phi->pressure_rhoRT(tau, delta);
+double WaterPropsIAPWS::pressureM_rhoRT() const {
+ double retn = m_phi->pressureM_rhoRT(tau, delta);
return retn;
}
@@ -284,7 +284,7 @@ double WaterPropsIAPWS::coeffThermExp(double temperature, double pressure) {
return retn;
}
-/**
+/*
* Returns the coefficient of isothermal compressibility
* of temperature and pressure.
* kappa = - d (ln V) / dP at constant T.
@@ -404,8 +404,7 @@ corr1(double temperature, double pressure, double &densLiq,
* p : Pascals : Newtons/m**2
*/
static int method = 1;
-double WaterPropsIAPWS::
-psat(double temperature) {
+double WaterPropsIAPWS::psat(double temperature) {
double densLiq = -1.0, densGas = -1.0, delGRT = 0.0;
double dp, pcorr;
double p = psat_est(temperature);
diff --git a/Cantera/src/thermo/WaterPropsIAPWS.h b/Cantera/src/thermo/WaterPropsIAPWS.h
index 66ba81a36..114a9d0b9 100644
--- a/Cantera/src/thermo/WaterPropsIAPWS.h
+++ b/Cantera/src/thermo/WaterPropsIAPWS.h
@@ -323,7 +323,7 @@ private:
double cp_R() const;
//! Return the current dimensionless pressure
- double pressure_rhoRT() const;
+ double pressureM_rhoRT() const;
protected:
diff --git a/Cantera/src/thermo/WaterPropsIAPWSphi.cpp b/Cantera/src/thermo/WaterPropsIAPWSphi.cpp
index e0d764e3f..128b6cf1e 100644
--- a/Cantera/src/thermo/WaterPropsIAPWSphi.cpp
+++ b/Cantera/src/thermo/WaterPropsIAPWSphi.cpp
@@ -631,7 +631,7 @@ double WaterPropsIAPWSphi::phi_d(double tau, double delta) {
*
* note: this is done so much, we have a seperate routine.
*/
-double WaterPropsIAPWSphi::pressure_rhoRT(double tau, double delta) {
+double WaterPropsIAPWSphi::pressureM_rhoRT(double tau, double delta) {
tdpolycalc(tau, delta);
double res = phiR_d();
double retn = 1.0 + delta * res;
diff --git a/Cantera/src/thermo/WaterPropsIAPWSphi.h b/Cantera/src/thermo/WaterPropsIAPWSphi.h
index a8e5653d5..a6f65d563 100644
--- a/Cantera/src/thermo/WaterPropsIAPWSphi.h
+++ b/Cantera/src/thermo/WaterPropsIAPWSphi.h
@@ -84,14 +84,14 @@ public:
//! Calculate the dimensionless pressure at tau and delta;
/*!
*
- * p/(rhoRT) = delta * phi_d() = 1.0 + delta phiR_d()
+ * pM/(rhoRT) = delta * phi_d() = 1.0 + delta phiR_d()
*
* @param tau Dimensionless temperature = T_c/T
* @param delta Dimensionless density = delta = rho / Rho_c
*
* note: this is done so much, we have a seperate routine.
*/
- double pressure_rhoRT(double tau, double delta);
+ double pressureM_rhoRT(double tau, double delta);
/**
* This program computes the reduced density, given the reduced pressure
diff --git a/Cantera/src/thermo/WaterSSTP.cpp b/Cantera/src/thermo/WaterSSTP.cpp
new file mode 100644
index 000000000..1d2b47388
--- /dev/null
+++ b/Cantera/src/thermo/WaterSSTP.cpp
@@ -0,0 +1,530 @@
+/**
+ * @file WaterTP.cpp
+ *
+ */
+/*
+ * Copywrite (2006) Sandia Corporation. Under the terms of
+ * Contract DE-AC04-94AL85000 with Sandia Corporation, the
+ * U.S. Government retains certain rights in this software.
+ */
+/*
+ * $Id$
+ */
+
+#include "xml.h"
+#include "WaterSSTP.h"
+#include "WaterPropsIAPWS.h"
+#include "importCTML.h"
+
+namespace Cantera {
+ /**
+ * Basic list of constructors and duplicators
+ */
+
+ WaterSSTP::WaterSSTP() :
+ SingleSpeciesTP(),
+ m_sub(0),
+ m_subflag(0),
+ m_mw(0.0),
+ EW_Offset(0.0),
+ SW_Offset(0.0),
+ m_verbose(0),
+ m_allowGasPhase(false)
+ {
+ constructPhase();
+ }
+
+
+ WaterSSTP::WaterSSTP(std::string inputFile, std::string id) :
+ SingleSpeciesTP(),
+ m_sub(0),
+ m_subflag(0),
+ m_mw(0.0),
+ EW_Offset(0.0),
+ SW_Offset(0.0),
+ m_verbose(0),
+ m_allowGasPhase(false)
+ {
+ constructPhaseFile(inputFile, id);
+ }
+
+
+ WaterSSTP::WaterSSTP(XML_Node& phaseRoot, std::string id) :
+ SingleSpeciesTP(),
+ m_sub(0),
+ m_subflag(0),
+ m_mw(0.0),
+ EW_Offset(0.0),
+ SW_Offset(0.0),
+ m_verbose(0),
+ m_allowGasPhase(false)
+ {
+ constructPhaseXML(phaseRoot, id) ;
+ }
+
+
+
+ WaterSSTP::WaterSSTP(const WaterSSTP &b) :
+ SingleSpeciesTP(b),
+ m_sub(0),
+ m_subflag(b.m_subflag),
+ m_mw(b.m_mw),
+ EW_Offset(b.EW_Offset),
+ SW_Offset(b.SW_Offset),
+ m_verbose(b.m_verbose),
+ m_allowGasPhase(b.m_allowGasPhase)
+ {
+ m_sub = new WaterPropsIAPWS(*(b.m_sub));
+ /*
+ * Use the assignment operator to do the brunt
+ * of the work for the copy construtor.
+ */
+ *this = b;
+ }
+
+ /*
+ * Assignment operator
+ */
+ WaterSSTP& WaterSSTP::operator=(const WaterSSTP&b) {
+ if (&b == this) return *this;
+ m_sub->operator=(*(b.m_sub));
+ m_subflag = b.m_subflag;
+ m_mw = b.m_mw;
+ m_verbose = b.m_verbose;
+ m_allowGasPhase = b.m_allowGasPhase;
+ return *this;
+ }
+
+
+ ThermoPhase *WaterSSTP::duplMyselfAsThermoPhase() {
+ WaterSSTP* wtp = new WaterSSTP(*this);
+ return (ThermoPhase *) wtp;
+ }
+
+ WaterSSTP::~WaterSSTP() {
+ delete m_sub;
+ }
+
+
+
+ void WaterSSTP::constructPhase() {
+ throw CanteraError("constructPhaseXML", "unimplemented");
+
+ }
+
+
+ /*
+ * @param infile XML file containing the description of the
+ * phase
+ *
+ * @param id Optional parameter identifying the name of the
+ * phase. If none is given, the first XML
+ * phase element will be used.
+ */
+ void WaterSSTP::constructPhaseXML(XML_Node& phaseNode, std::string id) {
+
+ /*
+ * Call the Cantera importPhase() function. This will import
+ * all of the species into the phase. This will also handle
+ * all of the solvent and solute standard states.
+ */
+ bool m_ok = importPhase(phaseNode, this);
+ if (!m_ok) {
+ throw CanteraError("initThermo","importPhase failed ");
+ }
+
+ }
+
+ /*
+ * constructPhaseFile
+ *
+ *
+ * This routine is a precursor to constructPhaseXML(XML_Node*)
+ * routine, which does most of the work.
+ *
+ * @param inputFile XML file containing the description of the
+ * phase
+ *
+ * @param id Optional parameter identifying the name of the
+ * phase. If none is given, the first XML
+ * phase element will be used.
+ */
+ void WaterSSTP::constructPhaseFile(std::string inputFile, std::string id) {
+
+ if (inputFile.size() == 0) {
+ throw CanteraError("WaterTp::initThermo",
+ "input file is null");
+ }
+ std::string path = findInputFile(inputFile);
+ std::ifstream fin(path.c_str());
+ if (!fin) {
+ throw CanteraError("WaterSSTP::initThermo","could not open "
+ +path+" for reading.");
+ }
+ /*
+ * The phase object automatically constructs an XML object.
+ * Use this object to store information.
+ */
+ XML_Node &phaseNode_XML = xml();
+ XML_Node *fxml = new XML_Node();
+ fxml->build(fin);
+ XML_Node *fxml_phase = findXMLPhase(fxml, id);
+ if (!fxml_phase) {
+ throw CanteraError("WaterSSTP::initThermo",
+ "ERROR: Can not find phase named " +
+ id + " in file named " + inputFile);
+ }
+ fxml_phase->copy(&phaseNode_XML);
+ constructPhaseXML(*fxml_phase, id);
+ delete fxml;
+ }
+
+
+
+ void WaterSSTP::initThermo() {
+ }
+
+ void WaterSSTP::
+ initThermoXML(XML_Node& phaseNode, std::string id) {
+ if (m_sub) delete m_sub;
+ m_sub = new WaterPropsIAPWS();
+ if (m_sub == 0) {
+ throw CanteraError("WaterSSTP::initThermo",
+ "could not create new substance object.");
+ }
+ /*
+ * Calculate the molecular weight. Note while there may
+ * be a very good calculated weight in the steam table
+ * class, using this weight may lead to codes exhibiting
+ * mass loss issues. We need to grab the elemental
+ * atomic weights used in the Element class and calculate
+ * a consistent H2O molecular weight based on that.
+ */
+ int nH = elementIndex("H");
+ if (nH < 0) {
+ throw CanteraError("WaterSSTP::initThermo",
+ "H not an element");
+ }
+ double mw_H = atomicWeight(nH);
+ int nO = elementIndex("O");
+ if (nO < 0) {
+ throw CanteraError("WaterSSTP::initThermo",
+ "O not an element");
+ }
+ double mw_O = atomicWeight(nO);
+ m_mw = 2.0 * mw_H + mw_O;
+ m_weight[0] = m_mw;
+ setMolecularWeight(0,m_mw);
+ double one = 1.0;
+ setMoleFractions(&one);
+
+ /*
+ * Set the baseline
+ */
+ doublereal T = 298.15;
+
+ doublereal presLow = 1.0E-2;
+ doublereal oneBar = 1.0E5;
+ doublereal dens = density();
+ doublereal dd = m_sub->density(T, presLow, WATER_GAS, dens);
+ setTemperature(T);
+ setDensity(dd);
+ SW_Offset = 0.0;
+ doublereal s = entropy_mole();
+ s -= GasConstant * log(oneBar/presLow);
+ if (s != 188.835E3) {
+ SW_Offset = 188.835E3 - s;
+ }
+ s = entropy_mole();
+ s -= GasConstant * log(oneBar/presLow);
+ printf("s = %g\n", s);
+
+ doublereal h = enthalpy_mole();
+ if (h != -241.826E6) {
+ EW_Offset = -241.826E6 - h;
+ }
+ h = enthalpy_mole();
+
+ printf("h = %g\n", h);
+
+
+ /*
+ * Set the initial state of the system to 298.15 K and
+ * 1 bar.
+ */
+ setTemperature(298.15);
+ double rho0 = m_sub->density(298.15, OneAtm, WATER_LIQUID);
+ setDensity(rho0);
+
+ /*
+ * We have to do something with the thermo function here.
+ */
+ if (m_spthermo) {
+ delete m_spthermo;
+ m_spthermo = 0;
+ }
+ }
+
+ void WaterSSTP::
+ setParametersFromXML(const XML_Node& eosdata) {
+ eosdata._require("model","PureFluid");
+ m_subflag = atoi(eosdata["fluid_type"].c_str());
+ if (m_subflag < 0)
+ throw CanteraError("WaterSSTP::setParametersFromXML",
+ "missing or negative substance flag");
+ }
+
+ /*
+ * Return the molar dimensionless enthalpy
+ */
+ void WaterSSTP::getEnthalpy_RT(doublereal* hrt) const {
+ double T = temperature();
+ double dens = density();
+ doublereal h = m_sub->enthalpy(T, dens);
+ *hrt = (h + EW_Offset)/(GasConstant*T);
+ }
+
+ /*
+ * Calculate the internal energy in mks units of
+ * J kmol-1
+ */
+ void WaterSSTP::getIntEnergy_RT(doublereal *ubar) const {
+ double T = temperature();
+ double dens = density();
+ doublereal u = m_sub->intEnergy(T, dens);
+ *ubar = (u + EW_Offset)/GasConstant;
+ }
+
+ /*
+ * Calculate the dimensionless entropy
+ */
+ void WaterSSTP::getEntropy_R(doublereal* sr) const {
+ double T = temperature();
+ double dens = density();
+ doublereal s = m_sub->entropy(T, dens);
+ sr[0] = (s + SW_Offset) / GasConstant;
+ }
+
+ /*
+ * Calculate the Gibbs free energy in mks units of
+ * J kmol-1 K-1.
+ */
+ void WaterSSTP::getGibbs_RT(doublereal *grt) const {
+ double T = temperature();
+ double dens = density();
+ doublereal g = m_sub->Gibbs(T, dens);
+ *grt = (g + EW_Offset - SW_Offset*T) / (GasConstant * T);
+ }
+
+ /*
+ * Calculate the Gibbs free energy in mks units of
+ * J kmol-1 K-1.
+ */
+ void WaterSSTP::getStandardChemPotentials(doublereal *gss) const {
+ double T = temperature();
+ double dens = density();
+ doublereal g = m_sub->Gibbs(T, dens);
+ *gss = (g + EW_Offset - SW_Offset*T);
+ }
+
+ void WaterSSTP::getCp_R(doublereal* cpr) const {
+ double T = temperature();
+ double dens = density();
+ doublereal cp = m_sub->cp(T, dens);
+ cpr[0] = cp / GasConstant;
+ }
+
+ /*
+ * Calculate the constant volume heat capacity
+ * in mks units of J kmol-1 K-1
+ */
+ doublereal WaterSSTP::
+ cv_mole() const {
+ double T = temperature();
+ double dens = density();
+ doublereal cv = m_sub->cv(T, dens);
+ return cv;
+ }
+
+ // @name Thermodynamic Values for the Species Reference State
+
+
+ void WaterSSTP::getEnthalpy_RT_ref(doublereal *hrt) const {
+ doublereal p = pressure();
+ double T = temperature();
+ double dens = density();
+ int waterState = WATER_GAS;
+ double rc = m_sub->Rhocrit();
+ if (dens > rc) {
+ waterState = WATER_LIQUID;
+ }
+ doublereal dd = m_sub->density(T, OneAtm, waterState, dens);
+ if (dd <= 0.0) {
+ throw CanteraError("setPressure", "error");
+ }
+ doublereal h = m_sub->enthalpy(T, dd);
+ *hrt = (h + EW_Offset) / (GasConstant * T);
+ dd = m_sub->density(T, p, waterState, dens);
+ }
+
+ void WaterSSTP::getGibbs_RT_ref(doublereal *grt) const {
+ doublereal p = pressure();
+ double T = temperature();
+ double dens = density();
+ int waterState = WATER_GAS;
+ double rc = m_sub->Rhocrit();
+ if (dens > rc) {
+ waterState = WATER_LIQUID;
+ }
+ doublereal dd = m_sub->density(T, OneAtm, waterState, dens);
+ if (dd <= 0.0) {
+ throw CanteraError("setPressure", "error");
+ }
+ m_sub->setState(T, dd);
+ doublereal g = m_sub->Gibbs(T, dd);
+ *grt = (g + EW_Offset - SW_Offset*T)/ (GasConstant * T);
+ dd = m_sub->density(T, p, waterState, dens);
+
+ }
+
+ void WaterSSTP::getGibbs_ref(doublereal *g) const {
+ getGibbs_RT_ref(g);
+ doublereal rt = _RT();
+ for (int k = 0; k < m_kk; k++) {
+ g[k] *= rt;
+ }
+ }
+
+ void WaterSSTP::getEntropy_R_ref(doublereal *sr) const {
+ doublereal p = pressure();
+ double T = temperature();
+ double dens = density();
+ int waterState = WATER_GAS;
+ double rc = m_sub->Rhocrit();
+ if (dens > rc) {
+ waterState = WATER_LIQUID;
+ }
+ doublereal dd = m_sub->density(T, OneAtm, waterState, dens);
+
+ if (dd <= 0.0) {
+ throw CanteraError("setPressure", "error");
+ }
+ m_sub->setState(T, dd);
+
+ doublereal s = m_sub->entropy(T, dd);
+ *sr = (s + SW_Offset)/ (GasConstant);
+ dd = m_sub->density(T, p, waterState, dens);
+
+ }
+
+ void WaterSSTP::getCp_R_ref(doublereal *cpr) const {
+ doublereal p = pressure();
+ double T = temperature();
+ double dens = density();
+ int waterState = WATER_GAS;
+ double rc = m_sub->Rhocrit();
+ if (dens > rc) {
+ waterState = WATER_LIQUID;
+ }
+ doublereal dd = m_sub->density(T, OneAtm, waterState, dens);
+ m_sub->setState(T, dd);
+ if (dd <= 0.0) {
+ throw CanteraError("setPressure", "error");
+ }
+ doublereal cp = m_sub->cp(T, dd);
+ *cpr = cp / (GasConstant);
+ dd = m_sub->density(T, p, waterState, dens);
+ }
+
+ void WaterSSTP::getStandardVolumes_ref(doublereal *vol) const {
+ doublereal p = pressure();
+ double T = temperature();
+ double dens = density();
+ int waterState = WATER_GAS;
+ double rc = m_sub->Rhocrit();
+ if (dens > rc) {
+ waterState = WATER_LIQUID;
+ }
+ doublereal dd = m_sub->density(T, OneAtm, waterState, dens);
+ if (dd <= 0.0) {
+ throw CanteraError("setPressure", "error");
+ }
+ *vol = meanMolecularWeight() /dd;
+ dd = m_sub->density(T, p, waterState, dens);
+ }
+
+ /*
+ * Calculate the pressure (Pascals), given the temperature and density
+ * Temperature: kelvin
+ * rho: density in kg m-3
+ */
+ doublereal WaterSSTP::
+ pressure() const {
+ double T = temperature();
+ double dens = density();
+ doublereal p = m_sub->pressure(T, dens);
+ return p;
+ }
+
+ void WaterSSTP::
+ setPressure(doublereal p) {
+ double T = temperature();
+ double dens = density();
+ int waterState = WATER_GAS;
+ double rc = m_sub->Rhocrit();
+ if (dens > rc) {
+ waterState = WATER_LIQUID;
+ }
+ doublereal dd = m_sub->density(T, p, waterState, dens);
+ if (dd <= 0.0) {
+ throw CanteraError("setPressure", "error");
+ }
+ setDensity(dd);
+ }
+
+
+ // critical temperature
+ doublereal WaterSSTP::critTemperature() const { return m_sub->Tcrit(); }
+
+ // critical pressure
+ doublereal WaterSSTP::critPressure() const { return m_sub->Pcrit(); }
+
+ // critical density
+ doublereal WaterSSTP::critDensity() const { return m_sub->Rhocrit(); }
+
+
+
+ void WaterSSTP::setTemperature(double temp) {
+ State::setTemperature(temp);
+ doublereal dd = density();
+ m_sub->setState(temp, dd);
+ }
+
+
+
+ // saturation pressure
+ doublereal WaterSSTP::satPressure(doublereal t){
+ doublereal pp = m_sub->psat(t);
+ double dens = density();
+ setTemperature(t);
+ setDensity(dens);
+ return pp;
+ }
+
+ // Return the fraction of vapor at the current conditions
+ doublereal WaterSSTP::vaporFraction() const {
+ if (temperature() >= m_sub->Tcrit()) {
+ double dens = density();
+ if (dens >= m_sub->Rhocrit()) {
+ return 0.0;
+ }
+ return 1.0;
+ }
+ /*
+ * If below tcrit we always return 0 from this class
+ */
+ return 0.0;
+ }
+
+
+}
diff --git a/Cantera/src/thermo/WaterTP.h b/Cantera/src/thermo/WaterSSTP.h
similarity index 58%
rename from Cantera/src/thermo/WaterTP.h
rename to Cantera/src/thermo/WaterSSTP.h
index 88e313419..13bc09c0a 100644
--- a/Cantera/src/thermo/WaterTP.h
+++ b/Cantera/src/thermo/WaterSSTP.h
@@ -1,7 +1,6 @@
/**
- * @file WaterTP.h
- *
- * Declares a ThermoPhase class consisting of
+ * @file WaterSSTP.h
+ * Declares a %ThermoPhase class consisting of
* pure water.
*/
/*
@@ -13,10 +12,10 @@
* $Id$
*/
-#ifndef CT_WATERTP_H
-#define CT_WATERTP_H
+#ifndef CT_WATERSSTP_H
+#define CT_WATERSSTP_H
-#include "ThermoPhase.h"
+#include "SingleSpeciesTP.h"
class WaterPropsIAPWS;
@@ -52,41 +51,38 @@ namespace Cantera {
* They assume u_liq(TP) = 0.0, s_liq(TP) = 0.0, where TP is the
* triple point conditions.
*
- * @todo
- * I should have made this inherit from SingleSpeciesTP!
- *
* @ingroup thermoprops
*
*/
- class WaterTP : public ThermoPhase {
+ class WaterSSTP : public SingleSpeciesTP {
public:
//! Base constructor
- WaterTP();
+ WaterSSTP();
//! Copy constructor
- WaterTP(const WaterTP &);
+ WaterSSTP(const WaterSSTP &);
//! Assignment operator
- WaterTP& operator=(const WaterTP&);
+ WaterSSTP& operator=(const WaterSSTP&);
//! Full constructor for a water phase
/*!
* @param inputFile String name of the input file
* @param id string id of the phase name
*/
- WaterTP(std::string inputFile, std::string id = "");
+ WaterSSTP(std::string inputFile, std::string id = "");
//! Full constructor for a water phase
/*!
* @param phaseRef XML node referencing the water phase.
* @param id string id of the phase name
*/
- WaterTP(XML_Node& phaseRef, std::string id = "");
+ WaterSSTP(XML_Node& phaseRef, std::string id = "");
//! Destructor
- virtual ~WaterTP();
+ virtual ~WaterSSTP();
//! Duplicator from a ThermoPhase object
ThermoPhase *duplMyselfAsThermoPhase();
@@ -103,11 +99,7 @@ namespace Cantera {
* @name Molar Thermodynamic Properties of the Solution --------------
* @{
*/
- virtual doublereal enthalpy_mole() const;
- virtual doublereal intEnergy_mole() const;
- virtual doublereal entropy_mole() const;
- virtual doublereal gibbs_mole() const;
- virtual doublereal cp_mole() const;
+
virtual doublereal cv_mole() const;
//@}
@@ -134,19 +126,158 @@ namespace Cantera {
//@{
- //! get the chemical potential of the water
- /*!
- * @param mu vector of chemical potentials.
- */
- virtual void getChemPotentials(doublereal* mu) const {
- mu[0] = gibbs_mole();
- }
-
//@}
/// @name Properties of the Standard State of the Species
// in the Solution --
//@{
+
+ //!Get the gibbs function for the species
+ //! standard states at the current T and P of the solution.
+ /*!
+ * @param grt Vector of length m_kk, which on return sr[k]
+ * will contain the
+ * standard state gibbs function for species k.
+ */
+ virtual void getStandardChemPotentials(doublereal* gss) const;
+
+ //!Get the nondimensional gibbs function for the species
+ //! standard states at the current T and P of the solution.
+ /*!
+ * @param grt Vector of length m_kk, which on return sr[k]
+ * will contain the nondimensional
+ * standard state gibbs function for species k.
+ */
+ virtual void getGibbs_RT(doublereal* grt) const;
+
+ //! Get the array of nondimensional Enthalpy functions for the standard state species
+ //! at the current T and P of the solution.
+ /*!
+ *
+ * @param hrt Vector of length m_kk, which on return hrt[k]
+ * will contain the nondimensional
+ * standard state enthalpy of species k.
+ */
+ void getEnthalpy_RT(doublereal* hrt) const;
+
+
+ //! Get the nondimensional Entropies for the species
+ //! standard states at the current T and P of the solution.
+ /*!
+ * @param sr Vector of length m_kk, which on return sr[k]
+ * will contain the nondimensional
+ * standard state entropy for species k.
+ */
+ void getEntropy_R(doublereal* sr) const;
+
+ //! Get the nondimensional heat capacity at constant pressure
+ //! function for the species standard states at the current T and P of the solution.
+ /*!
+ *
+ * @param cpr Vector of length m_kk, which on return cpr[k]
+ * will contain the nondimensional
+ * constant pressure heat capacity for species k.
+ */
+ virtual void getCp_R(doublereal* cpr) const;
+
+ //! Returns the vector of nondimensional
+ //! internal Energies of the standard state at the current
+ //! temperature and pressure of the solution for each species.
+ /*!
+ *
+ * @param urt Output vector of standard state nondimensional internal energies.
+ * Length: m_kk.
+ */
+ virtual void getIntEnergy_RT(doublereal *urt) const;
+
+ //@}
+ //! @name Thermodynamic Values for the Species Reference State
+ /*!
+ * All functions in this group need to be overrided, because
+ * the m_spthermo SpeciesThermo function is not adequate for
+ * the real equation of state.
+ *
+ */
+ //@{
+
+
+
+ //! Returns the vector of nondimensional
+ //! enthalpies of the reference state at the current temperature
+ //! of the solution and the reference pressure for the species.
+ /*!
+ * @param hrt Output vector containing the nondimensional reference state enthalpies
+ * Length: m_kk.
+ */
+ virtual void getEnthalpy_RT_ref(doublereal *hrt) const;
+
+ /*!
+ * Returns the vector of nondimensional
+ * enthalpies of the reference state at the current temperature
+ * of the solution and the reference pressure for the species.
+ *
+ * This function is resolved in this class. It is assumed that the m_spthermo species thermo
+ * pointer is populated and yields the reference state.
+ *
+ * @param grt Output vector containing the nondimensional reference state
+ * Gibbs Free energies. Length: m_kk.
+ */
+ virtual void getGibbs_RT_ref(doublereal *grt) const;
+
+
+ /*!
+ * Returns the vector of the
+ * gibbs function of the reference state at the current temperature
+ * of the solution and the reference pressure for the species.
+ * units = J/kmol
+ *
+ * This function is resolved in this class. It is assumed that the m_spthermo
+ * species thermo
+ * pointer is populated and yields the reference state.
+ *
+ * @param g Output vector containing the reference state
+ * Gibbs Free energies. Length: m_kk. Units: J/kmol.
+ */
+ virtual void getGibbs_ref(doublereal *g) const;
+
+ /*!
+ * Returns the vector of nondimensional
+ * entropies of the reference state at the current temperature
+ * of the solution and the reference pressure for each species.
+ *
+ * This function is resolved in this class. It is assumed that the m_spthermo species thermo
+ * pointer is populated and yields the reference state.
+ *
+ * @param er Output vector containing the nondimensional reference state
+ * entropies. Length: m_kk.
+ */
+ virtual void getEntropy_R_ref(doublereal *er) const;
+
+ /*!
+ * Returns the vector of nondimensional
+ * constant pressure heat capacities of the reference state
+ * at the current temperature of the solution
+ * and reference pressure for each species.
+ *
+ * This function is resolved in this class. It is assumed that the m_spthermo
+ * species thermo
+ * pointer is populated and yields the reference state.
+ *
+ * @param cprt Output vector of nondimensional reference state
+ * heat capacities at constant pressure for the species.
+ * Length: m_kk
+ */
+ virtual void getCp_R_ref(doublereal *cprt) const;
+
+ //! Get the molar volumes of the species reference states at the current
+ //! T and P_ref of the solution.
+ /*!
+ * units = m^3 / kmol
+ *
+ * @param vol Output vector containing the standard state volumes.
+ * Length: m_kk.
+ */
+ virtual void getStandardVolumes_ref(doublereal *vol) const;
/// critical temperature
virtual doublereal critTemperature() const;
@@ -282,6 +413,14 @@ namespace Cantera {
void check(doublereal v = 0.0) const;
void reportTPXError() const;
+ protected:
+ /**
+ * @internal
+ * This internal routine must be overwritten because
+ * it is not applicable.
+ */
+ void _updateThermo() const;
+
private:
mutable WaterPropsIAPWS *m_sub;
int m_subflag;
diff --git a/Cantera/src/thermo/WaterTP.cpp b/Cantera/src/thermo/WaterTP.cpp
deleted file mode 100644
index 1bc0d23d3..000000000
--- a/Cantera/src/thermo/WaterTP.cpp
+++ /dev/null
@@ -1,422 +0,0 @@
-/**
- * @file WaterTP.cpp
- *
- */
-/*
- * Copywrite (2006) Sandia Corporation. Under the terms of
- * Contract DE-AC04-94AL85000 with Sandia Corporation, the
- * U.S. Government retains certain rights in this software.
- */
-/*
- * $Id$
- */
-
-#include "xml.h"
-#include "WaterTP.h"
-#include "WaterPropsIAPWS.h"
-#include "importCTML.h"
-
-namespace Cantera {
- /**
- * Basic list of constructors and duplicators
- */
-
- WaterTP::WaterTP() :
- ThermoPhase(),
- m_sub(0),
- m_subflag(0),
- m_mw(0.0),
- EW_Offset(0.0),
- SW_Offset(0.0),
- m_verbose(0),
- m_allowGasPhase(false)
- {
- constructPhase();
- }
-
-
- WaterTP::WaterTP(std::string inputFile, std::string id) :
- ThermoPhase(),
- m_sub(0),
- m_subflag(0),
- m_mw(0.0),
- EW_Offset(0.0),
- SW_Offset(0.0),
- m_verbose(0),
- m_allowGasPhase(false)
- {
- constructPhaseFile(inputFile, id);
- }
-
-
- WaterTP::WaterTP(XML_Node& phaseRoot, std::string id) :
- ThermoPhase(),
- m_sub(0),
- m_subflag(0),
- m_mw(0.0),
- EW_Offset(0.0),
- SW_Offset(0.0),
- m_verbose(0),
- m_allowGasPhase(false)
- {
- constructPhaseXML(phaseRoot, id) ;
- }
-
-
-
- WaterTP::WaterTP(const WaterTP &b) :
- ThermoPhase(b),
- m_sub(0),
- m_subflag(b.m_subflag),
- m_mw(b.m_mw),
- EW_Offset(b.EW_Offset),
- SW_Offset(b.SW_Offset),
- m_verbose(b.m_verbose),
- m_allowGasPhase(b.m_allowGasPhase)
- {
- m_sub = new WaterPropsIAPWS(*(b.m_sub));
- /*
- * Use the assignment operator to do the brunt
- * of the work for the copy construtor.
- */
- *this = b;
- }
-
- /*
- * Assignment operator
- */
- WaterTP& WaterTP::operator=(const WaterTP&b) {
- if (&b == this) return *this;
- m_sub->operator=(*(b.m_sub));
- m_subflag = b.m_subflag;
- m_mw = b.m_mw;
- m_verbose = b.m_verbose;
- m_allowGasPhase = b.m_allowGasPhase;
- return *this;
- }
-
-
- ThermoPhase *WaterTP::duplMyselfAsThermoPhase() {
- WaterTP* wtp = new WaterTP(*this);
- return (ThermoPhase *) wtp;
- }
-
- WaterTP::~WaterTP() {
- delete m_sub;
- }
-
-
-
- void WaterTP::constructPhase() {
- throw CanteraError("constructPhaseXML", "unimplemented");
-
- }
-
-
- /*
- * @param infile XML file containing the description of the
- * phase
- *
- * @param id Optional parameter identifying the name of the
- * phase. If none is given, the first XML
- * phase element will be used.
- */
- void WaterTP::constructPhaseXML(XML_Node& phaseNode, std::string id) {
-
- /*
- * Call the Cantera importPhase() function. This will import
- * all of the species into the phase. This will also handle
- * all of the solvent and solute standard states.
- */
- bool m_ok = importPhase(phaseNode, this);
- if (!m_ok) {
- throw CanteraError("initThermo","importPhase failed ");
- }
-
- }
-
- /*
- * initThermo():
- *
- *
- * This routine is a precursor to constructPhaseXML(XML_Node*)
- * routine, which does most of the work.
- *
- * @param inputFile XML file containing the description of the
- * phase
- *
- * @param id Optional parameter identifying the name of the
- * phase. If none is given, the first XML
- * phase element will be used.
- */
- void WaterTP::constructPhaseFile(std::string inputFile, std::string id) {
-
- if (inputFile.size() == 0) {
- throw CanteraError("WaterTp::initThermo",
- "input file is null");
- }
- std::string path = findInputFile(inputFile);
- std::ifstream fin(path.c_str());
- if (!fin) {
- throw CanteraError("WaterTP::initThermo","could not open "
- +path+" for reading.");
- }
- /*
- * The phase object automatically constructs an XML object.
- * Use this object to store information.
- */
- XML_Node &phaseNode_XML = xml();
- XML_Node *fxml = new XML_Node();
- fxml->build(fin);
- XML_Node *fxml_phase = findXMLPhase(fxml, id);
- if (!fxml_phase) {
- throw CanteraError("WaterTP::initThermo",
- "ERROR: Can not find phase named " +
- id + " in file named " + inputFile);
- }
- fxml_phase->copy(&phaseNode_XML);
- constructPhaseXML(*fxml_phase, id);
- delete fxml;
- }
-
-
-
- void WaterTP::initThermo() {
- }
-
- void WaterTP::
- initThermoXML(XML_Node& phaseNode, std::string id) {
- if (m_sub) delete m_sub;
- m_sub = new WaterPropsIAPWS();
- if (m_sub == 0) {
- throw CanteraError("WaterTP::initThermo",
- "could not create new substance object.");
- }
- /*
- * Calculate the molecular weight. Note while there may
- * be a very good calculated weight in the steam table
- * class, using this weight may lead to codes exhibiting
- * mass loss issues. We need to grab the elemental
- * atomic weights used in the Element class and calculate
- * a consistent H2O molecular weight based on that.
- */
- int nH = elementIndex("H");
- if (nH < 0) {
- throw CanteraError("WaterTP::initThermo",
- "H not an element");
- }
- double mw_H = atomicWeight(nH);
- int nO = elementIndex("O");
- if (nO < 0) {
- throw CanteraError("WaterTP::initThermo",
- "O not an element");
- }
- double mw_O = atomicWeight(nO);
- m_mw = 2.0 * mw_H + mw_O;
- m_weight[0] = m_mw;
- setMolecularWeight(0,m_mw);
- double one = 1.0;
- setMoleFractions(&one);
-
- /*
- * Set the baseline
- */
- doublereal T = 298.15;
-
- doublereal presLow = 1.0E-2;
- doublereal oneBar = 1.0E5;
- doublereal dens = density();
- doublereal dd = m_sub->density(T, presLow, WATER_GAS, dens);
- setTemperature(T);
- setDensity(dd);
- SW_Offset = 0.0;
- doublereal s = entropy_mole();
- s -= GasConstant * log(oneBar/presLow);
- if (s != 188.835E3) {
- SW_Offset = 188.835E3 - s;
- }
- s = entropy_mole();
- s -= GasConstant * log(oneBar/presLow);
- printf("s = %g\n", s);
-
- doublereal h = enthalpy_mole();
- if (h != -241.826E6) {
- EW_Offset = -241.826E6 - h;
- }
- h = enthalpy_mole();
-
- printf("h = %g\n", h);
-
-
- /*
- * Set the initial state of the system to 298.15 K and
- * 1 bar.
- */
- setTemperature(298.15);
- double rho0 = m_sub->density(298.15, OneAtm, WATER_LIQUID);
- setDensity(rho0);
-
- /*
- * We have to do something with the thermo function here.
- */
- if (m_spthermo) {
- delete m_spthermo;
- m_spthermo = 0;
- }
- }
-
- void WaterTP::
- setParametersFromXML(const XML_Node& eosdata) {
- eosdata._require("model","PureFluid");
- m_subflag = atoi(eosdata["fluid_type"].c_str());
- if (m_subflag < 0)
- throw CanteraError("WaterTP::setParametersFromXML",
- "missing or negative substance flag");
- }
-
- /*
- * Return the molar enthalpy in units of J kmol-1
- */
- doublereal WaterTP::
- enthalpy_mole() const {
- double T = temperature();
- double dens = density();
- doublereal h = m_sub->enthalpy(T, dens);
- return (h + EW_Offset);
- }
-
- /*
- * Calculate the internal energy in mks units of
- * J kmol-1
- */
- doublereal WaterTP::
- intEnergy_mole() const {
- double T = temperature();
- double dens = density();
- doublereal u = m_sub->intEnergy(T, dens);
- return (u + EW_Offset);
- }
-
- /*
- * Calculate the entropy in mks units of
- * J kmol-1 K-1
- */
- doublereal WaterTP::
- entropy_mole() const {
- double T = temperature();
- double dens = density();
- doublereal s = m_sub->entropy(T, dens);
- return (s + SW_Offset);
- }
-
- /*
- * Calculate the Gibbs free energy in mks units of
- * J kmol-1 K-1.
- */
- doublereal WaterTP::
- gibbs_mole() const {
- double T = temperature();
- double dens = density();
- doublereal g = m_sub->Gibbs(T, dens);
- return (g + EW_Offset - SW_Offset*T);
- }
-
- /*
- * Calculate the constant pressure heat capacity
- * in mks units of J kmol-1 K-1
- */
- doublereal WaterTP::
- cp_mole() const {
- double T = temperature();
- double dens = density();
- doublereal cp = m_sub->cp(T, dens);
- return cp;
- }
-
- /*
- * Calculate the constant volume heat capacity
- * in mks units of J kmol-1 K-1
- */
- doublereal WaterTP::
- cv_mole() const {
- double T = temperature();
- double dens = density();
- doublereal cv = m_sub->cv(T, dens);
- return cv;
- }
-
- /*
- * Calculate the pressure (Pascals), given the temperature and density
- * Temperature: kelvin
- * rho: density in kg m-3
- */
- doublereal WaterTP::
- pressure() const {
- double T = temperature();
- double dens = density();
- doublereal p = m_sub->pressure(T, dens);
- return p;
- }
-
- void WaterTP::
- setPressure(doublereal p) {
- double T = temperature();
- double dens = density();
- int waterState = WATER_GAS;
- double rc = m_sub->Rhocrit();
- if (dens > rc) {
- waterState = WATER_LIQUID;
- }
- doublereal dd = m_sub->density(T, p, waterState, dens);
- if (dd <= 0.0) {
- throw CanteraError("setPressure", "error");
- }
- setDensity(dd);
- }
-
-
- // critical temperature
- doublereal WaterTP::critTemperature() const { return m_sub->Tcrit(); }
-
- // critical pressure
- doublereal WaterTP::critPressure() const { return m_sub->Pcrit(); }
-
- // critical density
- doublereal WaterTP::critDensity() const { return m_sub->Rhocrit(); }
-
-
-
- void WaterTP::setTemperature(double temp) {
- State::setTemperature(temp);
- doublereal dd = density();
- m_sub->setState(temp, dd);
- }
-
-
-
- // saturation pressure
- doublereal WaterTP::satPressure(doublereal t){
- doublereal pp = m_sub->psat(t);
- double dens = density();
- setTemperature(t);
- setDensity(dens);
- return pp;
- }
-
- // Return the fraction of vapor at the current conditions
- doublereal WaterTP::vaporFraction() const {
- if (temperature() >= m_sub->Tcrit()) {
- double dens = density();
- if (dens >= m_sub->Rhocrit()) {
- return 0.0;
- }
- return 1.0;
- }
- /*
- * If below tcrit we always return 0 from this class
- */
- return 0.0;
- }
-
-
-}
diff --git a/ext/f2c_libs/arith.h b/ext/f2c_libs/arith.h
index 9746f9c1c..76539f82b 100644
--- a/ext/f2c_libs/arith.h
+++ b/ext/f2c_libs/arith.h
@@ -1,3 +1,2 @@
-#define IEEE_8087
-#define Arith_Kind_ASL 1
-#define Double_Align
+#define IEEE_8087
+#define Arith_Kind_ASL 1
diff --git a/test_problems/cathermo/testIAPWSPres/output_blessed.txt b/test_problems/cathermo/testIAPWSPres/output_blessed.txt
index 69bfd8773..d60ad14e3 100644
--- a/test_problems/cathermo/testIAPWSPres/output_blessed.txt
+++ b/test_problems/cathermo/testIAPWSPres/output_blessed.txt
@@ -1,4 +1,4 @@
-pres = 182080
+pres = 10107
psat(273.16) = 611.655
dens (liquid) = 999.793 kg m-3
intEng (liquid) ~= 0.0 J/kmol (less than fabs(5.0E-7))
diff --git a/test_problems/cathermo/testIAPWSPres/runtest b/test_problems/cathermo/testIAPWSPres/runtest
index abc65a02b..a9b880b30 100755
--- a/test_problems/cathermo/testIAPWSPres/runtest
+++ b/test_problems/cathermo/testIAPWSPres/runtest
@@ -4,7 +4,7 @@
temp_success="1"
/bin/rm -f output.txt outputa.txt
-testName=testPress
+testName=testIAPWSPress
#################################################################
#
#################################################################
diff --git a/test_problems/cathermo/testWaterTP/.cvsignore b/test_problems/cathermo/testWaterTP/.cvsignore
index 0b0e7ecad..19060c353 100644
--- a/test_problems/cathermo/testWaterTP/.cvsignore
+++ b/test_problems/cathermo/testWaterTP/.cvsignore
@@ -8,4 +8,4 @@ csvCode.txt
diff_test.out
test.diff
test.out
-testWaterTP
+testWaterSSTP
diff --git a/test_problems/cathermo/testWaterTP/Makefile.in b/test_problems/cathermo/testWaterTP/Makefile.in
index bb3ddb68a..fd7e991ed 100644
--- a/test_problems/cathermo/testWaterTP/Makefile.in
+++ b/test_problems/cathermo/testWaterTP/Makefile.in
@@ -11,11 +11,11 @@
.SUFFIXES : .d
# the name of the executable program to be created
-PROG_NAME = testWaterTP
+PROG_NAME = testWaterSSTP
# the object files to be linked together. List those generated from Fortran
# and from C/C++ separately
-OBJS = testWaterTP.o
+OBJS = testWaterSSTP.o
# Location of the current build. Will assume that tests are run
# in the source directory tree location
diff --git a/test_problems/cathermo/testWaterTP/output_blessed.txt b/test_problems/cathermo/testWaterTP/output_blessed.txt
index e00e9c677..507913150 100644
--- a/test_problems/cathermo/testWaterTP/output_blessed.txt
+++ b/test_problems/cathermo/testWaterTP/output_blessed.txt
@@ -38,6 +38,46 @@ Liquid Densities:
400 2.46261 2.45769 937.486 0.0192166
500 26.4447 26.392 831.318 0.0216707
+Liquid 1bar or psat State: Partial Molar Quantities
+ T press psat Cpbar Sbar -(G0-H298)/T H0-H298 Volume
+ (K) (bar) (bar) (J/molK) (J/molK) (J/molK) (kJ/mol) m3/kmol
+ 273.19 1 0.00612989 76.0121 63.3157 70.2194 -1.88603 0.0180181
+ 298.15 1 0.0316993 75.3276 69.9224 69.9224 0 0.0180686
+ 300 1 0.0353681 75.3153 70.3884 69.9239 0.139344 0.0180775
+ 373.15 1.01621 1.01418 75.9465 86.857 71.6855 5.66125 0.0187982
+ 400 2.46261 2.45769 76.6642 92.1544 72.8766 7.71115 0.0192166
+ 500 26.4447 26.392 84.0126 109.805 78.4402 15.6825 0.0216707
+
+Liquid 1bar or psat State: Standard State Quantities
+ T press psat Cpbar Sbar -(G0-H298)/T H0-H298 Volume
+ (K) (bar) (bar) (J/molK) (J/molK) (J/molK) (kJ/mol) m3/kmol
+ 273.19 1 0.00612989 76.0121 63.3157 70.2194 -1.88603 0.0180181
+ 298.15 1 0.0316993 75.3276 69.9224 69.9224 0 0.0180686
+ 300 1 0.0353681 75.3153 70.3884 69.9239 0.139344 0.0180775
+ 373.15 1.01621 1.01418 75.9465 86.857 71.6855 5.66125 0.0187982
+ 400 2.46261 2.45769 76.6642 92.1544 72.8766 7.71115 0.0192166
+ 500 26.4447 26.392 84.0126 109.805 78.4402 15.6825 0.0216707
+
+Liquid 1bar or psat State: Reference State Quantities (Always 1 atm no matter what system pressure is)
+ T press psat Cpbar Sbar -(G0-H298)/T H0-H298 Volume
+ (K) (bar) (bar) (J/molK) (J/molK) (J/molK) (kJ/mol) m3/kmol
+ 273.19 1 0.00612989 76.0119 63.3157 70.2193 -1.886 0.0180181
+ 298.15 1 0.0316993 75.3275 69.9224 69.9224 2.21044e-05 0.0180686
+ 300 1 0.0353681 75.3153 70.3883 69.9238 0.139366 0.0180775
+ 373.15 1.01621 1.01418 75.9465 86.857 71.6855 5.66125 0.0187982
+ 400 2.46261 2.45769 76.6712 92.1569 72.8835 7.70936 0.0192181
+ 500 26.4447 26.392 84.4316 109.897 78.5506 15.6732 0.021734
+
+Liquid 1 atm: Standard State Quantities - Should agree with table above
+ T press psat Cpbar Sbar -(G0-H298)/T H0-H298 Volume
+ (K) (bar) (bar) (J/molK) (J/molK) (J/molK) (kJ/mol) m3/kmol
+ 273.19 1.01325 0.00612989 76.0119 63.3157 70.2193 -1.886 0.0180181
+ 298.15 1.01325 0.0316993 75.3275 69.9224 69.9224 2.21044e-05 0.0180686
+ 300 1.01325 0.0353681 75.3153 70.3883 69.9238 0.139366 0.0180775
+ 373.15 1.01325 1.01418 75.9465 86.857 71.6855 5.66125 0.0187982
+ 400 1.01325 2.45769 76.6712 92.1569 72.8835 7.70936 0.0192181
+ 500 1.01325 26.392 84.4316 109.897 78.5506 15.6732 0.021734
+
Table of increasing Enthalpy at 1 atm
diff --git a/test_problems/cathermo/testWaterTP/runtest b/test_problems/cathermo/testWaterTP/runtest
index ef6267d07..80a62b12e 100755
--- a/test_problems/cathermo/testWaterTP/runtest
+++ b/test_problems/cathermo/testWaterTP/runtest
@@ -11,7 +11,7 @@ testName=testWaterTP
CANTERA_DATA=${CANTERA_DATA:=../../../data/inputs}; export CANTERA_DATA
CANTERA_BIN=${CANTERA_BIN:=../../../bin}
-./testWaterTP > output.txt
+./testWaterSSTP > output.txt
retnStat=$?
if [ $retnStat != "0" ]
then
diff --git a/test_problems/cathermo/testWaterTP/testWaterTP.cpp b/test_problems/cathermo/testWaterTP/testWaterSSTP.cpp
similarity index 56%
rename from test_problems/cathermo/testWaterTP/testWaterTP.cpp
rename to test_problems/cathermo/testWaterTP/testWaterSSTP.cpp
index 07e67ea30..cec3e0fc6 100644
--- a/test_problems/cathermo/testWaterTP/testWaterTP.cpp
+++ b/test_problems/cathermo/testWaterTP/testWaterSSTP.cpp
@@ -3,7 +3,7 @@
*/
#include "stdio.h"
#include "math.h"
-#include "WaterTP.h"
+#include "WaterSSTP.h"
#include
using namespace std;
using namespace Cantera;
@@ -22,7 +22,7 @@ int main () {
double pres;
try {
- WaterTP *w = new WaterTP("waterTPphase.xml","");
+ WaterSSTP *w = new WaterSSTP("waterTPphase.xml","");
/*
@@ -34,6 +34,7 @@ int main () {
double presLow = 1.0E-2;
temp = 298.15;
double oneBar = 1.0E5;
+ double vol;
printf("Comparisons to NIST: (see http://webbook.nist.gov):\n\n");
@@ -55,7 +56,7 @@ int main () {
T[3] = 1000.;
double Cp0, delh0, delg0, g;
-
+ double Cp0_ss;
printf("\nIdeal Gas Standard State:\n");
printf (" T Cp0 S0 "
" -(G0-H298)/T H0-H298\n");
@@ -69,6 +70,14 @@ int main () {
g = w->gibbs_mole();
delg0 = (g - h298)/temp + GasConstant * log(oneBar/presLow);
Cp0 = w->cp_mole();
+ {
+ w->getCp_R(&Cp0_ss);
+ Cp0_ss *= GasConstant;
+ if (fabs(Cp0_ss - Cp0) > 1.0E-5) {
+ printf("Inconsistency!\n");
+ exit(-1);
+ }
+ }
s = w->entropy_mole();
s -= GasConstant * log(oneBar/presLow);
printf("%10g %10g %13g %13g %13g\n", temp, Cp0*1.0E-3, s*1.0E-3,
@@ -142,6 +151,120 @@ int main () {
}
+
+ printf("\nLiquid 1bar or psat State: Partial Molar Quantities\n");
+ printf (" T press psat Cpbar Sbar "
+ " -(G0-H298)/T H0-H298 Volume\n");
+ printf (" (K) (bar) (bar) (J/molK) (J/molK)"
+ " (J/molK) (kJ/mol) m3/kmol\n");
+
+ for (int i = 0; i < 6; i++) {
+ temp = T[i];
+ double psat = w->satPressure(temp);
+ double press = oneBar;
+ if (psat > press) {
+ press = psat*1.002;
+ }
+ w->setState_TP(temp, press);
+ w->getPartialMolarEnthalpies(&h);
+ delh0 = tvalue(h - h298l, 1.0E-6);
+ w->getChemPotentials(&g);
+ delg0 = (g - h298l)/temp;
+ w->getPartialMolarCp(&Cp0);
+ w->getPartialMolarEntropies(&s);
+ w->getPartialMolarVolumes(&vol);
+ printf("%10g %10g %12g %13g %13g %13g %13g %13g\n", temp, press*1.0E-5,
+ psat*1.0E-5,
+ Cp0*1.0E-3, s*1.0E-3,
+ -delg0*1.0E-3, delh0*1.0E-6, vol);
+ }
+
+ printf("\nLiquid 1bar or psat State: Standard State Quantities\n");
+ printf (" T press psat Cpbar Sbar "
+ " -(G0-H298)/T H0-H298 Volume\n");
+ printf (" (K) (bar) (bar) (J/molK) (J/molK)"
+ " (J/molK) (kJ/mol) m3/kmol\n");
+
+ for (int i = 0; i < 6; i++) {
+ temp = T[i];
+ double psat = w->satPressure(temp);
+ double press = oneBar;
+ if (psat > press) {
+ press = psat*1.002;
+ }
+ w->setState_TP(temp, press);
+ w->getEnthalpy_RT(&h);
+ h *= temp * GasConstant;
+ delh0 = tvalue(h - h298l, 1.0E-6);
+ w->getStandardChemPotentials(&g);
+ delg0 = (g - h298l)/temp;
+ w->getCp_R(&Cp0);
+ Cp0 *= GasConstant;
+ w->getEntropy_R(&s);
+ s *= GasConstant;
+ w->getStandardVolumes(&vol);
+ printf("%10g %10g %12g %13g %13g %13g %13g %13g\n", temp, press*1.0E-5,
+ psat*1.0E-5,
+ Cp0*1.0E-3, s*1.0E-3,
+ -delg0*1.0E-3, delh0*1.0E-6, vol);
+ }
+
+ printf("\nLiquid 1bar or psat State: Reference State Quantities (Always 1 atm no matter what system pressure is)\n");
+ printf (" T press psat Cpbar Sbar "
+ " -(G0-H298)/T H0-H298 Volume\n");
+ printf (" (K) (bar) (bar) (J/molK) (J/molK)"
+ " (J/molK) (kJ/mol) m3/kmol\n");
+
+ for (int i = 0; i < 6; i++) {
+ temp = T[i];
+ double psat = w->satPressure(temp);
+ double press = oneBar;
+ if (psat > press) {
+ press = psat*1.002;
+ }
+ w->setState_TP(temp, press);
+ w->getEnthalpy_RT_ref(&h);
+ h *= temp * GasConstant;
+ delh0 = tvalue(h - h298l, 1.0E-6);
+ w->getGibbs_ref(&g);
+ delg0 = (g - h298l)/temp;
+ w->getCp_R_ref(&Cp0);
+ Cp0 *= GasConstant;
+ w->getEntropy_R_ref(&s);
+ s *= GasConstant;
+ w->getStandardVolumes_ref(&vol);
+ printf("%10g %10g %12g %13g %13g %13g %13g %13g\n", temp, press*1.0E-5,
+ psat*1.0E-5,
+ Cp0*1.0E-3, s*1.0E-3,
+ -delg0*1.0E-3, delh0*1.0E-6, vol);
+ }
+
+ printf("\nLiquid 1 atm: Standard State Quantities - Should agree with table above\n");
+ printf (" T press psat Cpbar Sbar "
+ " -(G0-H298)/T H0-H298 Volume\n");
+ printf (" (K) (bar) (bar) (J/molK) (J/molK)"
+ " (J/molK) (kJ/mol) m3/kmol\n");
+
+ for (int i = 0; i < 6; i++) {
+ temp = T[i];
+ double psat = w->satPressure(temp);
+ double press = OneAtm;
+ w->setState_TP(temp, press);
+ w->getEnthalpy_RT(&h);
+ h *= temp * GasConstant;
+ delh0 = tvalue(h - h298l, 1.0E-6);
+ w->getStandardChemPotentials(&g);
+ delg0 = (g - h298l)/temp;
+ w->getCp_R(&Cp0);
+ Cp0 *= GasConstant;
+ w->getEntropy_R(&s);
+ s *= GasConstant;
+ w->getStandardVolumes(&vol);
+ printf("%10g %10g %12g %13g %13g %13g %13g %13g\n", temp, press*1.0E-5,
+ psat*1.0E-5,
+ Cp0*1.0E-3, s*1.0E-3,
+ -delg0*1.0E-3, delh0*1.0E-6, vol);
+ }
printf("\n\nTable of increasing Enthalpy at 1 atm\n\n");
double dens;
printf(" Enthalpy, Temperature, x_Vapor, Density, Entropy_mass, Gibbs_mass\n");