Cantera  3.1.0a1
WaterTransport.cpp
Go to the documentation of this file.
1 //! @file WaterTransport.cpp
2 
3 // This file is part of Cantera. See License.txt in the top-level directory or
4 // at https://cantera.org/license.txt for license and copyright information.
5 
10 
11 namespace {
12 
13 const double Tstar = 647.27;
14 const double rhoStar = 317.763; // kg / m3
15 const double presStar = 22.115E6; // Pa
16 const double muStar = 55.071E-6; //Pa s
17 
18 const double H[4] = {1., 0.978197, 0.579829, -0.202354};
19 const double Hij[6][7] = {
20  { 0.5132047, 0.2151778, -0.2818107, 0.1778064, -0.04176610, 0., 0.},
21  { 0.3205656, 0.7317883, -1.070786 , 0.4605040, 0., -0.01578386, 0.},
22  { 0., 1.241044 , -1.263184 , 0.2340379, 0., 0., 0.},
23  { 0., 1.476783 , 0., -0.4924179, 0.1600435, 0., -0.003629481},
24  {-0.7782567, 0.0 , 0., 0. , 0., 0., 0.},
25  { 0.1885447, 0.0 , 0., 0. , 0., 0., 0.},
26 };
27 
28 }
29 
30 namespace Cantera
31 {
32 
33 void WaterTransport::init(ThermoPhase* thermo, int mode, int log_level)
34 {
35  m_thermo = thermo;
36 }
37 
39 {
40  static const double TStar = 647.27; // Kelvin
41  double temp = m_thermo->temperature();
42  double dens = m_thermo->density();
43 
44  double rhobar = dens/rhoStar;
45  double tbar = temp / TStar;
46  double tbar2 = tbar * tbar;
47  double tbar3 = tbar2 * tbar;
48  double mu0bar = std::sqrt(tbar) / (H[0] + H[1]/tbar + H[2]/tbar2 + H[3]/tbar3);
49 
50  double tfac1 = 1.0 / tbar - 1.0;
51  double tfac2 = tfac1 * tfac1;
52  double tfac3 = tfac2 * tfac1;
53  double tfac4 = tfac3 * tfac1;
54  double tfac5 = tfac4 * tfac1;
55 
56  double rfac1 = rhobar - 1.0;
57  double rfac2 = rfac1 * rfac1;
58  double rfac3 = rfac2 * rfac1;
59  double rfac4 = rfac3 * rfac1;
60  double rfac5 = rfac4 * rfac1;
61  double rfac6 = rfac5 * rfac1;
62 
63  double sum = Hij[0][0] + Hij[1][0]*tfac1 + Hij[4][0]*tfac4 + Hij[5][0]*tfac5 +
64  Hij[0][1]*rfac1 + Hij[1][1]*tfac1*rfac1 + Hij[2][1]*tfac2*rfac1 + Hij[3][1]*tfac3*rfac1 +
65  Hij[0][2]*rfac2 + Hij[1][2]*tfac1*rfac2 + Hij[2][2]*tfac2*rfac2 +
66  Hij[0][3]*rfac3 + Hij[1][3]*tfac1*rfac3 + Hij[2][3]*tfac2*rfac3 + Hij[3][3]*tfac3*rfac3 +
67  Hij[0][4]*rfac4 + Hij[3][4]*tfac3*rfac4 +
68  Hij[1][5]*tfac1*rfac5 + Hij[3][6]*tfac3*rfac6;
69  double mu1bar = std::exp(rhobar * sum);
70 
71  // Apply the near-critical point corrections if necessary
72  double mu2bar = 1.0;
73  if (tbar >= 0.9970 && tbar <= 1.0082 && rhobar >= 0.755 && rhobar <= 1.290) {
74  double drhodp = m_thermo->isothermalCompressibility() * dens;
75  drhodp *= presStar / rhoStar;
76  double xsi = rhobar * drhodp;
77  if (xsi >= 21.93) {
78  mu2bar = 0.922 * std::pow(xsi, 0.0263);
79  }
80  }
81 
82  double mubar = mu0bar * mu1bar * mu2bar;
83  return mubar * muStar;
84 }
85 
87 {
88  static const double lambdastar = 0.4945;
89  static const double L[4] = {
90  1.0000,
91  6.978267,
92  2.599096,
93  -0.998254
94  };
95  static const double Lji[6][5] = {
96  { 1.3293046, 1.7018363, 5.2246158, 8.7127675, -1.8525999},
97  {-0.40452437, -2.2156845, -10.124111, -9.5000611, 0.93404690},
98  { 0.24409490, 1.6511057, 4.9874687, 4.3786606, 0.0},
99  { 0.018660751, -0.76736002, -0.27297694, -0.91783782, 0.0},
100  {-0.12961068, 0.37283344, -0.43083393, 0.0, 0.0},
101  { 0.044809953, -0.11203160, 0.13333849, 0.0, 0.0},
102  };
103 
104  double temp = m_thermo->temperature();
105  double dens = m_thermo->density();
106 
107  double rhobar = dens / rhoStar;
108  double tbar = temp / Tstar;
109  double tbar2 = tbar * tbar;
110  double tbar3 = tbar2 * tbar;
111  double lambda0bar = sqrt(tbar) / (L[0] + L[1]/tbar + L[2]/tbar2 + L[3]/tbar3);
112 
113  double tfac1 = 1.0 / tbar - 1.0;
114  double tfac2 = tfac1 * tfac1;
115  double tfac3 = tfac2 * tfac1;
116  double tfac4 = tfac3 * tfac1;
117  double tfac5 = tfac4 * tfac1;
118 
119  double rfac1 = rhobar - 1.0;
120  double rfac2 = rfac1 * rfac1;
121  double rfac3 = rfac2 * rfac1;
122  double rfac4 = rfac3 * rfac1;
123  double rfac5 = rfac4 * rfac1;
124  double rfac6 = rfac5 * rfac1;
125 
126  double sum = (Lji[0][0] + Lji[0][1]*tfac1 + Lji[0][2]*tfac2 + Lji[0][3]*tfac3 + Lji[0][4]*tfac4 +
127  Lji[1][0]*rfac1 + Lji[1][1]*tfac1*rfac1 + Lji[1][2]*tfac2*rfac1 + Lji[1][3]*tfac3*rfac1 + Lji[1][4]*tfac4*rfac1 +
128  Lji[2][0]*rfac2 + Lji[2][1]*tfac1*rfac2 + Lji[2][2]*tfac2*rfac2 + Lji[2][3]*tfac3*rfac2 +
129  Lji[3][0]*rfac3 + Lji[3][1]*tfac1*rfac3 + Lji[3][2]*tfac2*rfac3 + Lji[3][3]*tfac3*rfac3 +
130  Lji[4][0]*rfac4 + Lji[4][1]*tfac1*rfac4 + Lji[4][2]*tfac2*rfac4 +
131  Lji[5][0]*rfac5 + Lji[5][1]*tfac1*rfac5 + Lji[5][2]*tfac2*rfac5
132  );
133  double lambda1bar = exp(rhobar * sum);
134  double mu0bar = std::sqrt(tbar) / (H[0] + H[1]/tbar + H[2]/tbar2 + H[3]/tbar3);
135 
136  sum = (Hij[0][0] + Hij[1][0]*tfac1 + Hij[4][0]*tfac4 + Hij[5][0]*tfac5 +
137  Hij[0][1]*rfac1 + Hij[1][1]*tfac1*rfac1 + Hij[2][1]*tfac2*rfac1 + Hij[3][1]*tfac3*rfac1 +
138  Hij[0][2]*rfac2 + Hij[1][2]*tfac1*rfac2 + Hij[2][2]*tfac2*rfac2 +
139  Hij[0][3]*rfac3 + Hij[1][3]*tfac1*rfac3 + Hij[2][3]*tfac2*rfac3 + Hij[3][3]*tfac3*rfac3 +
140  Hij[0][4]*rfac4 + Hij[3][4]*tfac3*rfac4 +
141  Hij[1][5]*tfac1*rfac5 + Hij[3][6]*tfac3*rfac6
142  );
143  double mu1bar = std::exp(rhobar * sum);
144  double t2r2 = tbar2 / (rhobar * rhobar);
145  double kappa = m_thermo->isothermalCompressibility();
146  double xsi = rhobar * rhobar * kappa * presStar;
147  double xsipow = std::pow(xsi, 0.4678);
148  double temp2 = (tbar - 1.0) * (tbar - 1.0);
149  double dpdT_const_rho = m_thermo->thermalExpansionCoeff() / kappa;
150  dpdT_const_rho *= Tstar / presStar;
151  double lambda2bar = 0.0013848 / (mu0bar * mu1bar) * t2r2 * dpdT_const_rho * dpdT_const_rho *
152  xsipow * sqrt(rhobar) * exp(-18.66*temp2 - rfac4);
153  return (lambda0bar * lambda1bar + lambda2bar) * lambdastar;
154 }
155 
156 }
Implementation of a pressure dependent standard state virtual function for a Pure Water Phase (see Sp...
Header file for a derived class of ThermoPhase that handles variable pressure standard state methods ...
Declares a ThermoPhase class consisting of pure water (see Thermodynamic Properties and class WaterSS...
Header file defining class WaterTransport.
double temperature() const
Temperature (K).
Definition: Phase.h:562
virtual double density() const
Density (kg/m^3).
Definition: Phase.h:587
Base class for a phase with thermodynamic properties.
Definition: ThermoPhase.h:390
virtual double thermalExpansionCoeff() const
Return the volumetric thermal expansion coefficient. Units: 1/K.
Definition: ThermoPhase.h:580
virtual double isothermalCompressibility() const
Returns the isothermal compressibility. Units: 1/Pa.
Definition: ThermoPhase.h:569
ThermoPhase * m_thermo
pointer to the object representing the phase
Definition: Transport.h:420
ThermoPhase & thermo()
Phase object.
Definition: Transport.h:103
double thermalConductivity() override
Returns the thermal conductivity of water at the current conditions (W/m/K)
double viscosity() override
Returns the viscosity of water at the current conditions (kg/m/s)
void init(ThermoPhase *thermo, int mode=0, int log_level=0) override
Initialize a transport manager.
Namespace for the Cantera kernel.
Definition: AnyMap.cpp:564