[Thermo] Fix PDSS_HKFT initialization when one property is not given

Calculating one of G0, H0, or S0 requires the parent ThermoPhase object, so this
calculation has to be delayed until m_tp has been set.
This commit is contained in:
Ray Speth 2017-02-25 18:11:13 -05:00
parent 51f419fbad
commit 2c3512c22a
2 changed files with 28 additions and 32 deletions

View file

@ -27,10 +27,10 @@ PDSS_HKFT::PDSS_HKFT()
, m_densWaterSS(-1.0)
, m_born_coeff_j(-1.0)
, m_r_e_j(-1.0)
, m_deltaG_formation_tr_pr(0.0)
, m_deltaH_formation_tr_pr(0.0)
, m_deltaG_formation_tr_pr(NAN)
, m_deltaH_formation_tr_pr(NAN)
, m_Mu0_tr_pr(0.0)
, m_Entrop_tr_pr(0.0)
, m_Entrop_tr_pr(NAN)
, m_a1(0.0)
, m_a2(0.0)
, m_a3(0.0)
@ -239,6 +239,28 @@ void PDSS_HKFT::initThermo()
{
PDSS::initThermo();
// Ok, if we are missing one, then we construct its value from the other two.
// This code has been internally verified.
m_charge_j = m_tp->charge(m_spindex);
if (std::isnan(m_deltaH_formation_tr_pr)) {
convertDGFormation();
doublereal Hcalc = m_Mu0_tr_pr + 298.15 * (m_Entrop_tr_pr * toSI("cal/gmol"));
m_deltaH_formation_tr_pr = Hcalc / toSI("cal/gmol");
} else if (std::isnan(m_deltaG_formation_tr_pr)) {
doublereal DHjmol = m_deltaH_formation_tr_pr * toSI("cal/gmol");
m_Mu0_tr_pr = DHjmol - 298.15 * (m_Entrop_tr_pr * toSI("cal/gmol"));
m_deltaG_formation_tr_pr = m_Mu0_tr_pr / toSI("cal/gmol");
double tmp = m_Mu0_tr_pr;
convertDGFormation();
double totalSum = m_Mu0_tr_pr - tmp;
m_Mu0_tr_pr = tmp;
m_deltaG_formation_tr_pr = (m_Mu0_tr_pr - totalSum)/ toSI("cal/gmol");
} else if (std::isnan(m_Entrop_tr_pr)) {
convertDGFormation();
doublereal DHjmol = m_deltaH_formation_tr_pr * toSI("cal/gmol");
m_Entrop_tr_pr = (DHjmol - m_Mu0_tr_pr) / (298.15 * toSI("cal/gmol"));
}
m_waterSS = &dynamic_cast<PDSS_Water&>(*m_tp->providePDSS(0));
// Section to initialize m_Z_pr_tr and m_Y_pr_tr
@ -253,7 +275,6 @@ void PDSS_HKFT::initThermo()
m_waterProps.reset(new WaterProps(m_waterSS));
m_presR_bar = OneAtm / 1.0E5;
m_presR_bar = 1.0;
m_charge_j = m_tp->charge(m_spindex);
convertDGFormation();
// Ok, we have mu. Let's check it against the input value
@ -408,31 +429,6 @@ void PDSS_HKFT::setParametersFromXML(const XML_Node& speciesNode)
"Missing 2 or more of DG0_f_Pr_Tr, DH0_f_Pr_Tr, or S0_f_Pr_Tr fields. "
"Need to supply at least two of these fields");
}
// Ok, if we are missing one, then we construct its value from the other two.
// This code has been internally verified.
if (hasDHO == 0) {
m_charge_j = m_tp->charge(m_spindex);
convertDGFormation();
doublereal Hcalc = m_Mu0_tr_pr + 298.15 * (m_Entrop_tr_pr * toSI("cal/gmol"));
m_deltaH_formation_tr_pr = Hcalc / toSI("cal/gmol");
}
if (hasDGO == 0) {
doublereal DHjmol = m_deltaH_formation_tr_pr * toSI("cal/gmol");
m_Mu0_tr_pr = DHjmol - 298.15 * (m_Entrop_tr_pr * toSI("cal/gmol"));
m_deltaG_formation_tr_pr = m_Mu0_tr_pr / toSI("cal/gmol");
double tmp = m_Mu0_tr_pr;
m_charge_j = m_tp->charge(m_spindex);
convertDGFormation();
double totalSum = m_Mu0_tr_pr - tmp;
m_Mu0_tr_pr = tmp;
m_deltaG_formation_tr_pr = (m_Mu0_tr_pr - totalSum)/ toSI("cal/gmol");
}
if (hasSO == 0) {
m_charge_j = m_tp->charge(m_spindex);
convertDGFormation();
doublereal DHjmol = m_deltaH_formation_tr_pr * toSI("cal/gmol");
m_Entrop_tr_pr = (DHjmol - m_Mu0_tr_pr) / (298.15 * toSI("cal/gmol"));
}
}
doublereal PDSS_HKFT::deltaH() const

View file

@ -133,7 +133,7 @@
<charge> +1 </charge>
<thermo model="HKFT">
<HKFT Pref="1 atm" Tmax=" 640." Tmin=" 273.15">
<DG0_f_Pr_Tr units="cal/gmol"> -62591. </DG0_f_Pr_Tr>
<!-- <DG0_f_Pr_Tr units="cal/gmol"> -62591. </DG0_f_Pr_Tr> -->
<DH0_f_Pr_Tr units="cal/gmol"> -57433. </DH0_f_Pr_Tr>
<S0_Pr_Tr units="cal/gmol/K"> 13.96 </S0_Pr_Tr>
</HKFT>
@ -161,7 +161,7 @@
<thermo model="HKFT">
<HKFT Pref="1 atm" Tmax=" 623.15" Tmin=" 298.00">
<DG0_f_Pr_Tr units="cal/gmol"> -31379. </DG0_f_Pr_Tr>
<DH0_f_Pr_Tr units="cal/gmol"> -39933. </DH0_f_Pr_Tr>
<!-- <DH0_f_Pr_Tr units="cal/gmol"> -39933. </DH0_f_Pr_Tr> -->
<S0_Pr_Tr units="cal/gmol/K"> 13.56 </S0_Pr_Tr>
</HKFT>
</thermo>
@ -189,7 +189,7 @@
<HKFT Pref="1 atm" Tmax=" 623.15" Tmin=" 298.00">
<DG0_f_Pr_Tr units="cal/gmol"> 0.0 </DG0_f_Pr_Tr>
<DH0_f_Pr_Tr units="cal/gmol"> 0.0 </DH0_f_Pr_Tr>
<S0_Pr_Tr units="cal/gmol/K"> 0.0 </S0_Pr_Tr>
<!-- <S0_Pr_Tr units="cal/gmol/K"> 0.0 </S0_Pr_Tr> -->
</HKFT>
</thermo>
<standardState model="HKFT">