From f56bce1757ceeaa3e881003db845c7207ab26726 Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Sat, 13 Dec 2008 01:59:49 +0000 Subject: [PATCH] Added more output from activity coefficient calculation in HMWSoln. Added a capability to readh in and change the H298f value in each species. Needed this for analysis. --- Cantera/src/thermo/ConstCpPoly.cpp | 19 +++ Cantera/src/thermo/ConstCpPoly.h | 8 + Cantera/src/thermo/GeneralSpeciesThermo.cpp | 22 +++ Cantera/src/thermo/GeneralSpeciesThermo.h | 8 + Cantera/src/thermo/HMWSoln.cpp | 151 ++++++++++++++++-- Cantera/src/thermo/IdealGasPhase.h | 7 + Cantera/src/thermo/LatticePhase.h | 9 +- Cantera/src/thermo/LatticeSolidPhase.h | 9 ++ Cantera/src/thermo/NasaPoly1.h | 36 +++++ Cantera/src/thermo/NasaPoly2.h | 30 ++++ Cantera/src/thermo/NasaThermo.h | 47 ++++++ Cantera/src/thermo/ShomatePoly.h | 67 ++++++++ Cantera/src/thermo/ShomateThermo.h | 50 ++++++ Cantera/src/thermo/SimpleThermo.h | 12 ++ Cantera/src/thermo/SingleSpeciesTP.h | 9 +- Cantera/src/thermo/SpeciesThermo.h | 8 + .../src/thermo/SpeciesThermoInterpType.cpp | 14 ++ Cantera/src/thermo/SpeciesThermoInterpType.h | 7 + Cantera/src/thermo/SpeciesThermoMgr.h | 24 +++ Cantera/src/thermo/StoichSubstance.h | 7 + Cantera/src/thermo/SurfPhase.h | 8 + Cantera/src/thermo/ThermoPhase.h | 9 ++ Cantera/src/thermo/VPStandardStateTP.h | 9 +- 23 files changed, 557 insertions(+), 13 deletions(-) diff --git a/Cantera/src/thermo/ConstCpPoly.cpp b/Cantera/src/thermo/ConstCpPoly.cpp index e72f56681..c6406c7f9 100644 --- a/Cantera/src/thermo/ConstCpPoly.cpp +++ b/Cantera/src/thermo/ConstCpPoly.cpp @@ -137,6 +137,25 @@ namespace Cantera { m_logt0 = log(m_t0); } +#ifdef H298MODIFY_CAPABILITY + + doublereal ConstCpPoly::reportHf298(doublereal* const h298) const { + double temp = 298.15; + doublereal h = GasConstant * (m_h0_R + (temp - m_t0) * m_cp0_R); + if (h298) { + h298[m_index] = h; + } + return h; + } + + void ConstCpPoly::modifyOneHf298(const int k, const doublereal Hf298New) { + if (k != m_index) return; + doublereal hnow = reportHf298(); + doublereal delH = Hf298New - hnow; + m_h0_R += delH / GasConstant; + } + +#endif } diff --git a/Cantera/src/thermo/ConstCpPoly.h b/Cantera/src/thermo/ConstCpPoly.h index 8f09eed03..7577769c5 100644 --- a/Cantera/src/thermo/ConstCpPoly.h +++ b/Cantera/src/thermo/ConstCpPoly.h @@ -170,6 +170,14 @@ namespace Cantera { */ virtual void modifyParameters(doublereal* coeffs); +#ifdef H298MODIFY_CAPABILITY + + virtual doublereal reportHf298(doublereal* const h298 = 0) const; + + virtual void modifyOneHf298(const int k, const doublereal Hf298New); + +#endif + protected: //! Base temperature doublereal m_t0; diff --git a/Cantera/src/thermo/GeneralSpeciesThermo.cpp b/Cantera/src/thermo/GeneralSpeciesThermo.cpp index e8eb969aa..b96fd7bbe 100644 --- a/Cantera/src/thermo/GeneralSpeciesThermo.cpp +++ b/Cantera/src/thermo/GeneralSpeciesThermo.cpp @@ -352,4 +352,26 @@ namespace Cantera { return (m_sp[k]); } +#ifdef H298MODIFY_CAPABILITY + + doublereal GeneralSpeciesThermo::reportOneHf298(int k) const { + SpeciesThermoInterpType * sp_ptr = m_sp[k]; + doublereal h = -1.0; + if (sp_ptr) { + h = sp_ptr->reportHf298(0); + } + return h; + } + + void GeneralSpeciesThermo::modifyOneHf298(const int k, const doublereal Hf298New) { + SpeciesThermoInterpType * sp_ptr = m_sp[k]; + if (sp_ptr) { + sp_ptr->modifyOneHf298(k, Hf298New); + } + } + + +#endif + + } diff --git a/Cantera/src/thermo/GeneralSpeciesThermo.h b/Cantera/src/thermo/GeneralSpeciesThermo.h index 57f06c76d..0697d7e9d 100644 --- a/Cantera/src/thermo/GeneralSpeciesThermo.h +++ b/Cantera/src/thermo/GeneralSpeciesThermo.h @@ -218,6 +218,14 @@ namespace Cantera { */ virtual void modifyParams(int index, doublereal *c); +#ifdef H298MODIFY_CAPABILITY + + virtual doublereal reportOneHf298(int k) const; + + virtual void modifyOneHf298(const int k, const doublereal Hf298New); + +#endif + private: //! Provide the SpeciesthermoInterpType object /*! diff --git a/Cantera/src/thermo/HMWSoln.cpp b/Cantera/src/thermo/HMWSoln.cpp index d5a57c085..25a6c27a7 100644 --- a/Cantera/src/thermo/HMWSoln.cpp +++ b/Cantera/src/thermo/HMWSoln.cpp @@ -2605,7 +2605,7 @@ namespace Cantera { } #ifdef DEBUG_MODE if (m_debugCalc) { - printf(" Step 8: \n"); + printf(" Step 8: Summing in All Contributions to Activity Coefficients \n"); } #endif @@ -2616,9 +2616,21 @@ namespace Cantera { * -------- -> equations agree with my notes, Eqn. (118). * -> Equations agree with Pitzer, eqn.(63) */ - if (charge[i] > 0 ) { + if (charge[i] > 0.0 ) { + +#ifdef DEBUG_MODE + if (m_debugCalc) { + sni = speciesName(i); + printf(" Contributions to ln(ActCoeff_%s):\n", sni.c_str()); + } +#endif // species i is the cation (positive) to calc the actcoeff zsqF = charge[i]*charge[i]*F; +#ifdef DEBUG_MODE + if (m_debugCalc) { + printf(" Unary term: z*z*F = %10.5f\n", zsqF); + } +#endif sum1 = 0.0; sum2 = 0.0; sum3 = 0.0; @@ -2634,7 +2646,16 @@ namespace Cantera { if (charge[j] < 0.0) { // sum over all anions sum1 = sum1 + molality[j]* - (2.0*BMX[counterIJ]+molarcharge*CMX[counterIJ]); + (2.0*BMX[counterIJ] + molarcharge*CMX[counterIJ]); +#ifdef DEBUG_MODE + if (m_debugCalc) { + snj = speciesName(j) + ":"; + printf(" Bin term with %-13s 2 m_j BMX = %10.5f\n", snj.c_str(), + molality[j]*2.0*BMX[counterIJ]); + printf(" m_j Z CMX = %10.5f\n", + molality[j]* molarcharge*CMX[counterIJ]); + } +#endif if (j < m_kk-1) { /* * This term is the ternary interaction involving the @@ -2646,6 +2667,15 @@ namespace Cantera { if (charge[k] < 0.0) { n = k + j * m_kk + i * m_kk * m_kk; sum3 = sum3 + molality[j]*molality[k]*psi_ijk[n]; +#ifdef DEBUG_MODE + if (m_debugCalc) { + if (psi_ijk[n] != 0.0) { + snj = speciesName(j) + "," + speciesName(k) + ":"; + printf(" Psi term on %-16s m_j m_k psi_ijk = %10.5f\n", snj.c_str(), + molality[j]*molality[k]*psi_ijk[n]); + } + } +#endif } } } @@ -2654,13 +2684,33 @@ namespace Cantera { if (charge[j] > 0.0) { // sum over all cations - if (j != i) sum2 = sum2 + molality[j]*(2.0*Phi[counterIJ]); + if (j != i) { + sum2 = sum2 + molality[j]*(2.0*Phi[counterIJ]); +#ifdef DEBUG_MODE + if (m_debugCalc) { + if ((molality[j] * Phi[counterIJ])!= 0.0) { + snj = speciesName(j) + ":"; + printf(" Phi term with %-12s 2 m_j Phi_cc = %10.5f\n", snj.c_str(), + molality[j]*(2.0*Phi[counterIJ])); + } + } +#endif + } for (k = 1; k < m_kk; k++) { if (charge[k] < 0.0) { // two inner sums over anions n = k + j * m_kk + i * m_kk * m_kk; sum2 = sum2 + molality[j]*molality[k]*psi_ijk[n]; +#ifdef DEBUG_MODE + if (m_debugCalc) { + if (psi_ijk[n] != 0.0) { + snj = speciesName(j) + "," + speciesName(k) + ":"; + printf(" Psi term on %-16s m_j m_k psi_ijk = %10.5f\n", snj.c_str(), + molality[j]*molality[k]*psi_ijk[n]); + } + } +#endif /* * Find the counterIJ for the j,k interaction */ @@ -2668,6 +2718,15 @@ namespace Cantera { counterIJ2 = m_CounterIJ[n]; sum4 = sum4 + (fabs(charge[i])* molality[j]*molality[k]*CMX[counterIJ2]); +#ifdef DEBUG_MODE + if (m_debugCalc) { + if ((molality[j]*molality[k]*CMX[counterIJ2]) != 0.0) { + snj = speciesName(j) + "," + speciesName(k) + ":"; + printf(" Tern CMX term on %-16s abs(z_i) m_j m_k CMX = %10.5f\n", snj.c_str(), + fabs(charge[i])* molality[j]*molality[k]*CMX[counterIJ2]); + } + } +#endif } } } @@ -2677,6 +2736,15 @@ namespace Cantera { */ if (charge[j] == 0) { sum5 = sum5 + molality[j]*2.0*m_Lambda_nj(j,i); +#ifdef DEBUG_MODE + if (m_debugCalc) { + if ((molality[j]*2.0*m_Lambda_nj(j,i)) != 0.0) { + snj = speciesName(j) + ":"; + printf(" Lambda term with %-12s 2 m_j lam_ji = %10.5f\n", snj.c_str(), + molality[j]*2.0*m_Lambda_nj(j,i)); + } + } +#endif } } /* @@ -2688,10 +2756,8 @@ namespace Cantera { #ifdef DEBUG_MODE if (m_debugCalc) { sni = speciesName(i); - printf(" %-16s lngamma[i]=%10.6f gamma[i]=%10.6f \n", + printf(" Net %-16s lngamma[i] = %9.5f gamma[i]=%10.6f \n", sni.c_str(), m_lnActCoeffMolal[i], gamma[i]); - printf(" %12g %12g %12g %12g %12g %12g\n", - zsqF, sum1, sum2, sum3, sum4, sum5); } #endif } @@ -2702,8 +2768,21 @@ namespace Cantera { * -> Equations agree with Pitzer, eqn.(64) */ if (charge[i] < 0 ) { + +#ifdef DEBUG_MODE + if (m_debugCalc) { + sni = speciesName(i); + printf(" Contributions to ln(ActCoeff_%s):\n", sni.c_str()); + } +#endif + // species i is an anion (negative) zsqF = charge[i]*charge[i]*F; +#ifdef DEBUG_MODE + if (m_debugCalc) { + printf(" Unary term: z*z*F = %10.5f\n", zsqF); + } +#endif sum1 = 0.0; sum2 = 0.0; sum3 = 0.0; @@ -2722,12 +2801,30 @@ namespace Cantera { if (charge[j] > 0) { sum1 = sum1 + molality[j]* (2.0*BMX[counterIJ]+molarcharge*CMX[counterIJ]); +#ifdef DEBUG_MODE + if (m_debugCalc) { + snj = speciesName(j) + ":"; + printf(" Bin term with %-13s 2 m_j BMX = %10.5f\n", snj.c_str(), + molality[j]*2.0*BMX[counterIJ]); + printf(" m_j Z CMX = %10.5f\n", + molality[j]* molarcharge*CMX[counterIJ]); + } +#endif if (j < m_kk-1) { for (k = j+1; k < m_kk; k++) { // an inner sum over all cations if (charge[k] > 0) { n = k + j * m_kk + i * m_kk * m_kk; sum3 = sum3 + molality[j]*molality[k]*psi_ijk[n]; +#ifdef DEBUG_MODE + if (m_debugCalc) { + if (psi_ijk[n] != 0.0) { + snj = speciesName(j) + "," + speciesName(k) + ":"; + printf(" Psi term on %-16s m_j m_k psi_ijk = %10.5f\n", snj.c_str(), + molality[j]*molality[k]*psi_ijk[n]); + } + } +#endif } } } @@ -2740,12 +2837,30 @@ namespace Cantera { // sum over all anions if (j != i) { sum2 = sum2 + molality[j]*(2.0*Phi[counterIJ]); +#ifdef DEBUG_MODE + if (m_debugCalc) { + if ((molality[j] * Phi[counterIJ])!= 0.0) { + snj = speciesName(j) + ":"; + printf(" Phi term with %-12s 2 m_j Phi_aa = %10.5f\n", snj.c_str(), + molality[j]*(2.0*Phi[counterIJ])); + } + } +#endif } for (k = 1; k < m_kk; k++) { if (charge[k] > 0.0) { // two inner sums over cations n = k + j * m_kk + i * m_kk * m_kk; sum2 = sum2 + molality[j]*molality[k]*psi_ijk[n]; +#ifdef DEBUG_MODE + if (m_debugCalc) { + if (psi_ijk[n] != 0.0) { + snj = speciesName(j) + "," + speciesName(k) + ":"; + printf(" Psi term on %-16s m_j m_k psi_ijk = %10.5f\n", snj.c_str(), + molality[j]*molality[k]*psi_ijk[n]); + } + } +#endif /* * Find the counterIJ for the symmetric binary interaction */ @@ -2754,6 +2869,15 @@ namespace Cantera { sum4 = sum4 + (fabs(charge[i])* molality[j]*molality[k]*CMX[counterIJ2]); +#ifdef DEBUG_MODE + if (m_debugCalc) { + if ((molality[j]*molality[k]*CMX[counterIJ2]) != 0.0) { + snj = speciesName(j) + "," + speciesName(k) + ":"; + printf(" Tern CMX term on %-16s abs(z_i) m_j m_k CMX = %10.5f\n", snj.c_str(), + fabs(charge[i])* molality[j]*molality[k]*CMX[counterIJ2]); + } + } +#endif } } } @@ -2763,6 +2887,15 @@ namespace Cantera { */ if (charge[j] == 0.0) { sum5 = sum5 + molality[j]*2.0*m_Lambda_nj(j,i); +#ifdef DEBUG_MODE + if (m_debugCalc) { + if ((molality[j]*2.0*m_Lambda_nj(j,i)) != 0.0) { + snj = speciesName(j) + ":"; + printf(" Lambda term with %-12s 2 m_j lam_ji = %10.5f\n", snj.c_str(), + molality[j]*2.0*m_Lambda_nj(j,i)); + } + } +#endif } } m_lnActCoeffMolal[i] = zsqF + sum1 + sum2 + sum3 + sum4 + sum5; @@ -2770,10 +2903,8 @@ namespace Cantera { #ifdef DEBUG_MODE if (m_debugCalc) { sni = speciesName(i); - printf(" %-16s lngamma[i]=%10.6f gamma[i]=%10.6f\n", + printf(" Net %-16s lngamma[i] = %9.5f gamma[i]=%10.6f\n", sni.c_str(), m_lnActCoeffMolal[i], gamma[i]); - printf(" %12g %12g %12g %12g %12g %12g\n", - zsqF, sum1, sum2, sum3, sum4, sum5); } #endif } diff --git a/Cantera/src/thermo/IdealGasPhase.h b/Cantera/src/thermo/IdealGasPhase.h index 804049878..eadd648bd 100644 --- a/Cantera/src/thermo/IdealGasPhase.h +++ b/Cantera/src/thermo/IdealGasPhase.h @@ -712,6 +712,13 @@ namespace Cantera { */ virtual void getEnthalpy_RT_ref(doublereal *hrt) const; +#ifdef H298MODIFY_CAPABILITY + + virtual void modifyOneHf298SS(const int k, const doublereal Hf298New) { + m_spthermo->modifyOneHf298(k, Hf298New); + m_tlast += 0.0001234; + } +#endif //! Returns the vector of nondimensional //! Gibbs Free Energies of the reference state at the current temperature //! of the solution and the reference pressure for the species. diff --git a/Cantera/src/thermo/LatticePhase.h b/Cantera/src/thermo/LatticePhase.h index c6d2b890d..7c160560c 100644 --- a/Cantera/src/thermo/LatticePhase.h +++ b/Cantera/src/thermo/LatticePhase.h @@ -632,7 +632,14 @@ namespace Cantera { /// @name Thermodynamic Values for the Species Reference States //@{ - +#ifdef H298MODIFY_CAPABILITY + + virtual void modifyOneHf298SS(const int k, const doublereal Hf298New) { + m_spthermo->modifyOneHf298(k, Hf298New); + m_tlast += 0.0001234; + } +#endif + //! Returns the vector of nondimensional //! Enthalpies of the reference state at the current temperature //! of the solution and the reference pressure for the phase. diff --git a/Cantera/src/thermo/LatticeSolidPhase.h b/Cantera/src/thermo/LatticeSolidPhase.h index a08585043..d97116f3b 100644 --- a/Cantera/src/thermo/LatticeSolidPhase.h +++ b/Cantera/src/thermo/LatticeSolidPhase.h @@ -128,6 +128,15 @@ namespace Cantera { void setLatticeMoleFractions(int n, std::string x); + +#ifdef H298MODIFY_CAPABILITY + + virtual void modifyOneHf298SS(const int k, const doublereal Hf298New) { + m_spthermo->modifyOneHf298(k, Hf298New); + m_tlast += 0.0001234; + } +#endif + protected: int m_mm; diff --git a/Cantera/src/thermo/NasaPoly1.h b/Cantera/src/thermo/NasaPoly1.h index 976ce9629..60736ecca 100755 --- a/Cantera/src/thermo/NasaPoly1.h +++ b/Cantera/src/thermo/NasaPoly1.h @@ -269,6 +269,42 @@ namespace Cantera { } } +#ifdef H298MODIFY_CAPABILITY + + virtual doublereal reportHf298(doublereal* const h298 = 0) const { + double tt[6]; + double temp = 298.15; + tt[0] = temp; + tt[1] = temp * temp; + tt[2] = tt[1] * temp; + tt[3] = tt[2] * temp; + tt[4] = 1.0 / temp; + //tt[5] = std::log(temp); + doublereal ct0 = m_coeff[2]; // a0 + doublereal ct1 = m_coeff[3]*tt[0]; // a1 * T + doublereal ct2 = m_coeff[4]*tt[1]; // a2 * T^2 + doublereal ct3 = m_coeff[5]*tt[2]; // a3 * T^3 + doublereal ct4 = m_coeff[6]*tt[3]; // a4 * T^4 + + double h_RT = ct0 + 0.5*ct1 + OneThird*ct2 + 0.25*ct3 + 0.2*ct4 + + m_coeff[0]*tt[4]; // last t + + double h = h_RT * GasConstant * temp; + if (h298) { + h298[m_index] = h; + } + return h; + } + + virtual void modifyOneHf298(const int k, const doublereal Hf298New) { + if (k != m_index) return; + double hcurr = reportHf298(0); + double delH = Hf298New - hcurr; + m_coeff[0] += (delH) / GasConstant; + } + +#endif + protected: //! lowest valid temperature doublereal m_lowT; diff --git a/Cantera/src/thermo/NasaPoly2.h b/Cantera/src/thermo/NasaPoly2.h index 6456b236a..3e885f561 100644 --- a/Cantera/src/thermo/NasaPoly2.h +++ b/Cantera/src/thermo/NasaPoly2.h @@ -263,6 +263,36 @@ namespace Cantera { } } +#ifdef H298MODIFY_CAPABILITY + + doublereal reportHf298(doublereal* const h298 = 0) const { + double h; + if (298.15 <= m_midT) { + h = mnp_low->reportHf298(0); + } else { + h = mnp_high->reportHf298(0); + } + if (h298) { + h298[m_index] = h; + } + return h; + } + + void modifyOneHf298(const int k, const doublereal Hf298New) { + if (k != m_index) return; + + doublereal h298now = reportHf298(0); + doublereal delH = Hf298New - h298now; + double h = mnp_low->reportHf298(0); + double hnew = h + delH; + mnp_low->modifyOneHf298(k, hnew); + h = mnp_high->reportHf298(0); + hnew = h + delH; + mnp_high->modifyOneHf298(k, hnew); + } + +#endif + protected: //! lowest valid temperature doublereal m_lowT; diff --git a/Cantera/src/thermo/NasaThermo.h b/Cantera/src/thermo/NasaThermo.h index 4eca544ed..704ee5730 100755 --- a/Cantera/src/thermo/NasaThermo.h +++ b/Cantera/src/thermo/NasaThermo.h @@ -455,6 +455,51 @@ namespace Cantera { } } +#ifdef H298MODIFY_CAPABILITY + virtual doublereal reportOneHf298(const int k) const { + + int grp = m_group_map[k]; + int pos = m_posInGroup_map[k]; + const vector &mlg = m_low[grp-1]; + const NasaPoly1 *nlow = &(mlg[pos]); + doublereal tmid = nlow->maxTemp(); + double h; + if (298.15 <= tmid) { + h = nlow->reportHf298(0); + } else { + const vector &mhg = m_high[grp-1]; + const NasaPoly1 *nhigh = &(mhg[pos]); + h = nhigh->reportHf298(0); + } + return h; + } + + virtual void modifyOneHf298(const int k, const doublereal Hf298New) { + int grp = m_group_map[k]; + int pos = m_posInGroup_map[k]; + vector &mlg = m_low[grp-1]; + NasaPoly1 *nlow = &(mlg[pos]); + vector &mhg = m_high[grp-1]; + NasaPoly1 *nhigh = &(mhg[pos]); + doublereal tmid = nlow->maxTemp(); + + double hnow = reportOneHf298(k); + double delH = Hf298New - hnow; + if (298.15 <= tmid) { + nlow->modifyOneHf298(k, Hf298New); + double h = nhigh->reportHf298(0); + double hnew = h + delH; + nhigh->modifyOneHf298(k, hnew); + } else { + nhigh->modifyOneHf298(k, Hf298New); + double h = nlow->reportHf298(0); + double hnew = h + delH; + nlow->modifyOneHf298(k, hnew); + } + + } +#endif + protected: //! Vector of vector of NasaPoly1's for the high temp region. /*! @@ -568,6 +613,8 @@ namespace Cantera { + c[6]; } + + }; } diff --git a/Cantera/src/thermo/ShomatePoly.h b/Cantera/src/thermo/ShomatePoly.h index d496d07fe..2937d8f0b 100755 --- a/Cantera/src/thermo/ShomatePoly.h +++ b/Cantera/src/thermo/ShomatePoly.h @@ -289,6 +289,43 @@ namespace Cantera { std::copy(coeffs, coeffs + 7, m_coeff.begin()); } +#ifdef H298MODIFY_CAPABILITY + + + virtual doublereal reportHf298(doublereal* const h298 = 0) const { + + double tPoly[4]; + doublereal tt = 1.e-3*298.15; + tPoly[0] = tt; + tPoly[1] = tt * tt; + tPoly[2] = tPoly[1] * tt; + tPoly[3] = 1.0/tPoly[1]; + + doublereal A = m_coeff[0]; + doublereal Bt = m_coeff[1]*tPoly[0]; + doublereal Ct2 = m_coeff[2]*tPoly[1]; + doublereal Dt3 = m_coeff[3]*tPoly[2]; + doublereal Etm2 = m_coeff[4]*tPoly[3]; + doublereal F = m_coeff[5]; + + doublereal h = tPoly[0]*(A + 0.5*Bt + OneThird*Ct2 + 0.25*Dt3 - Etm2) + F; + + double hh = 1.e6 * h; + if (h298) { + h298[m_index] = 1.e6 * h; + } + return hh; + } + + + virtual void modifyOneHf298(const int k, const doublereal Hf298New) { + doublereal hnow = reportHf298(); + doublereal delH = Hf298New - hnow; + m_coeff[5] += delH / 1.0E6; + } + +#endif + protected: //! Minimum temperature for which the parameterization is valid (Kelvin) doublereal m_lowT; @@ -576,6 +613,36 @@ namespace Cantera { msp_high = new ShomatePoly(m_index, m_midT, m_highT, m_Pref, coeffs+8); } +#ifdef H298MODIFY_CAPABILITY + + virtual doublereal reportHf298(doublereal* const h298 = 0) const { + doublereal h; + if (298.15 <= m_midT) { + h = msp_low->reportHf298(h298); + } else { + h = msp_high->reportHf298(h298); + } + if (h298) { + h298[m_index] = h; + } + return h; + } + + virtual void modifyOneHf298(const int k, const doublereal Hf298New) { + if (k != m_index) return; + + doublereal h298now = reportHf298(0); + doublereal delH = Hf298New - h298now; + double h = msp_low->reportHf298(0); + double hnew = h + delH; + msp_low->modifyOneHf298(k, hnew); + h = msp_high->reportHf298(0); + hnew = h + delH; + msp_high->modifyOneHf298(k, hnew); + } + +#endif + protected: //! Minimum temperature the representation is valid(kelvin) doublereal m_lowT; diff --git a/Cantera/src/thermo/ShomateThermo.h b/Cantera/src/thermo/ShomateThermo.h index bc9d31957..c37822801 100755 --- a/Cantera/src/thermo/ShomateThermo.h +++ b/Cantera/src/thermo/ShomateThermo.h @@ -445,6 +445,56 @@ namespace Cantera { } } +#ifdef H298MODIFY_CAPABILITY + + virtual doublereal reportOneHf298(int k) const { + doublereal h; + doublereal t = 298.15; + + int grp = m_group_map[k]; + int pos = m_posInGroup_map[k]; + const vector &mlg = m_low[grp-1]; + const ShomatePoly *nlow = &(mlg[pos]); + + doublereal tmid = nlow->maxTemp(); + if (t <= tmid) { + h = nlow->reportHf298(); + } else { + const vector &mhg = m_high[grp-1]; + const ShomatePoly *nhigh = &(mhg[pos]); + h = nhigh->reportHf298(); + } + return h; + } + + virtual void modifyOneHf298(const int k, const doublereal Hf298New) { + + int grp = m_group_map[k]; + int pos = m_posInGroup_map[k]; + vector &mlg = m_low[grp-1]; + ShomatePoly *nlow = &(mlg[pos]); + vector &mhg = m_high[grp-1]; + ShomatePoly *nhigh = &(mhg[pos]); + doublereal tmid = nlow->maxTemp(); + + double hnow = reportOneHf298(k); + double delH = Hf298New - hnow; + if (298.15 <= tmid) { + nlow->modifyOneHf298(k, Hf298New); + double h = nhigh->reportHf298(0); + double hnew = h + delH; + nhigh->modifyOneHf298(k, hnew); + } else { + nhigh->modifyOneHf298(k, Hf298New); + double h = nlow->reportHf298(0); + double hnew = h + delH; + nlow->modifyOneHf298(k, hnew); + } + + } + + +#endif protected: //! Vector of vector of NasaPoly1's for the high temp region. diff --git a/Cantera/src/thermo/SimpleThermo.h b/Cantera/src/thermo/SimpleThermo.h index 406f1157e..2a548c739 100644 --- a/Cantera/src/thermo/SimpleThermo.h +++ b/Cantera/src/thermo/SimpleThermo.h @@ -361,6 +361,18 @@ namespace Cantera { m_cp0_R[loc] = c[3] / GasConstant; } +#ifdef H298MODIFY_CAPABILITY + + virtual doublereal reportOneHf298(int k) const { + throw CanteraError("reportHF298", "unimplemented"); + } + + virtual void modifyOneHf298(const int k, const doublereal Hf298New) { + throw CanteraError("reportHF298", "unimplemented"); + } + + +#endif protected: //! Mapping between the species index and the vector index where the coefficients are kept diff --git a/Cantera/src/thermo/SingleSpeciesTP.h b/Cantera/src/thermo/SingleSpeciesTP.h index 3eb58c95b..59e726bf2 100644 --- a/Cantera/src/thermo/SingleSpeciesTP.h +++ b/Cantera/src/thermo/SingleSpeciesTP.h @@ -393,7 +393,14 @@ namespace Cantera { /// equation of state. //@{ - + #ifdef H298MODIFY_CAPABILITY + + virtual void modifyOneHf298SS(const int k, const doublereal Hf298New) { + m_spthermo->modifyOneHf298(k, Hf298New); + m_tlast += 0.0001234; + } +#endif + /*! * Returns the vector of nondimensional * enthalpies of the reference state at the current temperature diff --git a/Cantera/src/thermo/SpeciesThermo.h b/Cantera/src/thermo/SpeciesThermo.h index b30e326a9..35b8c8edb 100755 --- a/Cantera/src/thermo/SpeciesThermo.h +++ b/Cantera/src/thermo/SpeciesThermo.h @@ -367,7 +367,15 @@ namespace Cantera { * parameters for the standard state. */ virtual void modifyParams(int index, doublereal *c) = 0; + +#ifdef H298MODIFY_CAPABILITY + virtual doublereal reportOneHf298(int k) const = 0; + + virtual void modifyOneHf298(const int k, const doublereal Hf298New) = 0; + + +#endif }; //@} } diff --git a/Cantera/src/thermo/SpeciesThermoInterpType.cpp b/Cantera/src/thermo/SpeciesThermoInterpType.cpp index e448f20ff..da2ecf503 100644 --- a/Cantera/src/thermo/SpeciesThermoInterpType.cpp +++ b/Cantera/src/thermo/SpeciesThermoInterpType.cpp @@ -29,6 +29,20 @@ namespace Cantera { updatePropertiesTemp(T, cp_R, h_RT, s_R); } +#ifdef H298MODIFY_CAPABILITY + + doublereal SpeciesThermoInterpType::reportHf298(doublereal* const h298) const { + throw CanteraError("SpeciesThermoInterpType::reportHf298", + "Not implemented"); + } + + void SpeciesThermoInterpType::modifyOneHf298(const int k, const doublereal Hf298New) { + throw CanteraError("SpeciesThermoInterpType::modifyOneHf298", + "Not implemented"); + } + +#endif + /***************************************************************************************************/ //! Constructor diff --git a/Cantera/src/thermo/SpeciesThermoInterpType.h b/Cantera/src/thermo/SpeciesThermoInterpType.h index 3a460d2e6..6e8f83b04 100644 --- a/Cantera/src/thermo/SpeciesThermoInterpType.h +++ b/Cantera/src/thermo/SpeciesThermoInterpType.h @@ -264,6 +264,13 @@ namespace Cantera { */ virtual void modifyParameters(doublereal* coeffs) {} +#ifdef H298MODIFY_CAPABILITY + + virtual doublereal reportHf298(doublereal* const h298 = 0) const; + + virtual void modifyOneHf298(const int k, const doublereal Hf298New); + +#endif }; //! Class for the thermoydnamic manager for an individual species' reference state diff --git a/Cantera/src/thermo/SpeciesThermoMgr.h b/Cantera/src/thermo/SpeciesThermoMgr.h index af982c7c9..0dc405c9b 100755 --- a/Cantera/src/thermo/SpeciesThermoMgr.h +++ b/Cantera/src/thermo/SpeciesThermoMgr.h @@ -316,6 +316,19 @@ namespace Cantera { */ virtual void modifyParams(int index, doublereal *c); + +#ifdef H298MODIFY_CAPABILITY + + virtual doublereal reportOneHf298(int k) const { + throw CanteraError("reportHF298", "unimplemented"); + } + + virtual void modifyOneHf298(const int k, const doublereal Hf298New) { + throw CanteraError("reportHF298", "unimplemented"); + } + +#endif + private: //! Thermo Type 1 @@ -465,6 +478,17 @@ namespace Cantera { */ virtual void modifyParams(int index, doublereal *c); +#ifdef H298MODIFY_CAPABILITY + + virtual doublereal reportOneHf298(int k) const { + throw CanteraError("reportHF298", "unimplemented"); + } + + virtual void modifyOneHf298(const int k, const doublereal Hf298New) { + throw CanteraError("reportHF298", "unimplemented"); + } + +#endif private: //! Vector of SPM objects. There are m_kk of them std::vector m_thermo; diff --git a/Cantera/src/thermo/StoichSubstance.h b/Cantera/src/thermo/StoichSubstance.h index b8d78065b..022cf5e67 100644 --- a/Cantera/src/thermo/StoichSubstance.h +++ b/Cantera/src/thermo/StoichSubstance.h @@ -329,6 +329,13 @@ namespace Cantera { */ virtual void getEnthalpy_RT_ref(doublereal *hrt) const; +#ifdef H298MODIFY_CAPABILITY + + virtual void modifyOneHf298SS(const int k, const doublereal Hf298New) { + m_spthermo->modifyOneHf298(k, Hf298New); + m_tlast += 0.0001234; + } +#endif /** * Returns the vector of nondimensional * enthalpies of the reference state at the current temperature diff --git a/Cantera/src/thermo/SurfPhase.h b/Cantera/src/thermo/SurfPhase.h index 4f90f431b..29a0e9f30 100644 --- a/Cantera/src/thermo/SurfPhase.h +++ b/Cantera/src/thermo/SurfPhase.h @@ -546,6 +546,14 @@ namespace Cantera { */ virtual void getEnthalpy_RT_ref(doublereal* hrt) const; +#ifdef H298MODIFY_CAPABILITY + + virtual void modifyOneHf298SS(const int k, const doublereal Hf298New) { + m_spthermo->modifyOneHf298(k, Hf298New); + m_tlast += 0.0001234; + } +#endif + //! Returns the vector of nondimensional //! entropies of the reference state at the current temperature //! of the solution and the reference pressure for each species. diff --git a/Cantera/src/thermo/ThermoPhase.h b/Cantera/src/thermo/ThermoPhase.h index c868652a5..40d0025ae 100755 --- a/Cantera/src/thermo/ThermoPhase.h +++ b/Cantera/src/thermo/ThermoPhase.h @@ -780,6 +780,15 @@ namespace Cantera { return m_spthermo->minTemp(k); } +#ifdef H298MODIFY_CAPABILITY + doublereal Hf298SS(const int k) const { + return (m_spthermo->reportOneHf298(k)); + } + + virtual void modifyOneHf298SS(const int k, const doublereal Hf298New) { + m_spthermo->modifyOneHf298(k, Hf298New); + } +#endif //! Maximum temperature for which the thermodynamic data for the species //! are valid. /*! diff --git a/Cantera/src/thermo/VPStandardStateTP.h b/Cantera/src/thermo/VPStandardStateTP.h index 95204725c..0e09a2925 100644 --- a/Cantera/src/thermo/VPStandardStateTP.h +++ b/Cantera/src/thermo/VPStandardStateTP.h @@ -353,7 +353,14 @@ namespace Cantera { * length = m_kk, units = dimensionless. */ virtual void getEnthalpy_RT_ref(doublereal *hrt) const; - + +#ifdef H298MODIFY_CAPABILITY + + void modifyOneHf298SS(const int k, const doublereal Hf298New) { + m_spthermo->modifyOneHf298(k, Hf298New); + m_Tlast_ss += 0.0001234; + } +#endif /*! * Returns the vector of nondimensional * Gibbs free energies of the reference state at the current temperature