30 HighPressureGasTransport::HighPressureGasTransport(
thermo_t* thermo)
39 doublereal Lprime_m = 0.0;
40 const doublereal c1 = 1./16.04;
52 m_thermo -> getPartialMolarVolumes(&V_k[0]);
55 for (
size_t i = 0; i <
m_nsp; i++) {
56 doublereal Tc_i = Tcrit_i(i);
57 doublereal Vc_i = Vcrit_i(i);
59 doublereal V_r = V_k[i]/Vc_i;
60 doublereal T_p = std::min(T_r,2.0);
61 doublereal V_p = std::max(0.5,std::min(V_r,2.0));
64 doublereal theta_p = 1.0 + (
m_w_ac[i] - 0.011)*(0.56553
65 - 0.86276*log(T_p) - 0.69852/T_p);
66 doublereal phi_p = (1.0 + (
m_w_ac[i] - 0.011)*(0.38560
67 - 1.1617*log(T_p)))*0.288/Zcrit_i(i);
68 doublereal f_fac = Tc_i*theta_p/190.4;
69 doublereal h_fac = 1000*Vc_i*phi_p/99.2;
70 doublereal T_0 =
m_temp/f_fac;
71 doublereal mu_0 = 1e-7*(2.90774e6/T_0 - 3.31287e6*pow(T_0,-2./3.)
72 + 1.60810e6*pow(T_0,-1./3.) - 4.33190e5 + 7.06248e4*pow(T_0,1./3.)
73 - 7.11662e3*pow(T_0,2./3.) + 4.32517e2*T_0 - 1.44591e1*pow(T_0,4./3.)
74 + 2.03712e-1*pow(T_0,5./3.));
75 doublereal H = sqrt(f_fac*16.04/
m_mw[i])*pow(h_fac,-2./3.);
76 doublereal mu_i = mu_0*H*
m_mw[i]*c1;
78 L_i_min = min(L_i_min,L_i[i]);
80 doublereal theta_s = 1 + (
m_w_ac[i] - 0.011)*(0.09057 - 0.86276*log(T_p)
81 + (0.31664 - 0.46568/T_p)*(V_p - 0.5));
82 doublereal phi_s = (1 + (
m_w_ac[i] - 0.011)*(0.39490*(V_p - 1.02355)
83 - 0.93281*(V_p - 0.75464)*log(T_p)))*0.288/Zcrit_i(i);
84 f_i[i] = Tc_i*theta_s/190.4;
85 h_i[i] = 1000*Vc_i*phi_s/99.2;
91 for (
size_t i = 0; i <
m_nsp; i++) {
92 for (
size_t j = 0; j <
m_nsp; j++) {
94 doublereal L_ij = 2*L_i[i]*L_i[j]/(L_i[i] + L_i[j] +
Tiny);
95 Lprime_m += molefracs[i]*molefracs[j]*L_ij;
97 doublereal f_ij = sqrt(f_i[i]*f_i[j]);
98 doublereal h_ij = 0.125*pow(pow(h_i[i],1./3.) + pow(h_i[j],1./3.),3.);
100 f_m += molefracs[i]*molefracs[j]*f_ij*h_ij;
101 h_m += molefracs[i]*molefracs[j]*h_ij;
102 mw_m += molefracs[i]*molefracs[j]*sqrt(mw_ij_inv*f_ij)*pow(h_ij,-4./3.);
107 mw_m = pow(mw_m,-2.)*f_m*pow(h_m,-8./3.);
110 doublereal T_0 =
m_temp/f_m;
111 doublereal mu_0 = 1e-7*(2.90774e6/T_0 - 3.31287e6*pow(T_0,-2./3.)
112 + 1.60810e6*pow(T_0,-1./3.) - 4.33190e5 + 7.06248e4
113 *pow(T_0,1./3.) - 7.11662e3*pow(T_0,2./3.) + 4.32517e2*T_0
114 - 1.44591e1*pow(T_0,4./3.) + 2.03712e-1*pow(T_0,5./3.));
115 doublereal L_1m = 1944*mu_0;
116 doublereal L_2m = (-2.5276e-4 + 3.3433e-4*pow(1.12 - log(T_0/1.680e2),2))*rho_0;
117 doublereal L_3m = exp(-7.19771 + 85.67822/T_0)*(exp((12.47183
118 - 984.6252*pow(T_0,-1.5))*pow(rho_0,0.1) + (rho_0/0.1617 - 1)
119 *sqrt(rho_0)*(0.3594685 + 69.79841/T_0 - 872.8833*pow(T_0,-2))) - 1.)*1e-3;
120 doublereal H_m = sqrt(f_m*16.04/mw_m)*pow(h_m,-2./3.);
121 doublereal Lstar_m = H_m*(L_1m + L_2m + L_3m);
122 return Lprime_m + Lstar_m;
133 throw CanteraError(
"HighPressureGasTransport::getThermalDiffCoeffs",
134 "Not yet implemented.");
152 throw CanteraError(
"HighPressureTransport::getBinaryDiffCoeffs()",
"ld is too small");
155 for (
size_t i = 0; i < nsp; i++) {
156 for (
size_t j = 0; j < nsp; j++) {
159 doublereal x_i = std::max(
Tiny, molefracs[i]);
160 doublereal x_j = std::max(
Tiny, molefracs[j]);
163 x_i = x_i/(x_i + x_j);
164 x_j = x_j/(x_i + x_j);
167 double Tr_ij =
m_temp/(x_i*Tcrit_i(i) + x_j*Tcrit_i(j));
177 P_corr_ij = setPcorr(Pr_ij, Tr_ij);
188 d[ld*j + i] = P_corr_ij*rp *
m_bdiff(i,j);
197 throw CanteraError(
"HighPressureTransport:getMultiDiffCoeffs()",
198 "Routine not yet implemented");
220 throw CanteraError(
"HighPressureTransport::getMultiDiffCoeffs()",
223 for (
size_t i = 0; i <
m_nsp; i++) {
224 for (
size_t j = 0; j <
m_nsp; j++) {
227 doublereal x_i = std::max(
Tiny, molefracs[i]);
228 doublereal x_j = std::max(
Tiny, molefracs[j]);
231 double Tr_ij =
m_temp/(x_i*Tcrit_i(i) + x_j*Tcrit_i(j));
238 P_corr_ij = setPcorr(Pr_ij, Tr_ij);
261 throw CanteraError(
"HighPressureGasTransport::getMultiDiffCoeffs",
262 "invert returned ierr = {}", ierr);
265 m_lmatrix_soln_ok =
false;
267 doublereal prefactor = 16.0 *
m_temp 270 for (
size_t i = 0; i <
m_nsp; i++) {
271 for (
size_t j = 0; j <
m_nsp; j++) {
272 double c = prefactor/
m_mw[j];
273 d[ld*j + i] = c*molefracs[i]*(m_Lmatrix(i,j) - m_Lmatrix(i,i));
282 double Pc_mix_n = 0.;
283 double Pc_mix_d = 0.;
285 double MW_H =
m_mw[0];
286 double MW_L =
m_mw[0];
287 doublereal FP_mix_o = 0;
288 doublereal FQ_mix_o = 0;
295 double x_H = molefracs[0];
296 for (
size_t i = 0; i <
m_nsp; i++) {
299 double Tc = Tcrit_i(i);
300 double Tr = tKelvin/Tc;
301 double Zc = Zcrit_i(i);
302 Tc_mix += Tc*molefracs[i];
303 Pc_mix_n += molefracs[i]*Zc;
304 Pc_mix_d += molefracs[i]*Vcrit_i(i);
307 if (
m_mw[i] > MW_H) {
310 }
else if (
m_mw[i] < MW_L) {
317 FP_mix_o += molefracs[i];
318 }
else if (mu_ri < 0.075) {
319 FP_mix_o += molefracs[i]*(1. + 30.55*pow(0.292 - Zc, 1.72));
320 }
else { FP_mix_o += molefracs[i]*(1. + 30.55*pow(0.292 - Zc, 1.72)
321 *fabs(0.96 + 0.1*(Tr - 0.7)));
330 if (spnames[i] ==
"He") {
331 FQ_mix_o += molefracs[i]*FQ_i(1.38,Tr,
m_mw[i]);
332 }
else if (spnames[i] ==
"H2") {
333 FQ_mix_o += molefracs[i]*(FQ_i(0.76,Tr,
m_mw[i]));
334 }
else if (spnames[i] ==
"D2") {
335 FQ_mix_o += molefracs[i]*(FQ_i(0.52,Tr,
m_mw[i]));
337 FQ_mix_o += molefracs[i];
341 double Tr_mix = tKelvin/Tc_mix;
342 double Pc_mix =
GasConstant*Tc_mix*Pc_mix_n/Pc_mix_d;
344 double ratio = MW_H/MW_L;
345 double ksi = pow(
GasConstant*Tc_mix*3.6277*pow(10.0,53.0)/(pow(MW_mix,3)
346 *pow(Pc_mix,4)),1.0/6.0);
348 if (ratio > 9 && x_H > 0.05 && x_H < 0.7) {
349 FQ_mix_o *= 1 - 0.01*pow(ratio,0.87);
353 double Z1m = (0.807*pow(Tr_mix,0.618) - 0.357*exp(-0.449*Tr_mix)
354 + 0.340*exp(-4.058*Tr_mix)+0.018)*FP_mix_o*FQ_mix_o;
359 if (Pr_mix < Pvp_mix/Pc_mix) {
360 doublereal alpha = 3.262 + 14.98*pow(Pr_mix,5.508);
361 doublereal beta = 1.390 + 5.746*Pr_mix;
362 Z2m = 0.600 + 0.760*pow(Pr_mix,alpha) + (0.6990*pow(Pr_mix,beta) -
365 throw CanteraError(
"HighPressureGasTransport::viscosity",
366 "State is outside the limits of the Lucas model, Tr <= 1");
368 }
else if ((Tr_mix > 1.0) && (Tr_mix < 40.0)) {
369 if ((Pr_mix > 0.0) && (Pr_mix <= 100.0)) {
370 doublereal a_fac = 0.001245*exp(5.1726*pow(Tr_mix,-0.3286))/Tr_mix;
371 doublereal b_fac = a_fac*(1.6553*Tr_mix - 1.2723);
372 doublereal c_fac = 0.4489*exp(3.0578*pow(Tr_mix,-37.7332))/Tr_mix;
373 doublereal d_fac = 1.7368*exp(2.2310*pow(Tr_mix,-7.6351))/Tr_mix;
374 doublereal f_fac = 0.9425*exp(-0.1853*pow(Tr_mix,0.4489));
376 Z2m = Z1m*(1 + a_fac*pow(Pr_mix,1.3088)/(b_fac*pow(Pr_mix,f_fac)
377 + pow(1+c_fac*pow(Pr_mix,d_fac),-1)));
379 throw CanteraError(
"HighPressureGasTransport::viscosity",
380 "State is outside the limits of the Lucas model, 1.0 < Tr < 40");
383 throw CanteraError(
"HighPressureGasTransport::viscosity",
384 "State is outside the limits of the Lucas model, Tr > 40");
388 doublereal Y = Z2m/Z1m;
391 return Z2m*(1 + (FP_mix_o - 1)*pow(Y,-3))*(1 + (FQ_mix_o - 1)
392 *(1/Y - 0.007*pow(log(Y),4)))/(ksi*FP_mix_o*FQ_mix_o);
396 doublereal HighPressureGasTransport::Tcrit_i(
size_t i)
407 doublereal HighPressureGasTransport::Pcrit_i(
size_t i)
418 doublereal HighPressureGasTransport::Vcrit_i(
size_t i)
429 doublereal HighPressureGasTransport::Zcrit_i(
size_t i)
440 vector_fp HighPressureGasTransport::store(
size_t i,
size_t nsp)
452 doublereal HighPressureGasTransport::FQ_i(doublereal Q, doublereal Tr, doublereal MW)
454 return 1.22*pow(Q,0.15)*(1 + 0.00385*pow(pow(Tr - 12.,2.),1./MW)
455 *fabs(Tr-12)/(Tr-12));
460 doublereal HighPressureGasTransport::setPcorr(doublereal Pr, doublereal Tr)
462 const static double Pr_lookup[17] = {0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.8, 1.0,
463 1.2, 1.4, 1.6, 1.8, 2.0, 2.5, 3.0, 4.0, 5.0};
464 const static double DP_Rt_lookup[17] = {1.01, 1.01, 1.01, 1.01, 1.01, 1.01,
465 1.01, 1.02, 1.02, 1.02, 1.02, 1.03, 1.03, 1.04, 1.05, 1.06, 1.07};
466 const static double A_ij_lookup[17] = {0.038042, 0.067433, 0.098317,
467 0.137610, 0.175081, 0.216376, 0.314051, 0.385736, 0.514553, 0.599184,
468 0.557725, 0.593007, 0.696001, 0.790770, 0.502100, 0.837452, 0.890390};
469 const static double B_ij_lookup[17] = {1.52267, 2.16794, 2.42910, 2.77605,
470 2.98256, 3.11384, 3.50264, 3.07773, 3.54744, 3.61216, 3.41882, 3.18415,
471 3.37660, 3.27984, 3.39031, 3.23513, 3.13001};
472 const static double C_ij_lookup[17] = {0., 0., 0., 0., 0., 0., 0., 0.141211,
473 0.278407, 0.372683, 0.504894, 0.678469, 0.665702, 0., 0.602907, 0., 0.};
474 const static double E_ij_lookup[17] = {1., 1., 1., 1., 1., 1., 1., 13.45454,
475 14., 10.00900, 8.57519, 10.37483, 11.21674, 1., 6.19043, 1., 1.};
482 frac = (Pr - Pr_lookup[0])/(Pr_lookup[1] - Pr_lookup[0]);
484 for (
int j = 1; j < 17; j++) {
485 if (Pr_lookup[j] > Pr) {
486 frac = (Pr - Pr_lookup[j-1])/(Pr_lookup[j] - Pr_lookup[j-1]);
498 doublereal P_corr_1 = DP_Rt_lookup[Pr_i]*(1.0 - A_ij_lookup[Pr_i]
499 *pow(Tr,-B_ij_lookup[Pr_i]))*(1-C_ij_lookup[Pr_i]
500 *pow(Tr,-E_ij_lookup[Pr_i]));
501 doublereal P_corr_2 = DP_Rt_lookup[Pr_i+1]*(1.0 - A_ij_lookup[Pr_i+1]
502 *pow(Tr,-B_ij_lookup[Pr_i+1]))*(1-C_ij_lookup[Pr_i+1]
503 *pow(Tr,-E_ij_lookup[Pr_i+1]));
504 return P_corr_1*(1.0-frac) + P_corr_2*frac;
doublereal molarVolume() const
Molar volume (m^3/kmol).
int invert(DenseMatrix &A, size_t nn)
invert A. A is overwritten with A^-1.
virtual void getThermalDiffCoeffs(doublereal *const dt)
Return the thermal diffusion coefficients (kg/m/s)
DenseMatrix m_dipole
The effective dipole moment for (i,j) collisions.
doublereal temperature() const
Temperature (K).
Various templated functions that carry out common vector operations (see Templated Utility Functions)...
Interface for class MultiTransport.
bool m_bindiff_ok
Update boolean for the binary diffusivities at unit pressure.
thermo_t * m_thermo
pointer to the object representing the phase
size_t nSpecies() const
Returns the number of species in the phase.
void update_T()
Update basic temperature-dependent quantities if the temperature has changed.
DenseMatrix m_bdiff
Matrix of binary diffusion coefficients at the reference pressure and the current temperature Size is...
virtual doublereal satPressure(doublereal t)
Return the saturation pressure given the temperature.
doublereal m_temp
Current value of the temperature at which the properties in this object are calculated (Kelvin)...
const std::vector< std::string > & speciesNames() const
Return a const reference to the vector of species names.
virtual void getCp_R_ref(doublereal *cprt) const
Returns the vector of nondimensional constant pressure heat capacities of the reference state at the ...
Header file defining class TransportFactory (see TransportFactory)
Base class for a phase with thermodynamic properties.
virtual doublereal critPressure() const
Critical pressure (Pa).
ThermoPhase object for the ideal gas equation of state - workhorse for Cantera (see Thermodynamic Pro...
virtual void updateDiff_T()
Update the binary diffusion coefficients.
const doublereal BigNumber
largest number to compare to inf.
void eval_L0000(const doublereal *const x)
Evaluate the L0000 matrices.
virtual doublereal critCompressibility() const
Critical compressibility (unitless).
virtual void getBinaryDiffCoeffs(const size_t ld, doublereal *const d)
virtual double thermalConductivity()
Returns the mixture thermal conductivity in W/m/K.
Base class for exceptions thrown by Cantera classes.
virtual doublereal critVolume() const
Critical volume (m3/kmol).
virtual void setMoleFractions(const doublereal *const x)
Set the mole fractions to the specified values.
void updateThermal_T()
Update the temperature-dependent terms needed to compute the thermal conductivity and thermal diffusi...
Interface for class HighPressureGasTransport.
void getMoleFractions(doublereal *const x) const
Get the species mole fraction vector.
virtual doublereal pressure() const
Return the thermodynamic pressure (Pa).
virtual void getMultiDiffCoeffs(const size_t ld, doublereal *const d)
Return the Multicomponent diffusion coefficients. Units: [m^2/s].
Class MultiTransport implements multicomponent transport properties for ideal gas mixtures...
virtual doublereal critTemperature() const
Critical temperature (K).
virtual doublereal viscosity()
Viscosity of the mixture (kg /m /s)
std::vector< double > vector_fp
Turn on the use of stl vectors for the basic array type within cantera Vector of doubles.
doublereal meanMolecularWeight() const
The mean molecular weight. Units: (kg/kmol)
vector_fp m_w_ac
Pitzer acentric factor.
const doublereal Tiny
Small number to compare differences of mole fractions against.
const doublereal GasConstant
Universal Gas Constant. [J/kmol/K].
Contains declarations for string manipulation functions within Cantera.
size_t m_nsp
Number of species.
Namespace for the Cantera kernel.
Class that holds the data that is read in from the XML file, and which is used for processing of the ...
void update_C()
Update basic concentration-dependent quantities if the concentrations have changed.
vector_fp m_mw
Local copy of the species molecular weights.