research/cantera
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Added a viscosity model into WaterTransport.cpp for pure water

Browse source 
Added a test problem for the viscosity model.
Updated WaterSSTP object to support the use of WaterTransport 
with the object.
This commit is contained in:
Harry Moffat 2009-09-15 21:46:44 +00:00
parent 5d50aeef0e
commit 5009a3817f
13 changed files with 4643 additions and 5208 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

1
Cantera/src/base/vec_functions.h
View file

@ -192,6 +192,7 @@ namespace Cantera {
}
}
}
#endif

7
Cantera/src/kinetics/InterfaceKinetics.cpp
View file

@ -453,8 +453,9 @@ namespace Cantera {
*/
void InterfaceKinetics::getFwdRateConstants(doublereal* kfwd) {
_update_rates_T();
_update_rates_C();
// _update_rates_T();
// _update_rates_C();
updateROP();
const vector_fp& rf = m_kdata->m_rfn;
@ -492,7 +493,7 @@ namespace Cantera {
}
/**
* Update the rates of progress of the reactions in the reaciton
* Update the rates of progress of the reactions in the reaction
* mechanism. This routine operates on internal data.
*/
void InterfaceKinetics::updateROP() {

18
Cantera/src/thermo/WaterSSTP.cpp
View file

@ -16,6 +16,7 @@
#include "WaterSSTP.h"
#include "WaterPropsIAPWS.h"
//#include "importCTML.h"
#include "WaterProps.h"
#include "ThermoFactory.h"
#include <cmath>
@ -29,6 +30,7 @@ namespace Cantera {
WaterSSTP::WaterSSTP() :
SingleSpeciesTP(),
m_sub(0),
m_waterProps(0),
m_mw(0.0),
EW_Offset(0.0),
SW_Offset(0.0),
@ -43,6 +45,7 @@ namespace Cantera {
WaterSSTP::WaterSSTP(std::string inputFile, std::string id) :
SingleSpeciesTP(),
m_sub(0),
m_waterProps(0),
m_mw(0.0),
EW_Offset(0.0),
SW_Offset(0.0),
@ -57,6 +60,7 @@ namespace Cantera {
WaterSSTP::WaterSSTP(XML_Node& phaseRoot, std::string id) :
SingleSpeciesTP(),
m_sub(0),
m_waterProps(0),
m_mw(0.0),
EW_Offset(0.0),
SW_Offset(0.0),
@ -72,6 +76,7 @@ namespace Cantera {
WaterSSTP::WaterSSTP(const WaterSSTP &b) :
SingleSpeciesTP(b),
m_sub(0),
m_waterProps(0),
m_mw(b.m_mw),
EW_Offset(b.EW_Offset),
SW_Offset(b.SW_Offset),
@ -80,6 +85,8 @@ namespace Cantera {
m_allowGasPhase(b.m_allowGasPhase)
{
m_sub = new WaterPropsIAPWS(*(b.m_sub));
m_waterProps = new WaterProps(m_sub);
/*
* Use the assignment operator to do the brunt
* of the work for the copy construtor.
@ -93,6 +100,13 @@ namespace Cantera {
WaterSSTP& WaterSSTP::operator=(const WaterSSTP&b) {
if (&b == this) return *this;
m_sub->operator=(*(b.m_sub));
if (!m_waterProps) {
m_waterProps = new WaterProps(m_sub);
}
m_waterProps->operator=(*(b.m_waterProps));
m_mw = b.m_mw;
m_verbose = b.m_verbose;
m_ready = b.m_ready;
@ -108,6 +122,7 @@ namespace Cantera {
WaterSSTP::~WaterSSTP() {
delete m_sub;
delete m_waterProps;
}
@ -268,6 +283,9 @@ namespace Cantera {
double rho0 = m_sub->density(298.15, OneAtm, WATER_LIQUID);
setDensity(rho0);
m_waterProps = new WaterProps(m_sub);
/*
* We have to do something with the thermo function here.
*/

22
Cantera/src/thermo/WaterSSTP.h
View file

@ -21,6 +21,7 @@
namespace Cantera {
class WaterPropsIAPWS;
class WaterProps;
//! Class for single-component water. This is designed to cover just the
//! liquid part of water.
/*!
@ -501,6 +502,17 @@ namespace Cantera {
*/
virtual void setParametersFromXML(const XML_Node& eosdata);
//! Get a pointer to a changeable WaterPropsIAPWS object
WaterPropsIAPWS *getWater() {
return m_sub;
}
//! Get a pointer to a changeable WaterPropsIAPWS object
WaterProps *getWaterProps() {
return m_waterProps;
}
protected:
/**
@ -515,6 +527,16 @@ namespace Cantera {
//! of water.
mutable WaterPropsIAPWS *m_sub;
//! Pointer to the WaterProps object
/*!
* This class is used to house several approximation
* routines for properties of water.
*
* This object owns m_waterProps, and the WaterPropsIAPWS object used by
* WaterProps is m_sub, which is defined above.
*/
WaterProps *m_waterProps;
//! Molecular weight of Water -> Cantera assumption
doublereal m_mw;

134
Cantera/src/transport/WaterTransport.cpp
View file

@ -1,11 +1,14 @@
#include "ct_defs.h"
#include "WaterPropsIAPWS.h"
#include "TransportBase.h"
#include "DenseMatrix.h"
#include "LiquidTransportParams.h"
#include "VPStandardStateTP.h"
#include "WaterTransport.h"
#include "PDSS_Water.h"
#include "WaterSSTP.h"
#include "WaterProps.h"
#include <iostream>
using namespace std;
@ -17,7 +20,7 @@ namespace Cantera {
WaterTransport::WaterTransport(thermo_t* thermo, int ndim) :
Transport(thermo, ndim)
{
initTP();
}
// Copy Constructor for the %WaterThermo object.
@ -28,7 +31,6 @@ namespace Cantera {
Transport(right.m_thermo, right.m_nDim)
{
*this = right;
}
// Assignment operator
@ -43,7 +45,11 @@ namespace Cantera {
return *this;
}
Transport::operator=(right);
// All pointers in this routine are shallow pointers. Therefore, it's
// ok just to reinitialize them
initTP();
return *this;
}
@ -66,96 +72,44 @@ namespace Cantera {
WaterTransport::~WaterTransport() {
}
// Routine to do some common initializations at the start of using
// this routine.
void WaterTransport::initTP() {
// The expectation is that we have a VPStandardStateTP derived object
VPStandardStateTP *vpthermo = dynamic_cast<VPStandardStateTP *>(m_thermo);
if (!vpthermo) {
WaterSSTP *wsstp = dynamic_cast<WaterSSTP *>(m_thermo);
if (!wsstp) {
throw CanteraError("WaterTransport::initTP()",
"Expectation is that ThermoPhase be a VPStandardStateTP");
} else {
const double H[4] = {1.,
0.978197,
0.579829,
-0.202354};
const double Hij[6][7] =
{
{ 0.5132047, 0.2151778, -0.2818107, 0.1778064, -0.04176610, 0., 0.},
{ 0.3205656, 0.7317883, -1.070786 , 0.4605040, 0., -0.01578386, 0.},
{ 0., 1.241044 , -1.263184 , 0.2340379, 0., 0., 0.},
{ 0., 1.476783 , 0., -0.4924179, 0.1600435, 0., -0.003629481},
{-0.7782567, 0.0 , 0., 0. , 0., 0., 0.},
{ 0.1885447, 0.0 , 0., 0. , 0., 0., 0.},
};
const double TStar = 647.27; // Kelvin
const double rhoStar = 317.763; // kg / m3
const double presStar = 22.115E6; // Pa
const double muStar = 55.071E-6; //Pa s
m_sub = wsstp->getWater();
AssertTrace(m_sub != 0);
// Get a pointer to a changeable WaterProps object
m_waterProps = wsstp->getWaterProps();
AssertTrace(m_waterProps != 0);
}
} else {
m_waterPDSS = dynamic_cast<PDSS_Water *>(vpthermo->providePDSS(0));
if (!m_waterPDSS) {
throw CanteraError("WaterTransport::initTP()",
"Expectation is that first species be water with a PDSS_Water object");
}
// Get a pointer to a changeable WaterPropsIAPWS object
m_sub = m_waterPDSS->getWater();
AssertTrace(m_sub != 0);
// Get a pointer to a changeable WaterProps object
m_waterProps = m_waterPDSS->getWaterProps();
AssertTrace(m_waterProps != 0);
}
}
double WaterTransport::viscosity() {
double temp = m_thermo->temperature();
double dens = m_thermo->density();
WaterPropsIAPWS *waterP = new WaterPropsIAPWS();
waterP->setState_TR(temp, dens);
double pressure = waterP->pressure();
printf("pressure = %g\n", pressure);
//dens = 18.02 * pressure / (GasConstant * temp);
//printf ("mod dens = %g\n", dens);
double rhobar = dens/rhoStar;
double tbar = temp / TStar;
// double pbar = pressure / presStar;
double tbar2 = tbar * tbar;
double tbar3 = tbar2 * tbar;
double mu0bar = sqrt(tbar) / (H[0] + H[1]/tbar + H[2]/tbar2 + H[3]/tbar3);
printf("mu0bar = %g\n", mu0bar);
printf("mu0 = %g\n", mu0bar * muStar);
//double tfac0 = 1.0;
double tfac1 = 1.0 / tbar - 1.0;
double tfac2 = tfac1 * tfac1;
double tfac3 = tfac2 * tfac1;
double tfac4 = tfac3 * tfac1;
double tfac5 = tfac4 * tfac1;
//double rfac0 = 1.0;
double rfac1 = rhobar - 1.0;
double rfac2 = rfac1 * rfac1;
double rfac3 = rfac2 * rfac1;
double rfac4 = rfac3 * rfac1;
double rfac5 = rfac4 * rfac1;
double rfac6 = rfac5 * rfac1;
double sum = (Hij[0][0] + Hij[1][0]*tfac1 + Hij[4][0]*tfac4 + Hij[5][0]*tfac5 +
Hij[0][1]*rfac1 + Hij[1][1]*tfac1*rfac1 + Hij[2][1]*tfac2*rfac1 + Hij[3][1]*tfac3*rfac1 +
Hij[0][2]*rfac2 + Hij[1][2]*tfac1*rfac2 + Hij[2][2]*tfac2*rfac2 +
Hij[0][3]*rfac3 + Hij[1][3]*tfac1*rfac3 + Hij[2][3]*tfac2*rfac3 + Hij[3][3]*tfac3*rfac3 +
Hij[0][4]*rfac4 + Hij[3][4]*tfac3*rfac4 +
Hij[1][5]*tfac1*rfac5 + Hij[3][6]*tfac3*rfac6
);
double mu1bar = exp(rhobar * sum);
double mu2bar = 1.0;
double mubar = mu0bar * mu1bar * mu2bar;
return mubar * muStar;
double visc = m_waterProps->viscosityWater();
return visc;
}
}

64
Cantera/src/transport/WaterTransport.h
View file

@ -40,6 +40,8 @@ namespace Cantera {
class TransportParams;
class WaterProps;
class PDSS_Water;
class WaterTransport : public Transport {
@ -83,28 +85,54 @@ namespace Cantera {
return cWaterTransport;
}
//! overloaded base class methods
//! Returns the viscosity of the solution
//! Returns the viscosity of water at the current conditions
//! (kg/m/s)
/*!
* The viscosity is computed using the Wilke mixture rule.
* \f[
* \mu = \sum_k \frac{\mu_k X_k}{\sum_j \Phi_{k,j} X_j}.
* \f]
* Here \f$ \mu_k \f$ is the viscosity of pure species \e k,
* and
* \f[
* \Phi_{k,j} = \frac{\left[1
* + \sqrt{\left(\frac{\mu_k}{\mu_j}\sqrt{\frac{M_j}{M_k}}\right)}\right]^2}
* {\sqrt{8}\sqrt{1 + M_k/M_j}}
* \f]
* @see updateViscosity_T();
* This function calculates the value of the viscosity of pure
* water at the current T and P.
*
* Controlling update boolean m_viscmix_ok
*/
* The formulas used are from the paper
*
* J. V. Sengers, J. T. R. Watson, "Improved International
* Formulations for the Viscosity and Thermal Conductivity of
* Water Substance", J. Phys. Chem. Ref. Data, 15, 1291 (1986).
*
* The formulation is accurate for all temperatures and pressures,
* for steam and for water, even near the critical point.
* Pressures above 500 MPa and temperature above 900 C are suspect.
*/
virtual doublereal viscosity();
private:
//! Routine to do some common initializations at the start of using
//! this routine.
void initTP();
//! Pointer to the WaterPropsIAPWS object, which does the actual calculations
//! for the real equation of state
/*!
* This object owns m_sub
*/
mutable WaterPropsIAPWS *m_sub;
//! Pointer to the WaterProps object
/*!
* This class is used to house several approximation
* routines for properties of water.
*
* This object owns m_waterProps, and the WaterPropsIAPWS object used by
* WaterProps is m_sub, which is defined above.
*/
WaterProps *m_waterProps;
//! Pressure dependent standard state object for water
/*!
* We assume that species 0 is water, with a PDSS_Water object.
*/
PDSS_Water *m_waterPDSS;
};
}

9486
configure vendored
View file

File diff suppressed because one or more lines are too long

1
configure.in
View file

@ -1927,6 +1927,7 @@ AC_OUTPUT(Makefile \
test_problems/cathermo/DH_graph_bdotak/Makefile \
test_problems/cathermo/HMW_dupl_test/Makefile \
test_problems/cathermo/VPissp/Makefile \
test_problems/cathermo/wtWater/Makefile \
test_problems/VCSnonideal/Makefile \
test_problems/VPsilane_test/Makefile \
test_problems/VPsilane_test/runtest \

4
test_problems/cathermo/Makefile.in
View file

@ -33,6 +33,7 @@ ifeq ($(test_electrolytes),1)
cd DH_graph_Pitzer; @MAKE@ all
cd HMW_dupl_test; @MAKE@ all
cd VPissp; @MAKE@ all
cd wtWater; @MAKE@ all
endif
test:
@ -61,6 +62,7 @@ ifeq ($(test_electrolytes),1)
cd DH_graph_NM; @MAKE@ -s test
cd DH_graph_Pitzer; @MAKE@ -s test
cd HMW_dupl_test; @MAKE@ -s test
cd wtWater; @MAKE@ -s test
endif
#
# Have to remove .depends before calling make, because
@ -91,6 +93,7 @@ clean:
cd DH_graph_Pitzer; $(RM) .depends ; @MAKE@ clean
cd HMW_dupl_test; $(RM) .depends ; @MAKE@ clean
cd VPissp; $(RM) .depends ; @MAKE@ clean
cd wtWater; $(RM) .depends ; @MAKE@ clean
depends:
ifeq ($(test_issp),1)
@ -118,5 +121,6 @@ ifeq ($(test_electrolytes),1)
cd DH_graph_Pitzer; @MAKE@ depends
cd HMW_dupl_test; @MAKE@ depends
cd VPissp; @MAKE@ depends
cd wtWater; @MAKE@ depends
endif

4
test_problems/cathermo/wtWater/Makefile.in
View file

@ -11,11 +11,11 @@
.SUFFIXES : .d
# the name of the executable program to be created
PROG_NAME = testWaterPDSS
PROG_NAME = wtWater
# the object files to be linked together. List those generated from Fortran
# and from C/C++ separately
OBJS = testWaterPDSS.o
OBJS = wtWater.o
# additional flags to be passed to the linker. If your program
# requires other external libraries, put them here

12
test_problems/cathermo/wtWater/output_blessed.txt Normal file
View file

@ -0,0 +1,12 @@
-------------------------------------------------------------------------
T(C) MPa Phase Visc Visc(paper)
10-6 kg/m/s
-------------------------------------------------------------------------
25 0.1 L 890.496 890.5
100 0.1 L 281.807 281.9
100 10 L 284.457 284.5
250 5 L 106.405 106.4
250 50 L 117.43 117.5
350 17.5 L 66.9916 67.0
400 15 SC 24.9278 24.93
-------------------------------------------------------------------------

4
test_problems/cathermo/wtWater/runtest
View file

@ -4,14 +4,14 @@
temp_success="1"
/bin/rm -f output.txt outputa.txt
testName=testWaterPDSS
testName=wtWater
#################################################################
#
#################################################################
CANTERA_DATA=${CANTERA_DATA:=../../../data/inputs}; export CANTERA_DATA
CANTERA_BIN=${CANTERA_BIN:=../../../bin}
./testWaterPDSS > output.txt
./wtWater > output.txt
retnStat=$?
if [ $retnStat != "0" ]
then

94
test_problems/cathermo/wtWater/wtWater.cpp
View file

@ -28,30 +28,80 @@ double tvalue(double val, double atol = 1.0E-9) {
int main () {
double pres;
try {
try {
WaterSSTP * w = new WaterSSTP("waterTPphase.xml", "");
WaterSSTP * w = new WaterSSTP("waterTPphase.xml", "");
WaterTransport *wtTran = new WaterTransport(w, 3);
//double T = 273.15 + 150.0;
double T = 273.15 + 100.0;
w->setState_TP(T, 1.0E7);
double visc = wtTran->viscosity();
printf("visc = %g\n", visc);
delete w;
} catch (CanteraError) {
showErrors();
Cantera::appdelete();
return -1;
}
WaterTransport *wtTran = new WaterTransport(w, 3);
printf("-------------------------------------------------------------------------\n");
printf(" T(C) MPa Phase Visc Visc(paper) \n");
printf(" 10-6 kg/m/s \n");
printf("-------------------------------------------------------------------------\n");
double T = 273.15 + 25.0;
double pres = 1.0E5;
w->setState_TP(T, pres);
double visc = wtTran->viscosity();
printf("%8g %10.3g L %13.6g 890.5\n",
T - 273.15, pres * 1.0E-6, visc * 1.0E6);
return 0;
T = 273.15 + 100.0;
pres = 1.0E5;
w->setState_TP(T, pres);
visc = wtTran->viscosity();
printf("%8g %10.3g L %13.6g 281.9\n",
T - 273.15, pres * 1.0E-6, visc * 1.0E6);
T = 273.15 + 100.0;
pres = 1.0E7;
w->setState_TP(T, pres);
visc = wtTran->viscosity();
printf("%8g %10.3g L %13.6g 284.5\n",
T - 273.15, pres * 1.0E-6, visc * 1.0E6);
T = 273.15 + 250.0;
pres = 5.0E6;
w->setState_TP(T, pres);
visc = wtTran->viscosity();
printf("%8g %10.3g L %13.6g 106.4\n",
T - 273.15, pres * 1.0E-6, visc * 1.0E6);
T = 273.15 + 250.0;
pres = 5.0E7;
w->setState_TP(T, pres);
visc = wtTran->viscosity();
printf("%8g %10.3g L %13.6g 117.5\n",
T - 273.15, pres * 1.0E-6, visc * 1.0E6);
T = 273.15 + 350.0;
pres = 1.75E7;
w->setState_TP(T, pres);
visc = wtTran->viscosity();
printf("%8g %10.3g L %13.6g 67.0\n",
T - 273.15, pres * 1.0E-6, visc * 1.0E6);
T = 273.15 + 400.0;
pres = 1.50E7;
w->setState_TP(T, pres);
visc = wtTran->viscosity();
printf("%8g %10.3g SC %13.6g 24.93\n",
T - 273.15, pres * 1.0E-6, visc * 1.0E6);
printf("-------------------------------------------------------------------------\n");
delete w;
} catch (CanteraError) {
showErrors();
Cantera::appdelete();
return -1;
}
return 0;
}