From c5f28afd19c37c1931cfbe4b8582d39ad31ed0a7 Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Fri, 21 Dec 2007 02:58:41 +0000 Subject: [PATCH] Made the report function a virtual function within the ThermoPhase class. This means that the function can be modified depending on the inheritance. This allows for the printout of molalities and pH's for liquid phases where appropriate. Took out a duplicate function, report(), in the Cantera namespace. There was one in cxxutils.cpp and one in phasereport.cpp. --- Cantera/cxx/include/GRI30.h | 3 +- Cantera/cxx/include/IdealGasMix.h | 3 +- Cantera/cxx/include/Metal.h | 2 +- Cantera/cxx/include/PureFluid.h | 4 +- Cantera/cxx/include/importPhase.h | 2 + Cantera/cxx/src/Makefile.in | 2 +- Cantera/cxx/src/cxxutils.cpp | 101 -------- Cantera/src/thermo/MolalityVPSSTP.cpp | 131 +++++++++- Cantera/src/thermo/MolalityVPSSTP.h | 7 + Cantera/src/thermo/PureFluidPhase.cpp | 120 +++++++++ Cantera/src/thermo/PureFluidPhase.h | 11 +- Cantera/src/thermo/ThermoPhase.cpp | 114 +++++++++ Cantera/src/thermo/ThermoPhase.h | 7 + Cantera/src/thermo/phasereport.cpp | 185 +++----------- test_problems/VCSnonideal/NaCl_equil/README | 3 +- .../VCSnonideal/NaCl_equil/good_dout.txt | 29 ++- .../VCSnonideal/NaCl_equil/good_out.txt | 237 +++++++++++++++++- test_problems/VCSnonideal/NaCl_equil/runtest | 2 +- 18 files changed, 683 insertions(+), 280 deletions(-) diff --git a/Cantera/cxx/include/GRI30.h b/Cantera/cxx/include/GRI30.h index ea9f30b99..9edd9da54 100644 --- a/Cantera/cxx/include/GRI30.h +++ b/Cantera/cxx/include/GRI30.h @@ -9,6 +9,7 @@ #include "kernel/IdealGasPhase.h" #include "kernel/GRI_30_Kinetics.h" #include "kernel/importKinetics.h" +#include "kernel/stringUtils.h" namespace Cantera { @@ -29,7 +30,7 @@ namespace Cantera { bool operator!() { return !m_ok;} bool ready() const { return m_ok; } friend std::ostream& operator<<(std::ostream& s, GRI30& mix) { - std::string r = report(mix, true); + std::string r = Cantera::report(mix, true); s << r; return s; } diff --git a/Cantera/cxx/include/IdealGasMix.h b/Cantera/cxx/include/IdealGasMix.h index f0593f444..006a30ac9 100644 --- a/Cantera/cxx/include/IdealGasMix.h +++ b/Cantera/cxx/include/IdealGasMix.h @@ -6,6 +6,7 @@ #include "kernel/IdealGasPhase.h" #include "kernel/GasKinetics.h" #include "kernel/importKinetics.h" +#include "kernel/stringUtils.h" namespace Cantera { @@ -35,7 +36,7 @@ namespace Cantera { bool operator!() { return !m_ok;} bool ready() const { return m_ok; } friend std::ostream& operator<<(std::ostream& s, IdealGasMix& mix) { - std::string r = report(mix, true); + std::string r = Cantera::report(mix, true); s << r; return s; } diff --git a/Cantera/cxx/include/Metal.h b/Cantera/cxx/include/Metal.h index cca906e7e..3beaefb06 100644 --- a/Cantera/cxx/include/Metal.h +++ b/Cantera/cxx/include/Metal.h @@ -27,7 +27,7 @@ namespace Cantera { bool ready() const { return m_ok; } //friend std::ostream& operator<<(std::ostream& s, IdealGasMix& mix) { - // std::string r = report(mix, true); + // std::string r = Cantera::report(mix, true); // s << r; // return s; diff --git a/Cantera/cxx/include/PureFluid.h b/Cantera/cxx/include/PureFluid.h index 051caeb40..ce7ec4c9e 100644 --- a/Cantera/cxx/include/PureFluid.h +++ b/Cantera/cxx/include/PureFluid.h @@ -8,6 +8,8 @@ #include "kernel/PureFluidPhase.h" #include "kinetics.h" +#include "kernel/stringUtils.h" + namespace Cantera { @@ -36,7 +38,7 @@ namespace Cantera { bool operator!() { return !m_ok;} bool ready() const { return m_ok; } friend std::ostream& operator<<(std::ostream& s, PureFluid& mix) { - std::string r = report(mix, true); + std::string r = Cantera::report(mix, true); s << r; return s; } diff --git a/Cantera/cxx/include/importPhase.h b/Cantera/cxx/include/importPhase.h index 94e86f4be..2d8d5ce30 100644 --- a/Cantera/cxx/include/importPhase.h +++ b/Cantera/cxx/include/importPhase.h @@ -13,6 +13,8 @@ namespace Cantera { ThermoPhase* importPhase(std::string infile, std::string id=""); + // -> this is a duplicate of a src/thermo/phasereport function + // We'll leave it here so that these are available externally std::string report(const ThermoPhase& th, bool show_thermo); std::string formatCompList(const Phase& mix, int xyc); } diff --git a/Cantera/cxx/src/Makefile.in b/Cantera/cxx/src/Makefile.in index 15136f6be..6b6fc24d3 100644 --- a/Cantera/cxx/src/Makefile.in +++ b/Cantera/cxx/src/Makefile.in @@ -15,7 +15,7 @@ SUFFIXES= .cpp .d .o PIC_FLAG=@PIC@ CXX_FLAGS = @CXXFLAGS@ $(CXX_OPT) $(PIC_FLAG) -OBJS = cxxutils.o importPhase.o +OBJS = importPhase.o DEPENDS = $(OBJS:.o=.d) diff --git a/Cantera/cxx/src/cxxutils.cpp b/Cantera/cxx/src/cxxutils.cpp index 70eb60ecd..cedaddbaf 100644 --- a/Cantera/cxx/src/cxxutils.cpp +++ b/Cantera/cxx/src/cxxutils.cpp @@ -11,107 +11,6 @@ using namespace std; namespace Cantera { - /** - * Format a summary of the mixture state for output. - */ - std::string report(const ThermoPhase& th, bool show_thermo) { - - try { - char p[200]; - string s = ""; - - sprintf(p, "\n temperature %12.6g K\n", th.temperature()); - s += p; - sprintf(p, " pressure %12.6g Pa\n", th.pressure()); - s += p; - sprintf(p, " density %12.6g kg/m^3\n", th.density()); - s += p; - sprintf(p, " mean mol. weight %12.6g amu\n", th.meanMolecularWeight()); - s += p; - - if (show_thermo) { - sprintf(p, "\n"); - s += p; - sprintf(p, " 1 kg 1 kmol\n"); - s += p; - sprintf(p, " ----------- ------------\n"); - s += p; - sprintf(p, " enthalpy %12.6g %12.4g J\n", - th.enthalpy_mass(), th.enthalpy_mole()); - s += p; - sprintf(p, " internal energy %12.6g %12.4g J\n", - th.intEnergy_mass(), th.intEnergy_mole()); - s += p; - sprintf(p, " entropy %12.6g %12.4g J/K\n", - th.entropy_mass(), th.entropy_mole()); - s += p; - sprintf(p, " Gibbs function %12.6g %12.4g J\n", - th.gibbs_mass(), th.gibbs_mole()); - s += p; - sprintf(p, " heat capacity c_p %12.6g %12.4g J/K\n", - th.cp_mass(), th.cp_mole()); - s += p; - sprintf(p, " heat capacity c_v %12.6g %12.4g J/K\n", - th.cv_mass(), th.cv_mole()); - s += p; - } - - int kk = th.nSpecies(); - array_fp x(kk); - array_fp y(kk); - th.getMoleFractions(DATA_PTR(x)); - th.getMassFractions(DATA_PTR(y)); - - int k; - - sprintf(p, "\n X Y \n"); - s += p; - sprintf(p, " ------------- ------------\n"); - s += p; - for (k = 0; k < kk; k++) { - sprintf(p, "%18s %12.6e %12.6e\n", - th.speciesName(k).c_str(), x[k], y[k]); - s += p; - } - return s; - } - catch (CanteraError) { - return std::string(""); - } - } - - /** - * Format a composition list for output. - */ - std::string formatCompList(const Phase& mix, int xyc) { - - const doublereal Threshold = 1.e-20; - - char p[200]; - std::string s = ""; - int kk = mix.nSpecies(); - array_fp zz(kk); - switch (xyc) { - case 0: mix.getMoleFractions(DATA_PTR(zz)); break; - case 1: mix.getMassFractions(DATA_PTR(zz)); break; - case 2: mix.getConcentrations(DATA_PTR(zz)); break; - default: return "error: xyc must be 0, 1, or 2"; - } - - doublereal z; - int k; - for (k = 0; k < kk; k++) { - z = fabs(zz[k]); - if (z < Threshold) zz[k] = 0.0; - } - - for (k = 0; k < kk; k++) { - sprintf(p, "%18s\t %12.6e\n", mix.speciesName(k).c_str(), - zz[k]); - s += p; - } - return s; - } } diff --git a/Cantera/src/thermo/MolalityVPSSTP.cpp b/Cantera/src/thermo/MolalityVPSSTP.cpp index cbfe5df3d..a220e8ffe 100644 --- a/Cantera/src/thermo/MolalityVPSSTP.cpp +++ b/Cantera/src/thermo/MolalityVPSSTP.cpp @@ -595,9 +595,134 @@ namespace Cantera { VPStandardStateTP::initThermoXML(phaseNode, id); } + /** + * Format a summary of the mixture state for output. + */ + std::string MolalityVPSSTP::report(bool show_thermo) const { + + char p[800]; + string s = ""; + try { + if (name() != "") { + sprintf(p, " \n %s:\n", name().c_str()); + s += p; + } + sprintf(p, " \n temperature %12.6g K\n", temperature()); + s += p; + sprintf(p, " pressure %12.6g Pa\n", pressure()); + s += p; + sprintf(p, " density %12.6g kg/m^3\n", density()); + s += p; + sprintf(p, " mean mol. weight %12.6g amu\n", meanMolecularWeight()); + s += p; + + doublereal phi = electricPotential(); + sprintf(p, " potential %12.6g V\n", phi); + s += p; + + int kk = nSpecies(); + array_fp x(kk); + array_fp molal(kk); + array_fp mu(kk); + array_fp muss(kk); + array_fp acMolal(kk); + array_fp actMolal(kk); + getMoleFractions(&x[0]); + getMolalities(&molal[0]); + getChemPotentials(&mu[0]); + getStandardChemPotentials(&muss[0]); + getMolalityActivityCoefficients(&acMolal[0]); + getActivities(&actMolal[0]); + + int iHp = speciesIndex("H+"); + if (iHp >= 0) { + double pH = -log(actMolal[iHp]) / log(10.0); + sprintf(p, " pH %12.4g \n", pH); + s += p; + } + + if (show_thermo) { + sprintf(p, " \n"); + s += p; + sprintf(p, " 1 kg 1 kmol\n"); + s += p; + sprintf(p, " ----------- ------------\n"); + s += p; + sprintf(p, " enthalpy %12.6g %12.4g J\n", + enthalpy_mass(), enthalpy_mole()); + s += p; + sprintf(p, " internal energy %12.6g %12.4g J\n", + intEnergy_mass(), intEnergy_mole()); + s += p; + sprintf(p, " entropy %12.6g %12.4g J/K\n", + entropy_mass(), entropy_mole()); + s += p; + sprintf(p, " Gibbs function %12.6g %12.4g J\n", + gibbs_mass(), gibbs_mole()); + s += p; + sprintf(p, " heat capacity c_p %12.6g %12.4g J/K\n", + cp_mass(), cp_mole()); + s += p; + try { + sprintf(p, " heat capacity c_v %12.6g %12.4g J/K\n", + cv_mass(), cv_mole()); + s += p; + } + catch(CanteraError) { + sprintf(p, " heat capacity c_v \n"); + s += p; + } + } + + + //doublereal rt = GasConstant * temperature(); + int k; + + + sprintf(p, " \n"); + s += p; + if (show_thermo) { + sprintf(p, " X " + " Molalities Chem.Pot. ChemPotSS ActCoeffMolal\n"); + s += p; + sprintf(p, " " + " (J/kmol) (J/kmol) \n"); + s += p; + sprintf(p, " ------------- " + " ------------ ------------ ------------ ------------\n"); + s += p; + for (k = 0; k < kk; k++) { + if (x[k] > SmallNumber) { + sprintf(p, "%18s %12.6g %12.6g %12.6g %12.6g %12.6g\n", + speciesName(k).c_str(), x[k], molal[k], mu[k], muss[k], acMolal[k]); + } + else { + sprintf(p, "%18s %12.6g %12.6g N/A %12.6g %12.6g \n", + speciesName(k).c_str(), x[k], molal[k], muss[k], acMolal[k]); + } + s += p; + } + } + else { + sprintf(p, " X" + "Molalities\n"); + s += p; + sprintf(p, " -------------" + " ------------\n"); + s += p; + for (k = 0; k < kk; k++) { + sprintf(p, "%18s %12.6g %12.6g\n", + speciesName(k).c_str(), x[k], molal[k]); + s += p; + } + } + } catch (CanteraError) { + ; + } + return s; + } + + } - - - diff --git a/Cantera/src/thermo/MolalityVPSSTP.h b/Cantera/src/thermo/MolalityVPSSTP.h index 32f0893f4..58056f713 100644 --- a/Cantera/src/thermo/MolalityVPSSTP.h +++ b/Cantera/src/thermo/MolalityVPSSTP.h @@ -744,6 +744,13 @@ namespace Cantera { */ void setState_TPM(doublereal t, doublereal p, const std::string& m); + //! returns a summary of the state of the phase as a string + /*! + * @param show_thermo If true, extra information is printed out + * about the thermodynamic state of the system. + */ + virtual std::string report(bool show_thermo = true) const; + private: void initLengths(); diff --git a/Cantera/src/thermo/PureFluidPhase.cpp b/Cantera/src/thermo/PureFluidPhase.cpp index 6f8de086f..80034d987 100644 --- a/Cantera/src/thermo/PureFluidPhase.cpp +++ b/Cantera/src/thermo/PureFluidPhase.cpp @@ -307,6 +307,126 @@ namespace Cantera { check(); } + + /** + * Format a summary of the mixture state for output. + */ + std::string PureFluidPhase::report(bool show_thermo) const { + + + char p[800]; + string s = ""; + try { + if (name() != "") { + sprintf(p, " \n %s:\n", name().c_str()); + s += p; + } + sprintf(p, " \n temperature %12.6g K\n", temperature()); + s += p; + sprintf(p, " pressure %12.6g Pa\n", pressure()); + s += p; + sprintf(p, " density %12.6g kg/m^3\n", density()); + s += p; + sprintf(p, " mean mol. weight %12.6g amu\n", meanMolecularWeight()); + s += p; + + if (eosType() == cPureFluid) { + double xx = ((PureFluidPhase *) (this))->vaporFraction(); + sprintf(p, " vapor fraction %12.6g \n", + xx); //th.vaporFraction()); + s += p; + } + + doublereal phi = electricPotential(); + if (phi != 0.0) { + sprintf(p, " potential %12.6g V\n", phi); + s += p; + } + if (show_thermo) { + sprintf(p, " \n"); + s += p; + sprintf(p, " 1 kg 1 kmol\n"); + s += p; + sprintf(p, " ----------- ------------\n"); + s += p; + sprintf(p, " enthalpy %12.6g %12.4g J\n", + enthalpy_mass(), enthalpy_mole()); + s += p; + sprintf(p, " internal energy %12.6g %12.4g J\n", + intEnergy_mass(), intEnergy_mole()); + s += p; + sprintf(p, " entropy %12.6g %12.4g J/K\n", + entropy_mass(), entropy_mole()); + s += p; + sprintf(p, " Gibbs function %12.6g %12.4g J\n", + gibbs_mass(), gibbs_mole()); + s += p; + sprintf(p, " heat capacity c_p %12.6g %12.4g J/K\n", + cp_mass(), cp_mole()); + s += p; + try { + sprintf(p, " heat capacity c_v %12.6g %12.4g J/K\n", + cv_mass(), cv_mole()); + s += p; + } + catch(CanteraError) { + sprintf(p, " heat capacity c_v \n"); + s += p; + } + } + + int kk = nSpecies(); + array_fp x(kk); + array_fp y(kk); + array_fp mu(kk); + getMoleFractions(&x[0]); + getMassFractions(&y[0]); + getChemPotentials(&mu[0]); + doublereal rt = GasConstant * temperature(); + int k; + //if (th.nSpecies() > 1) { + + if (show_thermo) { + sprintf(p, " \n X " + " Y Chem. Pot. / RT \n"); + s += p; + sprintf(p, " ------------- " + "------------ ------------\n"); + s += p; + for (k = 0; k < kk; k++) { + if (x[k] > SmallNumber) { + sprintf(p, "%18s %12.6g %12.6g %12.6g\n", + speciesName(k).c_str(), x[k], y[k], mu[k]/rt); + } + else { + sprintf(p, "%18s %12.6g %12.6g \n", + speciesName(k).c_str(), x[k], y[k]); + } + s += p; + } + } + else { + sprintf(p, " \n X" + "Y\n"); + s += p; + sprintf(p, " -------------" + " ------------\n"); + s += p; + for (k = 0; k < kk; k++) { + sprintf(p, "%18s %12.6g %12.6g\n", + speciesName(k).c_str(), x[k], y[k]); + s += p; + } + } + } + //} + catch (CanteraError) { + ; + } + return s; + } + + } #endif // WITH_PURE_FLUIDS diff --git a/Cantera/src/thermo/PureFluidPhase.h b/Cantera/src/thermo/PureFluidPhase.h index 73723f444..8512a9022 100644 --- a/Cantera/src/thermo/PureFluidPhase.h +++ b/Cantera/src/thermo/PureFluidPhase.h @@ -1,6 +1,7 @@ /** * @file PureFluidPhase.h - * Header for a ThermoPhase object for a pure fluid phase consisting of gas, liquid, mixed-gas-liquid + * Header for a ThermoPhase object for a pure fluid phase consisting of + * gas, liquid, mixed-gas-liquid * and supercrit fluid (see \ref thermoprops * and class \link Cantera::PureFluidPhase PureFluidPhase\endlink). * @@ -294,6 +295,14 @@ namespace Cantera { */ virtual void setParametersFromXML(const XML_Node& eosdata); + + //! returns a summary of the state of the phase as a string + /*! + * @param show_thermo If true, extra information is printed out + * about the thermodynamic state of the system. + */ + virtual std::string report(bool show_thermo = true) const; + protected: //! Main call to the tpx level to set the state of the system diff --git a/Cantera/src/thermo/ThermoPhase.cpp b/Cantera/src/thermo/ThermoPhase.cpp index 980c943f2..afec8b3a0 100644 --- a/Cantera/src/thermo/ThermoPhase.cpp +++ b/Cantera/src/thermo/ThermoPhase.cpp @@ -838,4 +838,118 @@ namespace Cantera { return (m_hasElementPotentials); } + /** + * Format a summary of the mixture state for output. + */ + std::string ThermoPhase::report(bool show_thermo) const { + + + char p[800]; + string s = ""; + try { + if (name() != "") { + sprintf(p, " \n %s:\n", name().c_str()); + s += p; + } + sprintf(p, " \n temperature %12.6g K\n", temperature()); + s += p; + sprintf(p, " pressure %12.6g Pa\n", pressure()); + s += p; + sprintf(p, " density %12.6g kg/m^3\n", density()); + s += p; + sprintf(p, " mean mol. weight %12.6g amu\n", meanMolecularWeight()); + s += p; + + doublereal phi = electricPotential(); + if (phi != 0.0) { + sprintf(p, " potential %12.6g V\n", phi); + s += p; + } + if (show_thermo) { + sprintf(p, " \n"); + s += p; + sprintf(p, " 1 kg 1 kmol\n"); + s += p; + sprintf(p, " ----------- ------------\n"); + s += p; + sprintf(p, " enthalpy %12.6g %12.4g J\n", + enthalpy_mass(), enthalpy_mole()); + s += p; + sprintf(p, " internal energy %12.6g %12.4g J\n", + intEnergy_mass(), intEnergy_mole()); + s += p; + sprintf(p, " entropy %12.6g %12.4g J/K\n", + entropy_mass(), entropy_mole()); + s += p; + sprintf(p, " Gibbs function %12.6g %12.4g J\n", + gibbs_mass(), gibbs_mole()); + s += p; + sprintf(p, " heat capacity c_p %12.6g %12.4g J/K\n", + cp_mass(), cp_mole()); + s += p; + try { + sprintf(p, " heat capacity c_v %12.6g %12.4g J/K\n", + cv_mass(), cv_mole()); + s += p; + } + catch(CanteraError) { + sprintf(p, " heat capacity c_v \n"); + s += p; + } + } + + int kk = nSpecies(); + array_fp x(kk); + array_fp y(kk); + array_fp mu(kk); + getMoleFractions(&x[0]); + getMassFractions(&y[0]); + getChemPotentials(&mu[0]); + doublereal rt = GasConstant * temperature(); + int k; + //if (th.nSpecies() > 1) { + + if (show_thermo) { + sprintf(p, " \n X " + " Y Chem. Pot. / RT \n"); + s += p; + sprintf(p, " ------------- " + "------------ ------------\n"); + s += p; + for (k = 0; k < kk; k++) { + if (x[k] > SmallNumber) { + sprintf(p, "%18s %12.6g %12.6g %12.6g\n", + speciesName(k).c_str(), x[k], y[k], mu[k]/rt); + } + else { + sprintf(p, "%18s %12.6g %12.6g \n", + speciesName(k).c_str(), x[k], y[k]); + } + s += p; + } + } + else { + sprintf(p, " \n X" + "Y\n"); + s += p; + sprintf(p, " -------------" + " ------------\n"); + s += p; + for (k = 0; k < kk; k++) { + sprintf(p, "%18s %12.6g %12.6g\n", + speciesName(k).c_str(), x[k], y[k]); + s += p; + } + } + } + //} + catch (CanteraError) { + ; + } + return s; + } + + + + } diff --git a/Cantera/src/thermo/ThermoPhase.h b/Cantera/src/thermo/ThermoPhase.h index a96edcdf7..bdfb17e55 100755 --- a/Cantera/src/thermo/ThermoPhase.h +++ b/Cantera/src/thermo/ThermoPhase.h @@ -1501,6 +1501,13 @@ namespace Cantera { return m_chargeNeutralityNecessary; } + + //! returns a summary of the state of the phase as a string + /*! + * @param show_thermo If true, extra information is printed out + * about the thermodynamic state of the system. + */ + virtual std::string report(bool show_thermo = true) const; protected: diff --git a/Cantera/src/thermo/phasereport.cpp b/Cantera/src/thermo/phasereport.cpp index a944313ac..8eb90710c 100644 --- a/Cantera/src/thermo/phasereport.cpp +++ b/Cantera/src/thermo/phasereport.cpp @@ -14,161 +14,50 @@ using namespace std; namespace Cantera { - /** - * Format a summary of the mixture state for output. - */ - string report(const ThermoPhase& th, bool show_thermo) { - char p[200]; - string s = ""; - try { - if (th.name() != "") { - sprintf(p, " \n %s:\n", th.name().c_str()); - s += p; - } - sprintf(p, " \n temperature %12.6g K\n", th.temperature()); - s += p; - sprintf(p, " pressure %12.6g Pa\n", th.pressure()); - s += p; - sprintf(p, " density %12.6g kg/m^3\n", th.density()); - s += p; - sprintf(p, " mean mol. weight %12.6g amu\n", th.meanMolecularWeight()); - s += p; -#ifdef WITH_PURE_FLUIDS - if (th.eosType() == cPureFluid) { - double xx = ((PureFluidPhase*)(&th))->vaporFraction(); - // if (th.temperature() < th.critTemperature()) { - sprintf(p, " vapor fraction %12.6g \n", - xx); //th.vaporFraction()); - s += p; - //} - } -#endif - doublereal phi = th.electricPotential(); - if (phi != 0.0) { - sprintf(p, " potential %12.6g V\n", phi); - s += p; - } - if (show_thermo) { - sprintf(p, " \n"); - s += p; - sprintf(p, " 1 kg 1 kmol\n"); - s += p; - sprintf(p, " ----------- ------------\n"); - s += p; - sprintf(p, " enthalpy %12.6g %12.4g J\n", - th.enthalpy_mass(), th.enthalpy_mole()); - s += p; - sprintf(p, " internal energy %12.6g %12.4g J\n", - th.intEnergy_mass(), th.intEnergy_mole()); - s += p; - sprintf(p, " entropy %12.6g %12.4g J/K\n", - th.entropy_mass(), th.entropy_mole()); - s += p; - sprintf(p, " Gibbs function %12.6g %12.4g J\n", - th.gibbs_mass(), th.gibbs_mole()); - s += p; - sprintf(p, " heat capacity c_p %12.6g %12.4g J/K\n", - th.cp_mass(), th.cp_mole()); - s += p; - try { - sprintf(p, " heat capacity c_v %12.6g %12.4g J/K\n", - th.cv_mass(), th.cv_mole()); - s += p; - } - catch(CanteraError) { - sprintf(p, " heat capacity c_v \n"); - s += p; - } - } + /** + * Format a summary of the mixture state for output. + */ + std::string report(const ThermoPhase& th, bool show_thermo) { + return th.report(show_thermo); + } - int kk = th.nSpecies(); - array_fp x(kk); - array_fp y(kk); - array_fp mu(kk); - th.getMoleFractions(&x[0]); - th.getMassFractions(&y[0]); - th.getChemPotentials(&mu[0]); - doublereal rt = GasConstant * th.temperature(); - int k; - //if (th.nSpecies() > 1) { + void writephase(const ThermoPhase& th, bool show_thermo) { + string s = report(th, show_thermo); + writelog(s+"\n"); + } - if (show_thermo) { - sprintf(p, " \n X " - " Y Chem. Pot. / RT \n"); - s += p; - sprintf(p, " ------------- " - "------------ ------------\n"); - s += p; - for (k = 0; k < kk; k++) { - if (x[k] > SmallNumber) { - sprintf(p, "%18s %12.6g %12.6g %12.6g\n", - th.speciesName(k).c_str(), x[k], y[k], mu[k]/rt); - } - else { - sprintf(p, "%18s %12.6g %12.6g \n", - th.speciesName(k).c_str(), x[k], y[k]); - } - s += p; - } - } - else { - sprintf(p, " \n X" - "Y\n"); - s += p; - sprintf(p, " -------------" - " ------------\n"); - s += p; - for (k = 0; k < kk; k++) { - sprintf(p, "%18s %12.6g %12.6g\n", - th.speciesName(k).c_str(), x[k], y[k]); - s += p; - } - } - } - //} - catch (CanteraError) { - ; - } - return s; + /** + * Format a composition list for output. + */ + std::string formatCompList(const Phase& mix, int xyc) { + + const doublereal Threshold = 1.e-20; + + char p[200]; + string s = ""; + int kk = mix.nSpecies(); + array_fp zz(kk); + switch (xyc) { + case 0: mix.getMoleFractions(&zz[0]); break; + case 1: mix.getMassFractions(&zz[0]); break; + case 2: mix.getConcentrations(&zz[0]); break; + default: return "error: xyc must be 0, 1, or 2"; } - void writephase(const ThermoPhase& th, bool show_thermo) { - string s = report(th, show_thermo); - writelog(s+"\n"); + doublereal z; + int k; + for (k = 0; k < kk; k++) { + z = fabs(zz[k]); + if (z < Threshold) zz[k] = 0.0; } - /** - * Format a composition list for output. - */ - string formatCompList(const Phase& mix, int xyc) { - - const doublereal Threshold = 1.e-20; - - char p[200]; - string s = ""; - int kk = mix.nSpecies(); - array_fp zz(kk); - switch (xyc) { - case 0: mix.getMoleFractions(&zz[0]); break; - case 1: mix.getMassFractions(&zz[0]); break; - case 2: mix.getConcentrations(&zz[0]); break; - default: return "error: xyc must be 0, 1, or 2"; - } - - doublereal z; - int k; - for (k = 0; k < kk; k++) { - z = fabs(zz[k]); - if (z < Threshold) zz[k] = 0.0; - } - - for (k = 0; k < kk; k++) { - sprintf(p, "%18s\t %12.6e\n", mix.speciesName(k).c_str(), - zz[k]); - s += p; - } - return s; + for (k = 0; k < kk; k++) { + sprintf(p, "%18s\t %12.6e\n", mix.speciesName(k).c_str(), + zz[k]); + s += p; } + return s; + } } diff --git a/test_problems/VCSnonideal/NaCl_equil/README b/test_problems/VCSnonideal/NaCl_equil/README index 612af0d9a..fb6dec2aa 100644 --- a/test_problems/VCSnonideal/NaCl_equil/README +++ b/test_problems/VCSnonideal/NaCl_equil/README @@ -12,6 +12,8 @@ The problem can be divided up into two parts: estimating the Gibbs reaction delta for Na+ + Cl- = NaCl(solid) + (so that Del(Gf) = G(NaCl(solid)) - G(Na+) - G(Cl-)) + and estimating the activity coefficients for the electrolytes at the solubility limit. @@ -39,7 +41,6 @@ From Codata key values for Thermodynamics: - In addition, the relative humidity of the salt solution may be compared to the humidity above a pure water solution in order to understand the effects of the lowering of the water activity. diff --git a/test_problems/VCSnonideal/NaCl_equil/good_dout.txt b/test_problems/VCSnonideal/NaCl_equil/good_dout.txt index 3028e0dfa..9f1553c30 100644 --- a/test_problems/VCSnonideal/NaCl_equil/good_dout.txt +++ b/test_problems/VCSnonideal/NaCl_equil/good_dout.txt @@ -2814,16 +2814,16 @@ Chemical Potentials of the Species: (dimensionless) Counters: Iterations Time (seconds) - vcs_basopt: 3 0.00000E+00 - vcs_TP: 36 1.75000E+00 + vcs_basopt: 3 2.00000E-02 + vcs_TP: 36 1.74000E+00 -------------------------------------------------------------------------------- -------------------------------------------------------------------------------- TCounters: Num_Calls Total_Its Total_Time (seconds) - vcs_basopt: 3 3 0.00000E+00 - vcs_TP: 1 36 1.75000E+00 + vcs_basopt: 3 3 2.00000E-02 + vcs_TP: 1 36 1.74000E+00 vcs_inest: 0 0.00000E+00 - vcs_TotalTime: 1.77000E+00 + vcs_TotalTime: 1.75000E+00 Results from vcs: @@ -2847,7 +2847,7 @@ N2 4.000e+00 9.763e-01 -5.719e+07 OH 3.747e-07 9.146e-08 -2.670e+08 NaCl(S) 4.788e+00 1.000e+00 -4.326e+08 ------------------------------------------------------------- -Total time = 1.780000e+00 seconds +Total time = 1.770000e+00 seconds *************** NaCl_electrolyte ***************** Moles: 2.32742 @@ -2857,6 +2857,8 @@ Moles: 2.32742 pressure 101325 Pa density 1216.41 kg/m^3 mean mol. weight 20.0596 amu + potential 0 V + pH 6.59 1 kg 1 kmol ----------- ------------ @@ -2867,13 +2869,14 @@ Moles: 2.32742 heat capacity c_p 3053.16 6.125e+04 J/K heat capacity c_v - X Y Chem. Pot. / RT - ------------- ------------ ------------ - H2O(L) 0.817565 0.734244 -124.003 - Cl- 0.0912175 0.161217 -72.3775 - H+ 7.91271e-10 3.97376e-11 -15.1739 - Na+ 0.0912175 0.104539 -102.139 - OH- 7.91271e-10 6.70891e-10 -108.829 + X Molalities Chem.Pot. ChemPotSS ActCoeffMolal + (J/kmol) (J/kmol) + ------------- ------------ ------------ ------------ ------------ + H2O(L) 0.817565 55.5084 -3.07399e+08 -3.06686e+08 0.917303 + Cl- 0.0912175 6.1932 -1.79421e+08 -1.83974e+08 1.01322 + H+ 7.91271e-10 5.37232e-08 -3.76154e+07 0 4.78537 + Na+ 0.0912175 6.1932 -2.53199e+08 -2.57752e+08 1.01322 + OH- 7.91271e-10 5.37232e-08 -2.69784e+08 -2.26784e+08 0.545262 *************** air ***************** Moles: 4.09718 diff --git a/test_problems/VCSnonideal/NaCl_equil/good_out.txt b/test_problems/VCSnonideal/NaCl_equil/good_out.txt index 7b89012c2..7e6f8f4cc 100644 --- a/test_problems/VCSnonideal/NaCl_equil/good_out.txt +++ b/test_problems/VCSnonideal/NaCl_equil/good_out.txt @@ -1,3 +1,67 @@ +Unknown Cantera EOS to VCSnonideal: 45012 + +================================================================================ +================ Cantera_to_vprob: START OF PROBLEM STATEMENT ==================== +================================================================================ + Phase IDs of species + species phaseID phaseName Initial_Estimated_gMols + H2O(L) 0 NaCl_electrolyte 2000 + Cl- 0 NaCl_electrolyte 0 + H+ 0 NaCl_electrolyte 0 + Na+ 0 NaCl_electrolyte 0 + OH- 0 NaCl_electrolyte 0 + O2 1 air 0 + H2 1 air 0 + CO2 1 air 0 + H2O 1 air 0 + NaCl 1 air 0 + N2 1 air 4000 + OH 1 air 0 + NaCl(S) 2 NaCl(S) 5000 + +-------------------------------------------------------------------------------- + Information about phases + PhaseName PhaseNum SingSpec GasPhase EqnState NumSpec TMolesInert Tmoles(gmol) +NaCl_electrolyte 0 0 0 UnkType: -1 5 0.000000e+00 2.000000e+03 + air 1 0 1 Ideal Gas 7 0.000000e+00 4.000000e+03 + NaCl(S) 2 1 0 Stoich Sub 1 0.000000e+00 5.000000e+03 + +================================================================================ +================ Cantera_to_vprob: END OF PROBLEM STATEMENT ==================== +================================================================================ + + +================================================================================ +==================== Cantera_to_vprob: START OF PROBLEM STATEMENT ==================== +================================================================================ + + Phase IDs of species + species phaseID phaseName Initial_Estimated_gMols + H2O(L) 0 NaCl_electrolyte 2000 + Cl- 0 NaCl_electrolyte 0 + H+ 0 NaCl_electrolyte 0 + Na+ 0 NaCl_electrolyte 0 + OH- 0 NaCl_electrolyte 0 + O2 1 air 0 + H2 1 air 0 + CO2 1 air 0 + H2O 1 air 0 + NaCl 1 air 0 + N2 1 air 4000 + OH 1 air 0 + NaCl(S) 2 NaCl(S) 5000 + +-------------------------------------------------------------------------------- + Information about phases + PhaseName PhaseNum SingSpec GasPhase EqnState NumSpec TMolesInert Tmoles(gmol) +NaCl_electrolyte 0 0 0 UnkType: -1 5 0.000000e+00 2.000000e+03 + air 1 0 1 Ideal Gas 7 0.000000e+00 4.000000e+03 + NaCl(S) 2 1 0 Stoich Sub 1 0.000000e+00 5.000000e+03 + +================================================================================ +==================== Cantera_to_vprob: END OF PROBLEM STATEMENT ==================== +================================================================================ + ================================================================================ ==================== VCS_PROB: PROBLEM STATEMENT =============================== @@ -63,6 +127,162 @@ Chemical Potentials: (J/kmol) ==================== VCS_PROB: END OF PROBLEM STATEMENT ======================== ================================================================================ +VCS CALCULATION METHOD + + MultiPhase Object + + + 13 SPECIES 11 ELEMENTS 7 COMPONENTS + 5 PHASE1 SPECIES 7 PHASE2 SPECIES 1 SINGLE SPECIES PHASES + + PRESSURE 101325.000 ATM + TEMPERATURE 298.150 K + PHASE1 INERTS 0.000 + PHASE2 INERTS 0.000 + + ELEMENTAL ABUNDANCES CORRECT FROM ESTIMATE Type + + O 2.000000000000E+03 2.000000000000E+03 0 + H 4.000000000000E+03 4.000000000000E+03 0 + C 0.000000000000E+00 0.000000000000E+00 0 + N 8.000000000000E+03 8.000000000000E+03 0 + Na 5.000000000000E+03 5.000000000000E+03 0 + Cl 5.000000000000E+03 5.000000000000E+03 0 + cn 0.000000000000E+00 0.000000000000E+00 2 + Fe 0.000000000000E+00 0.000000000000E+00 0 + E 0.000000000000E+00 0.000000000000E+00 1 + Si 0.000000000000E+00 0.000000000000E+00 0 + Ca 0.000000000000E+00 0.000000000000E+00 0 + + USER ESTIMATE OF EQUILIBRIUM + Stan. Chem. Pot. in J/kmol + + SPECIES FORMULA VECTOR STAN_CHEM_POT EQUILIBRIUM_EST. Species_Type + + O H C N Na Cl cn Fe E Si Ca SI(I) + NaCl(S) 0 0 0 0 1 1 0 0 0 0 0 0 -4.32620E+08 5.00000E+03 Mol_Num + N2 0 0 0 2 0 0 0 0 0 0 0 2 -5.71282E+07 4.00000E+03 Mol_Num + H2O(L) 1 2 0 0 0 0 0 0 -0 0 0 1 -3.06686E+08 2.00000E+03 Mol_Num + H+ 0 1 0 0 0 0 1 0 -1 0 0 1 0.00000E+00 0.00000E+00 Mol_Num + Na+ 0 0 0 0 1 0 1 0 -1 0 0 1 -2.57752E+08 0.00000E+00 Mol_Num + OH 1 1 0 0 0 0 0 0 0 0 0 2 -2.26793E+08 0.00000E+00 Mol_Num + CO2 2 0 1 0 0 0 0 0 0 0 0 2 -4.57249E+08 0.00000E+00 Mol_Num + H2 0 2 0 0 0 0 0 0 0 0 0 2 -3.89624E+07 0.00000E+00 Mol_Num + H2O 1 2 0 0 0 0 0 0 0 0 0 2 -2.98124E+08 0.00000E+00 Mol_Num + NaCl 0 0 0 0 1 1 0 0 0 0 0 2 -2.49904E+08 0.00000E+00 Mol_Num + Cl- 0 0 0 0 0 1 -1 0 1 0 0 1 -1.83974E+08 0.00000E+00 Mol_Num + O2 2 0 0 0 0 0 0 0 0 0 0 2 -6.11650E+07 0.00000E+00 Mol_Num + OH- 1 1 0 0 0 0 -1 0 1 0 0 1 -2.26784E+08 0.00000E+00 Mol_Num + + + + +-------------------------------------------------------------------------------- +-------------------------------------------------------------------------------- + VCS_TP REPORT +-------------------------------------------------------------------------------- +-------------------------------------------------------------------------------- + Temperature = 3e+02 Kelvin + Pressure = 1.0132e+05 Atmos + Volume = 1.0041e+05 cm**3 + + +-------------------------------------------------------------------------------- + Species Equilibrium moles Mole Fraction ChemPot/RT SpecUnkType +-------------------------------------------------------------------------------- + NaCl(S) 4.7876981E+03 1.0000000E+00 -1.7452E+02 0 + N2 4.0000000E+03 9.7628143E-01 -2.3069E+01 0 + H2O(L) 1.9028209E+03 8.1756495E-01 -1.2400E+02 0 + Na+ 2.1230193E+02 9.1217525E-02 -1.0214E+02 0 + H+ 1.8416226E-06 7.9127051E-10 -1.5174E+01 0 + OH 3.7473450E-04 9.1461582E-08 -1.0769E+02 0 + CO2 0.0000000E+00 0.0000000E+00 -5.1513E+02 0 + Cl- 2.1230193E+02 9.1217525E-02 -7.2377E+01 MolNum + H2O 9.7178674E+01 2.3718434E-02 -1.2400E+02 MolNum + H2 1.8736725E-04 4.5730791E-08 -3.2618E+01 MolNum + OH- 1.8416226E-06 7.9127051E-10 -1.0883E+02 MolNum + NaCl 3.9996395E-29 9.7619345E-33 -1.7452E+02 MolNum + O2 8.9465959E-66 2.1835989E-69 -1.8277E+02 MolNum +-------------------------------------------------------------------------------- + + +------------------------------------------------------------------------------------------------------------------- + |ComponentID| 0 1 2 3 4 5 6 | | + | Components| NaCl(S) N2 H2O(L) Na+ H+ OH CO2 | | + NonComponent | Moles | 4.79e+03 4e+03 1.9e+03 212 1.84e-06 0.000375 0 | DG/RT Rxn | +------------------------------------------------------------------------------------------------------------------- + 7 O2 | 8.95e-66 | 0.00 0.00 2.00 0.00 0.00 -4.00 0.00 | 0 | + 8 H2O | 97.2 | 0.00 0.00 -1.00 0.00 0.00 0.00 0.00 | -3.15e-10 | + 9 NaCl | 4e-29 | -1.00 0.00 0.00 0.00 0.00 0.00 0.00 | -5.03e-11 | + 10 Cl- | 212 | -1.00 0.00 0.00 1.00 0.00 0.00 0.00 | 3.93e-10 | + 11 H2 | 0.000187 | 0.00 0.00 -2.00 0.00 0.00 2.00 0.00 | -6.93e-10 | + 12 OH- | 1.84e-06 | 0.00 0.00 -1.00 0.00 1.00 0.00 0.00 | -1.57e-09 | +------------------------------------------------------------------------------------------------------------------- + + + + +------------------------------------------------------------------------------------------------------------------------------------------------------------------------ + | ElementID | 0 1 2 3 4 5 6 7 8 9 10 | | + | Element | O H C N Na Cl cn_NaCl_el Fe E Si Ca | | + PhaseName | MolTarget | 2e+03 4e+03 0 8e+03 5e+03 5e+03 0 0 0 0 0 | Gibbs Total | +------------------------------------------------------------------------------------------------------------------------------------------------------------------------ + 0 NaCl_electro | 2.327e+03 | 1.9e+03 3.81e+03 0 0 212 212 -5.56e-15 0 5.56e-15 0 0 | -2.73006234953E+05 | + 1 air | 4.097e+03 | 97.2 194 0 8e+03 4e-29 4e-29 0 0 0 0 0 | -1.04327434124E+05 | + 2 NaCl(S) | 4.788e+03 | 0 0 0 0 4.79e+03 4.79e+03 0 0 0 0 0 | -8.35533676806E+05 | +------------------------------------------------------------------------------------------------------------------------------------------------------------------------ + TOTAL | 1.121e+04 | 2e+03 4e+03 0 8e+03 5e+03 5e+03 -5.56e-15 0 5.56e-15 0 0 | -1.21286734588E+06 | +------------------------------------------------------------------------------------------------------------------------------------------------------------------------ + + + Total Dimensionless Gibbs Free Energy = G/RT = -1.2128673E+06 + +Elemental Abundances: Actual Target Type ElActive + O 2.000000000000E+03 2.000000000000E+03 0 1 + H 4.000000000000E+03 4.000000000000E+03 0 1 + C 0.000000000000E+00 0.000000000000E+00 0 1 + N 8.000000000000E+03 8.000000000000E+03 0 1 + Na 5.000000000000E+03 5.000000000000E+03 0 1 + Cl 5.000000000000E+03 5.000000000000E+03 0 1 + cn -5.561230390582E-15 0.000000000000E+00 2 1 + Fe 0.000000000000E+00 0.000000000000E+00 0 1 + E 5.561230390582E-15 0.000000000000E+00 1 0 + Si 0.000000000000E+00 0.000000000000E+00 0 1 + Ca 0.000000000000E+00 0.000000000000E+00 0 1 + + +--------------------------------------------------------------------------------------------- +Chemical Potentials of the Species: (dimensionless) + (RT = 2.47896e+06 J/kmol) + Name TMoles StandStateChemPot ln(AC) ln(X_i) | F z_i phi | ChemPot | (-lnMnaught) +------------------------------------------------------------------------------------------------------------------- + NaCl(S) 4.7876981E+03 -1.7451679E+02 0.0000000E+00 0.0000000E+00 | 0.0000000E+00 | -1.7452E+02 | + N2 4.0000000E+03 -2.3045225E+01 0.0000000E+00 -2.4004385E-02 | 0.0000000E+00 | -2.3069E+01 | + H2O(L) 1.9028209E+03 -1.2371551E+02 -8.6317426E-02 -2.0142493E-01 | 0.0000000E+00 | -1.2400E+02 | + Na+ 2.1230193E+02 -1.0397591E+02 2.1455966E-01 -2.3945082E+00 | 0.0000000E+00 | -1.0214E+02 | ( 4.0165350E+00) + H+ 1.8416226E-06 0.0000000E+00 1.7669883E+00 -2.0957381E+01 | 0.0000000E+00 | -1.5174E+01 | ( 4.0165350E+00) + OH 3.7473450E-04 -9.1487045E+01 0.0000000E+00 -1.6207347E+01 | 0.0000000E+00 | -1.0769E+02 | + CO2 0.0000000E+00 -1.8445182E+02 0.0000000E+00 -3.3067997E+02 | 0.0000000E+00 | -5.1513E+02 | + Cl- 2.1230193E+02 -7.4214051E+01 2.1455966E-01 -2.3945082E+00 |-0.0000000E+00 | -7.2377E+01 | ( 4.0165350E+00) + H2O 9.7178674E+01 -1.2026175E+02 0.0000000E+00 -3.7415027E+00 | 0.0000000E+00 | -1.2400E+02 | + H2 1.8736725E-04 -1.5717224E+01 0.0000000E+00 -1.6900494E+01 | 0.0000000E+00 | -3.2618E+01 | + OH- 1.8416226E-06 -9.1483483E+01 -4.0506365E-01 -2.0957381E+01 |-0.0000000E+00 | -1.0883E+02 | ( 4.0165350E+00) + NaCl 3.9996395E-29 -1.0080997E+02 0.0000000E+00 -7.3706817E+01 | 0.0000000E+00 | -1.7452E+02 | + O2 8.9465959E-66 -2.4673669E+01 0.0000000E+00 -1.5809740E+02 | 0.0000000E+00 | -1.8277E+02 | +------------------------------------------------------------------------------------------------------------------- + + +Counters: Iterations Time (seconds) + vcs_basopt: 3 0.00000E+00 + vcs_TP: 36 1.66000E+00 +-------------------------------------------------------------------------------- +-------------------------------------------------------------------------------- + +TCounters: Num_Calls Total_Its Total_Time (seconds) + vcs_basopt: 3 3 0.00000E+00 + vcs_TP: 1 36 1.66000E+00 + vcs_inest: 0 0.00000E+00 + vcs_TotalTime: 1.68000E+00 Results from vcs: @@ -95,6 +315,8 @@ Moles: 2.32742 pressure 101325 Pa density 1216.41 kg/m^3 mean mol. weight 20.0596 amu + potential 0 V + pH 6.59 1 kg 1 kmol ----------- ------------ @@ -105,13 +327,14 @@ Moles: 2.32742 heat capacity c_p 3053.16 6.125e+04 J/K heat capacity c_v - X Y Chem. Pot. / RT - ------------- ------------ ------------ - H2O(L) 0.817565 0.734244 -124.003 - Cl- 0.0912175 0.161217 -72.3775 - H+ 7.91271e-10 3.97376e-11 -15.1739 - Na+ 0.0912175 0.104539 -102.139 - OH- 7.91271e-10 6.70891e-10 -108.829 + X Molalities Chem.Pot. ChemPotSS ActCoeffMolal + (J/kmol) (J/kmol) + ------------- ------------ ------------ ------------ ------------ + H2O(L) 0.817565 55.5084 -3.07399e+08 -3.06686e+08 0.917303 + Cl- 0.0912175 6.1932 -1.79421e+08 -1.83974e+08 1.01322 + H+ 7.91271e-10 5.37232e-08 -3.76154e+07 0 4.78537 + Na+ 0.0912175 6.1932 -2.53199e+08 -2.57752e+08 1.01322 + OH- 7.91271e-10 5.37232e-08 -2.69784e+08 -2.26784e+08 0.545262 *************** air ***************** Moles: 4.09718 diff --git a/test_problems/VCSnonideal/NaCl_equil/runtest b/test_problems/VCSnonideal/NaCl_equil/runtest index 75d1cafd1..0ce58caeb 100755 --- a/test_problems/VCSnonideal/NaCl_equil/runtest +++ b/test_problems/VCSnonideal/NaCl_equil/runtest @@ -17,7 +17,7 @@ testName=NaCl_equil ################################################################# MPEQUIL_EXE=${MPEQUIL_EXE:=nacl_equil} -$MPEQUIL_EXE -d 2 > out.txt 2>err_out.txt +$MPEQUIL_EXE -d 3 > out.txt 2>err_out.txt retnStat=$? if test $retnStat != "0" then