Changed a setState function to a setState_TR function for consistency
with the rest of Cantera Added an analytical derivative for dpdT for the water object. Started adding states for the water object that refer to unstable conditions within the spinodal curve.
This commit is contained in:
parent
ebe950fc4c
commit
0f57869240
7 changed files with 88 additions and 51 deletions
|
|
@ -346,7 +346,7 @@ namespace Cantera {
|
|||
doublereal T = m_temp;
|
||||
doublereal dens0 = m_sub->density(T, m_p0);
|
||||
doublereal h = m_sub->enthalpy(T, dens0);
|
||||
m_sub->setState(m_temp, m_dens);
|
||||
m_sub->setState_TR(m_temp, m_dens);
|
||||
return ((h + EW_Offset - SW_Offset*T)/(T * GasConstant));
|
||||
}
|
||||
|
||||
|
|
@ -356,7 +356,7 @@ namespace Cantera {
|
|||
doublereal T = m_temp;
|
||||
doublereal dens0 = m_sub->density(T, m_p0);
|
||||
doublereal h = m_sub->enthalpy(T, dens0);
|
||||
m_sub->setState(m_temp, m_dens);
|
||||
m_sub->setState_TR(m_temp, m_dens);
|
||||
return ((h + EW_Offset)/(T * GasConstant));
|
||||
}
|
||||
|
||||
|
|
@ -365,7 +365,7 @@ namespace Cantera {
|
|||
doublereal T = m_temp;
|
||||
doublereal dens0 = m_sub->density(T, m_p0);
|
||||
doublereal s = m_sub->entropy(T, dens0);
|
||||
m_sub->setState(m_temp, m_dens);
|
||||
m_sub->setState_TR(m_temp, m_dens);
|
||||
return ((s + SW_Offset)/GasConstant);
|
||||
}
|
||||
|
||||
|
|
@ -374,7 +374,7 @@ namespace Cantera {
|
|||
doublereal T = m_temp;
|
||||
doublereal dens0 = m_sub->density(T, m_p0);
|
||||
doublereal cp = m_sub->cp(T, dens0);
|
||||
m_sub->setState(m_temp, m_dens);
|
||||
m_sub->setState_TR(m_temp, m_dens);
|
||||
return (cp/GasConstant);
|
||||
}
|
||||
|
||||
|
|
@ -383,7 +383,7 @@ namespace Cantera {
|
|||
doublereal T = m_temp;
|
||||
doublereal dens0 = m_sub->density(T, m_p0);
|
||||
doublereal mv = m_sub->molarVolume(T, dens0);
|
||||
m_sub->setState(m_temp, m_dens);
|
||||
m_sub->setState_TR(m_temp, m_dens);
|
||||
return (mv);
|
||||
}
|
||||
|
||||
|
|
@ -465,7 +465,7 @@ namespace Cantera {
|
|||
|
||||
void PDSS_Water::setDensity(doublereal dens) {
|
||||
m_dens = dens;
|
||||
m_sub->setState(m_temp, m_dens);
|
||||
m_sub->setState_TR(m_temp, m_dens);
|
||||
}
|
||||
|
||||
doublereal PDSS_Water::density() const {
|
||||
|
|
@ -475,7 +475,7 @@ namespace Cantera {
|
|||
void PDSS_Water::setTemperature(doublereal temp) {
|
||||
m_temp = temp;
|
||||
doublereal dd = m_dens;
|
||||
m_sub->setState(temp, dd);
|
||||
m_sub->setState_TR(temp, dd);
|
||||
}
|
||||
|
||||
void PDSS_Water::setState_TP(doublereal temp, doublereal pres) {
|
||||
|
|
|
|||
|
|
@ -440,7 +440,7 @@ namespace Cantera {
|
|||
skip = true;
|
||||
else
|
||||
throw CanteraError("importPhase",
|
||||
"duplicate species: "+name);
|
||||
"duplicate species: \"" + name + "\"");
|
||||
}
|
||||
if (!skip) {
|
||||
declared[name] = true;
|
||||
|
|
@ -453,8 +453,8 @@ namespace Cantera {
|
|||
++k;
|
||||
}
|
||||
else {
|
||||
throw CanteraError("importPhase","no data for species "
|
||||
+name);
|
||||
throw CanteraError("importPhase","no data for species, \""
|
||||
+ name + "\"");
|
||||
}
|
||||
}
|
||||
}
|
||||
|
|
|
|||
|
|
@ -14,6 +14,8 @@
|
|||
*/
|
||||
|
||||
#include "WaterPropsIAPWS.h"
|
||||
#include "ctexceptions.h"
|
||||
#include "stringUtils.h"
|
||||
#include <math.h>
|
||||
#include <stdio.h>
|
||||
#include <stdlib.h>
|
||||
|
|
@ -95,7 +97,7 @@ void WaterPropsIAPWS::calcDim(double temperature, double rho) {
|
|||
* J kmol-1 K-1.
|
||||
*/
|
||||
double WaterPropsIAPWS::helmholtzFE(double temperature, double rho) {
|
||||
setState(temperature, rho);
|
||||
setState_TR(temperature, rho);
|
||||
double retn = m_phi->phi(tau, delta);
|
||||
double RT = Rgas * temperature;
|
||||
return (retn * RT);
|
||||
|
|
@ -155,14 +157,20 @@ density(double temperature, double pressure, int phase, double rhoguess) {
|
|||
if (temperature > T_c) {
|
||||
rhoguess = pressure * M_water / (Rgas * temperature);
|
||||
} else {
|
||||
if (phase != WATER_LIQUID) {
|
||||
if (phase == WATER_GAS || phase == WATER_SUPERCRIT) {
|
||||
rhoguess = pressure * M_water / (Rgas * temperature);
|
||||
} else {
|
||||
} else if (phase == WATER_LIQUID) {
|
||||
/*
|
||||
* Provide a guess about the liquid density that is
|
||||
* relatively high -> convergnce from above seems robust.
|
||||
*/
|
||||
rhoguess = 1000.;
|
||||
} else if (phase == WATER_UNSTABLELIQUID || phase == WATER_UNSTABLEGAS) {
|
||||
throw Cantera::CanteraError("WaterPropsIAPWS::density",
|
||||
"Unstable Branch finder is untested");
|
||||
} else {
|
||||
throw Cantera::CanteraError("WaterPropsIAPWS::density",
|
||||
"unknown state: " + Cantera::int2str(phase));
|
||||
}
|
||||
}
|
||||
} else {
|
||||
|
|
@ -176,7 +184,7 @@ density(double temperature, double pressure, int phase, double rhoguess) {
|
|||
}
|
||||
double p_red = pressure * M_water / (Rgas * temperature * Rho_c);
|
||||
deltaGuess = rhoguess / Rho_c;
|
||||
setState(temperature, rhoguess);
|
||||
setState_TR(temperature, rhoguess);
|
||||
double delta_retn = m_phi->dfind(p_red, tau, deltaGuess);
|
||||
double density_retn;
|
||||
if (delta_retn >0.0) {
|
||||
|
|
@ -190,7 +198,7 @@ density(double temperature, double pressure, int phase, double rhoguess) {
|
|||
* Set the internal state -> this may be
|
||||
* a duplication. However, let's just be sure.
|
||||
*/
|
||||
setState(temperature, density_retn);
|
||||
setState_TR(temperature, density_retn);
|
||||
|
||||
|
||||
} else {
|
||||
|
|
@ -302,12 +310,17 @@ double WaterPropsIAPWS::isothermalCompressibility() const {
|
|||
return (1.0 / (dens * dpdrho));
|
||||
}
|
||||
|
||||
double WaterPropsIAPWS:: coeffPresExp() const {
|
||||
double retn = m_phi->dimdpdT(tau, delta);
|
||||
return (retn);
|
||||
}
|
||||
|
||||
/*
|
||||
* Calculate the Gibbs free energy in mks units of
|
||||
* J kmol-1 K-1.
|
||||
*/
|
||||
double WaterPropsIAPWS::Gibbs(double temperature, double rho) {
|
||||
setState(temperature, rho);
|
||||
setState_TR(temperature, rho);
|
||||
double gRT = m_phi->gibbs_RT();
|
||||
return (gRT * Rgas * temperature);
|
||||
}
|
||||
|
|
@ -331,7 +344,7 @@ corr(double temperature, double pressure, double &densLiq,
|
|||
printf("error liq\n");
|
||||
exit(-1);
|
||||
}
|
||||
setState(temperature, densLiq);
|
||||
setState_TR(temperature, densLiq);
|
||||
double gibbsLiqRT = m_phi->gibbs_RT();
|
||||
|
||||
densGas = density(temperature, pressure, WATER_GAS, densGas);
|
||||
|
|
@ -339,7 +352,7 @@ corr(double temperature, double pressure, double &densLiq,
|
|||
printf("error gas\n");
|
||||
exit(-1);
|
||||
}
|
||||
setState(temperature, densGas);
|
||||
setState_TR(temperature, densGas);
|
||||
double gibbsGasRT = m_phi->gibbs_RT();
|
||||
|
||||
delGRT = gibbsLiqRT - gibbsGasRT;
|
||||
|
|
@ -350,11 +363,11 @@ corr1(double temperature, double pressure, double &densLiq,
|
|||
double &densGas, double &pcorr) {
|
||||
|
||||
densLiq = density(temperature, pressure, WATER_LIQUID, densLiq);
|
||||
setState(temperature, densLiq);
|
||||
setState_TR(temperature, densLiq);
|
||||
double prL = m_phi->phiR();
|
||||
|
||||
densGas = density(temperature, pressure, WATER_GAS, densGas);
|
||||
setState(temperature, densGas);
|
||||
setState_TR(temperature, densGas);
|
||||
double prG = m_phi->phiR();
|
||||
|
||||
double rhs = (prL - prG) + log(densLiq/densGas);
|
||||
|
|
@ -405,7 +418,7 @@ int WaterPropsIAPWS::phaseState() const {
|
|||
* Sets the internal state of the object to the
|
||||
* specified temperature and density.
|
||||
*/
|
||||
void WaterPropsIAPWS::setState(double temperature, double rho) {
|
||||
void WaterPropsIAPWS::setState_TR(double temperature, double rho) {
|
||||
calcDim(temperature, rho);
|
||||
m_phi->tdpolycalc(tau, delta);
|
||||
}
|
||||
|
|
@ -417,7 +430,7 @@ void WaterPropsIAPWS::setState(double temperature, double rho) {
|
|||
*/
|
||||
double WaterPropsIAPWS::
|
||||
enthalpy(double temperature, double rho) {
|
||||
setState(temperature, rho);
|
||||
setState_TR(temperature, rho);
|
||||
double hRT = m_phi->enthalpy_RT();
|
||||
return (hRT * Rgas * temperature);
|
||||
}
|
||||
|
|
@ -435,7 +448,7 @@ enthalpy() const {
|
|||
*/
|
||||
double WaterPropsIAPWS::
|
||||
intEnergy(double temperature, double rho) {
|
||||
setState(temperature, rho);
|
||||
setState_TR(temperature, rho);
|
||||
double uRT = m_phi->intEnergy_RT();
|
||||
return (uRT * Rgas * temperature);
|
||||
}
|
||||
|
|
@ -452,7 +465,7 @@ intEnergy() const{
|
|||
*/
|
||||
double WaterPropsIAPWS::
|
||||
entropy(double temperature, double rho) {
|
||||
setState(temperature, rho);
|
||||
setState_TR(temperature, rho);
|
||||
double sR = m_phi->entropy_R();
|
||||
return (sR * Rgas);
|
||||
}
|
||||
|
|
@ -471,7 +484,7 @@ double WaterPropsIAPWS::entropy() const {
|
|||
* J kmol-1 K-1.
|
||||
*/
|
||||
double WaterPropsIAPWS::cv(double temperature, double rho) {
|
||||
setState(temperature, rho);
|
||||
setState_TR(temperature, rho);
|
||||
double cvR = m_phi->cv_R();
|
||||
return (cvR * Rgas);
|
||||
}
|
||||
|
|
@ -480,27 +493,23 @@ double WaterPropsIAPWS::cv(double temperature, double rho) {
|
|||
* Calculate heat capacity at constant pressure
|
||||
* J kmol-1 K-1.
|
||||
*/
|
||||
double WaterPropsIAPWS::
|
||||
cp(double temperature, double rho) {
|
||||
setState(temperature, rho);
|
||||
double WaterPropsIAPWS::cp(double temperature, double rho) {
|
||||
setState_TR(temperature, rho);
|
||||
double cpR = m_phi->cp_R();
|
||||
return (cpR * Rgas);
|
||||
}
|
||||
|
||||
double WaterPropsIAPWS::
|
||||
cp() const {
|
||||
double WaterPropsIAPWS::cp() const {
|
||||
double cpR = m_phi->cp_R();
|
||||
return (cpR * Rgas);
|
||||
}
|
||||
|
||||
double WaterPropsIAPWS::
|
||||
molarVolume(double temperature, double rho) {
|
||||
setState(temperature, rho);
|
||||
double WaterPropsIAPWS::molarVolume(double temperature, double rho) {
|
||||
setState_TR(temperature, rho);
|
||||
return (M_water / rho);
|
||||
}
|
||||
|
||||
double WaterPropsIAPWS::
|
||||
molarVolume() const {
|
||||
double WaterPropsIAPWS::molarVolume() const {
|
||||
double rho = delta * Rho_c;
|
||||
return (M_water / rho);
|
||||
}
|
||||
|
|
|
|||
|
|
@ -22,12 +22,20 @@
|
|||
* @name Names for the phase regions
|
||||
*
|
||||
* These constants are defined and used in the interface
|
||||
* to describe desired phases.
|
||||
* to describe the location of where we are in (T,rho) space.
|
||||
*
|
||||
* WATER_UNSTABLELIQUID indicates that we are in the unstable region, inside the
|
||||
* spinodal curve where dpdrho < 0.0 amonst other properties. The difference
|
||||
* between WATER_UNSTABLELIQUID and WATER_UNSTABLEGAS is that
|
||||
* for WATER_UNSTABLELIQUID d2pdrho2 > 0 and dpdrho < 0.0
|
||||
* for WATER_UNSTABLEGAS d2pdrho2 < 0 and dpdrho < 0.0
|
||||
*/
|
||||
//@{
|
||||
#define WATER_GAS 0
|
||||
#define WATER_LIQUID 1
|
||||
#define WATER_SUPERCRIT 2
|
||||
#define WATER_UNSTABLELIQUID 3
|
||||
#define WATER_UNSTABLEGAS 4
|
||||
//@}
|
||||
|
||||
//! Class for calculating the equation of state of water.
|
||||
|
|
@ -155,7 +163,7 @@ public:
|
|||
* @param temperature temperature (kelvin)
|
||||
* @param rho density (kg m-3)
|
||||
*/
|
||||
void setState(double temperature, double rho);
|
||||
void setState_TR(double temperature, double rho);
|
||||
|
||||
|
||||
//! Calculate the Helmholtz free energy in mks units of J kmol-1 K-1.
|
||||
|
|
@ -180,7 +188,6 @@ public:
|
|||
//! using the last temperature and density
|
||||
double Gibbs() const;
|
||||
|
||||
|
||||
//! Calculate the enthalpy in mks units of J kmol-1
|
||||
/*!
|
||||
* @param temperature temperature (kelvin)
|
||||
|
|
@ -311,6 +318,18 @@ public:
|
|||
*/
|
||||
double coeffThermExp(double temperature, double pressure);
|
||||
|
||||
//! Returns the isochoric pressure-temperature coefficient
|
||||
/*!
|
||||
*
|
||||
* beta = M / (rho * Rgas) (d (pressure) / dT) at constant rho
|
||||
*
|
||||
* Note for ideal gases this is equal to one.
|
||||
*
|
||||
* beta = delta (phi0_d() + phiR_d())
|
||||
* - tau delta (phi0_dt() + phiR_dt())
|
||||
*/
|
||||
double coeffPresExp() const;
|
||||
|
||||
//! Returns the coefficient of isothermal compressibility for the
|
||||
//! state of the object
|
||||
/*!
|
||||
|
|
|
|||
|
|
@ -784,7 +784,13 @@ double WaterPropsIAPWSphi::dimdpdrho(double tau, double delta) {
|
|||
return retn;
|
||||
}
|
||||
|
||||
|
||||
double WaterPropsIAPWSphi::dimdpdT(double tau, double delta) {
|
||||
tdpolycalc(tau, delta);
|
||||
double res1 = phiR_d();
|
||||
double res2 = phiR_dt();
|
||||
double retn = (1.0 + delta * res1) - tau * delta * (res2);
|
||||
return retn;
|
||||
}
|
||||
|
||||
/*
|
||||
* Calculate d_phi0/d(tau)
|
||||
|
|
|
|||
|
|
@ -73,13 +73,6 @@ public:
|
|||
*/
|
||||
double phi_tt(double tau, double delta);
|
||||
|
||||
//! Second derivative of phi wrt tau, then delta
|
||||
/*!
|
||||
* @param tau Dimensionless temperature = T_c/T
|
||||
* @param delta Dimensionless density = delta = rho / Rho_c
|
||||
*/
|
||||
double phi_dt(double tau, double delta);
|
||||
|
||||
//! Internal check # 1
|
||||
void check1();
|
||||
|
||||
|
|
@ -108,6 +101,16 @@ public:
|
|||
*/
|
||||
double dimdpdrho(double tau, double delta);
|
||||
|
||||
//! Dimensionless derivative of p wrt T at constant rho
|
||||
/*!
|
||||
* dp/dT * M/(Rho R) = (1.0 + delta phiR_d()
|
||||
* - tau delta (phiR_dt())
|
||||
*
|
||||
* @param tau Dimensionless temperature = T_c/T
|
||||
* @param delta Dimensionless density = delta = rho / Rho_c
|
||||
*/
|
||||
double dimdpdT(double tau, double delta);
|
||||
|
||||
/**
|
||||
* This program computes the reduced density, given the reduced pressure
|
||||
* and the reduced temperature, tau. It takes an initial guess, deltaGuess.
|
||||
|
|
|
|||
|
|
@ -398,7 +398,7 @@ namespace Cantera {
|
|||
if (dd <= 0.0) {
|
||||
throw CanteraError("setPressure", "error");
|
||||
}
|
||||
m_sub->setState(T, dd);
|
||||
m_sub->setState_TR(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);
|
||||
|
|
@ -427,7 +427,7 @@ namespace Cantera {
|
|||
if (dd <= 0.0) {
|
||||
throw CanteraError("setPressure", "error");
|
||||
}
|
||||
m_sub->setState(T, dd);
|
||||
m_sub->setState_TR(T, dd);
|
||||
|
||||
doublereal s = m_sub->entropy(T, dd);
|
||||
*sr = (s + SW_Offset)/ (GasConstant);
|
||||
|
|
@ -445,7 +445,7 @@ namespace Cantera {
|
|||
waterState = WATER_LIQUID;
|
||||
}
|
||||
doublereal dd = m_sub->density(T, OneAtm, waterState, dens);
|
||||
m_sub->setState(T, dd);
|
||||
m_sub->setState_TR(T, dd);
|
||||
if (dd <= 0.0) {
|
||||
throw CanteraError("setPressure", "error");
|
||||
}
|
||||
|
|
@ -539,7 +539,7 @@ namespace Cantera {
|
|||
void WaterSSTP::setTemperature(double temp) {
|
||||
State::setTemperature(temp);
|
||||
doublereal dd = density();
|
||||
m_sub->setState(temp, dd);
|
||||
m_sub->setState_TR(temp, dd);
|
||||
}
|
||||
|
||||
|
||||
|
|
@ -549,7 +549,7 @@ namespace Cantera {
|
|||
doublereal tsave = temperature();
|
||||
doublereal dsave = density();
|
||||
doublereal pp = m_sub->psat(t);
|
||||
m_sub->setState(tsave, dsave);
|
||||
m_sub->setState_TR(tsave, dsave);
|
||||
return pp;
|
||||
}
|
||||
|
||||
|
|
|
|||
Loading…
Add table
Reference in a new issue