Cantera  3.0.0
Loading...
Searching...
No Matches
HighPressureGasTransport.cpp
Go to the documentation of this file.
1/**
2 * @file HighPressureGasTransport.cpp
3 * Implementation file for class HighPressureGasTransport
4 **/
5
6// This file is part of Cantera. See License.txt in the top-level directory or
7// at https://cantera.org/license.txt for license and copyright information.
8
16
17using namespace std;
18
19namespace Cantera
20{
21
23: MultiTransport(thermo)
24{
25}
26
28{
29 // Method of Ely and Hanley:
30 update_T();
31 double Lprime_m = 0.0;
32 const double c1 = 1./16.04;
33 size_t nsp = m_thermo->nSpecies();
34 vector<double> molefracs(nsp);
35 m_thermo->getMoleFractions(&molefracs[0]);
36 vector<double> cp_0_R(nsp);
37 m_thermo->getCp_R_ref(&cp_0_R[0]);
38
39 vector<double> L_i(nsp);
40 vector<double> f_i(nsp);
41 vector<double> h_i(nsp);
42 vector<double> V_k(nsp);
43
44 m_thermo -> getPartialMolarVolumes(&V_k[0]);
45 double L_i_min = BigNumber;
46
47 for (size_t i = 0; i < m_nsp; i++) {
48 double Tc_i = Tcrit_i(i);
49 double Vc_i = Vcrit_i(i);
50 double T_r = m_thermo->temperature()/Tc_i;
51 double V_r = V_k[i]/Vc_i;
52 double T_p = std::min(T_r,2.0);
53 double V_p = std::max(0.5,std::min(V_r,2.0));
54
55 // Calculate variables for density-independent component:
56 double theta_p = 1.0 + (m_w_ac[i] - 0.011)*(0.56553
57 - 0.86276*log(T_p) - 0.69852/T_p);
58 double phi_p = (1.0 + (m_w_ac[i] - 0.011)*(0.38560
59 - 1.1617*log(T_p)))*0.288/Zcrit_i(i);
60 double f_fac = Tc_i*theta_p/190.4;
61 double h_fac = 1000*Vc_i*phi_p/99.2;
62 double T_0 = m_temp/f_fac;
63 double mu_0 = 1e-7*(2.90774e6/T_0 - 3.31287e6*pow(T_0,-2./3.)
64 + 1.60810e6*pow(T_0,-1./3.) - 4.33190e5 + 7.06248e4*pow(T_0,1./3.)
65 - 7.11662e3*pow(T_0,2./3.) + 4.32517e2*T_0 - 1.44591e1*pow(T_0,4./3.)
66 + 2.03712e-1*pow(T_0,5./3.));
67 double H = sqrt(f_fac*16.04/m_mw[i])*pow(h_fac,-2./3.);
68 double mu_i = mu_0*H*m_mw[i]*c1;
69 L_i[i] = mu_i*1.32*GasConstant*(cp_0_R[i] - 2.5)/m_mw[i];
70 L_i_min = min(L_i_min,L_i[i]);
71 // Calculate variables for density-dependent component:
72 double theta_s = 1 + (m_w_ac[i] - 0.011)*(0.09057 - 0.86276*log(T_p)
73 + (0.31664 - 0.46568/T_p)*(V_p - 0.5));
74 double phi_s = (1 + (m_w_ac[i] - 0.011)*(0.39490*(V_p - 1.02355)
75 - 0.93281*(V_p - 0.75464)*log(T_p)))*0.288/Zcrit_i(i);
76 f_i[i] = Tc_i*theta_s/190.4;
77 h_i[i] = 1000*Vc_i*phi_s/99.2;
78 }
79
80 double h_m = 0;
81 double f_m = 0;
82 double mw_m = 0;
83 for (size_t i = 0; i < m_nsp; i++) {
84 for (size_t j = 0; j < m_nsp; j++) {
85 // Density-independent component:
86 double L_ij = 2*L_i[i]*L_i[j]/(L_i[i] + L_i[j] + Tiny);
87 Lprime_m += molefracs[i]*molefracs[j]*L_ij;
88 // Additional variables for density-dependent component:
89 double f_ij = sqrt(f_i[i]*f_i[j]);
90 double h_ij = 0.125*pow(pow(h_i[i],1./3.) + pow(h_i[j],1./3.),3.);
91 double mw_ij_inv = (m_mw[i] + m_mw[j])/(2*m_mw[i]*m_mw[j]);
92 f_m += molefracs[i]*molefracs[j]*f_ij*h_ij;
93 h_m += molefracs[i]*molefracs[j]*h_ij;
94 mw_m += molefracs[i]*molefracs[j]*sqrt(mw_ij_inv*f_ij)*pow(h_ij,-4./3.);
95 }
96 }
97
98 f_m = f_m/h_m;
99 mw_m = pow(mw_m,-2.)*f_m*pow(h_m,-8./3.);
100
101 double rho_0 = 16.04*h_m/(1000*m_thermo->molarVolume());
102 double T_0 = m_temp/f_m;
103 double mu_0 = 1e-7*(2.90774e6/T_0 - 3.31287e6*pow(T_0,-2./3.)
104 + 1.60810e6*pow(T_0,-1./3.) - 4.33190e5 + 7.06248e4
105 *pow(T_0,1./3.) - 7.11662e3*pow(T_0,2./3.) + 4.32517e2*T_0
106 - 1.44591e1*pow(T_0,4./3.) + 2.03712e-1*pow(T_0,5./3.));
107 double L_1m = 1944*mu_0;
108 double L_2m = (-2.5276e-4 + 3.3433e-4*pow(1.12 - log(T_0/1.680e2),2))*rho_0;
109 double L_3m = exp(-7.19771 + 85.67822/T_0)*(exp((12.47183
110 - 984.6252*pow(T_0,-1.5))*pow(rho_0,0.1) + (rho_0/0.1617 - 1)
111 *sqrt(rho_0)*(0.3594685 + 69.79841/T_0 - 872.8833*pow(T_0,-2))) - 1.)*1e-3;
112 double H_m = sqrt(f_m*16.04/mw_m)*pow(h_m,-2./3.);
113 double Lstar_m = H_m*(L_1m + L_2m + L_3m);
114 return Lprime_m + Lstar_m;
115}
116
118{
119 // Method for MultiTransport class:
120 // solveLMatrixEquation();
121 // const double c = 1.6/GasConstant;
122 // for (size_t k = 0; k < m_nsp; k++) {
123 // dt[k] = c * m_mw[k] * m_molefracs[k] * m_a[k];
124 // }
125 throw NotImplementedError("HighPressureGasTransport::getThermalDiffCoeffs");
126}
127
128void HighPressureGasTransport::getBinaryDiffCoeffs(const size_t ld, double* const d)
129{
130 vector<double> PcP(5);
131 size_t nsp = m_thermo->nSpecies();
132 vector<double> molefracs(nsp);
133 m_thermo->getMoleFractions(&molefracs[0]);
134
135 update_T();
136 // Evaluate the binary diffusion coefficients from the polynomial fits.
137 // This should perhaps be preceded by a check to see whether any of T, P, or
138 // C have changed.
139 //if (!m_bindiff_ok) {
140 updateDiff_T();
141 //}
142 if (ld < nsp) {
143 throw CanteraError("HighPressureGasTransport::getBinaryDiffCoeffs",
144 "ld is too small");
145 }
146 double rp = 1.0/m_thermo->pressure();
147 for (size_t i = 0; i < nsp; i++) {
148 for (size_t j = 0; j < nsp; j++) {
149 // Add an offset to avoid a condition where x_i and x_j both equal
150 // zero (this would lead to Pr_ij = Inf):
151 double x_i = std::max(Tiny, molefracs[i]);
152 double x_j = std::max(Tiny, molefracs[j]);
153
154 // Weight mole fractions of i and j so that X_i + X_j = 1.0:
155 x_i = x_i/(x_i + x_j);
156 x_j = x_j/(x_i + x_j);
157
158 //Calculate Tr and Pr based on mole-fraction-weighted crit constants:
159 double Tr_ij = m_temp/(x_i*Tcrit_i(i) + x_j*Tcrit_i(j));
160 double Pr_ij = m_thermo->pressure()/(x_i*Pcrit_i(i) + x_j*Pcrit_i(j));
161
162 double P_corr_ij;
163 if (Pr_ij < 0.1) {
164 // If pressure is low enough, no correction is needed:
165 P_corr_ij = 1;
166 }else {
167 // Otherwise, calculate the parameters for Takahashi correlation
168 // by interpolating on Pr_ij:
169 P_corr_ij = setPcorr(Pr_ij, Tr_ij);
170
171 // If the reduced temperature is too low, the correction factor
172 // P_corr_ij will be < 0:
173 if (P_corr_ij<0) {
174 P_corr_ij = Tiny;
175 }
176 }
177
178 // Multiply the standard low-pressure binary diffusion coefficient
179 // (m_bdiff) by the Takahashi correction factor P_corr_ij:
180 d[ld*j + i] = P_corr_ij*rp * m_bdiff(i,j);
181 }
182 }
183}
184
185void HighPressureGasTransport::getMultiDiffCoeffs(const size_t ld, double* const d)
186{
187 // Not currently implemented. m_Lmatrix inversion returns NaN. Needs to be
188 // fixed. --SCD - 2-28-2014
189 throw NotImplementedError("HighPressureGasTransport:getMultiDiffCoeffs");
190 // Calculate the multi-component Stefan-Maxwell diffusion coefficients,
191 // based on the Takahashi-correlation-corrected binary diffusion coefficients.
192
193 // update the mole fractions
194 update_C();
195
196 // update the binary diffusion coefficients
197 update_T();
199
200 // Correct the binary diffusion coefficients for high-pressure effects; this
201 // is basically the same routine used in 'getBinaryDiffCoeffs,' above:
202 size_t nsp = m_thermo->nSpecies();
203 vector<double> molefracs(nsp);
204 m_thermo->getMoleFractions(&molefracs[0]);
205 update_T();
206 // Evaluate the binary diffusion coefficients from the polynomial fits -
207 // this should perhaps be preceded by a check for changes in T, P, or C.
208 updateDiff_T();
209
210 if (ld < m_nsp) {
211 throw CanteraError("HighPressureGasTransport::getMultiDiffCoeffs",
212 "ld is too small");
213 }
214 for (size_t i = 0; i < m_nsp; i++) {
215 for (size_t j = 0; j < m_nsp; j++) {
216 // Add an offset to avoid a condition where x_i and x_j both equal
217 // zero (this would lead to Pr_ij = Inf):
218 double x_i = std::max(Tiny, molefracs[i]);
219 double x_j = std::max(Tiny, molefracs[j]);
220 x_i = x_i/(x_i+x_j);
221 x_j = x_j/(x_i+x_j);
222 double Tr_ij = m_temp/(x_i*Tcrit_i(i) + x_j*Tcrit_i(j));
223 double Pr_ij = m_thermo->pressure()/(x_i*Pcrit_i(i) + x_j*Pcrit_i(j));
224
225 double P_corr_ij;
226 if (Pr_ij < 0.1) {
227 P_corr_ij = 1;
228 }else {
229 P_corr_ij = setPcorr(Pr_ij, Tr_ij);
230 if (P_corr_ij<0) {
231 P_corr_ij = Tiny;
232 }
233 }
234
235 m_bdiff(i,j) *= P_corr_ij;
236 }
237 }
238 m_bindiff_ok = false; // m_bdiff is overwritten by the above routine.
239
240 // Having corrected m_bdiff for pressure and concentration effects, the
241 // routine now proceeds the same as in the low-pressure case:
242
243 // evaluate L0000 if the temperature or concentrations have
244 // changed since it was last evaluated.
245 if (!m_l0000_ok) {
246 eval_L0000(molefracs.data());
247 }
248
249 // invert L00,00
250 int ierr = invert(m_Lmatrix, m_nsp);
251 if (ierr != 0) {
252 throw CanteraError("HighPressureGasTransport::getMultiDiffCoeffs",
253 "invert returned ierr = {}", ierr);
254 }
255 m_l0000_ok = false; // matrix is overwritten by inverse
256 m_lmatrix_soln_ok = false;
257
258 double prefactor = 16.0 * m_temp
260
261 for (size_t i = 0; i < m_nsp; i++) {
262 for (size_t j = 0; j < m_nsp; j++) {
263 double c = prefactor/m_mw[j];
264 d[ld*j + i] = c*molefracs[i]*(m_Lmatrix(i,j) - m_Lmatrix(i,i));
265 }
266 }
267}
268
270{
271 // Calculate the high-pressure mixture viscosity, based on the Lucas method.
272 double Tc_mix = 0.;
273 double Pc_mix_n = 0.;
274 double Pc_mix_d = 0.;
275 double MW_mix = m_thermo->meanMolecularWeight();
276 double MW_H = m_mw[0];
277 double MW_L = m_mw[0];
278 double FP_mix_o = 0;
279 double FQ_mix_o = 0;
280 double tKelvin = m_thermo->temperature();
281 double Pvp_mix = m_thermo->satPressure(tKelvin);
282 size_t nsp = m_thermo->nSpecies();
283 vector<double> molefracs(nsp);
284 m_thermo->getMoleFractions(&molefracs[0]);
285
286 double x_H = molefracs[0];
287 for (size_t i = 0; i < m_nsp; i++) {
288 // Calculate pure-species critical constants and add their contribution
289 // to the mole-fraction-weighted mixture averages:
290 double Tc = Tcrit_i(i);
291 double Tr = tKelvin/Tc;
292 double Zc = Zcrit_i(i);
293 Tc_mix += Tc*molefracs[i];
294 Pc_mix_n += molefracs[i]*Zc; //numerator
295 Pc_mix_d += molefracs[i]*Vcrit_i(i); //denominator
296
297 // Need to calculate ratio of heaviest to lightest species:
298 if (m_mw[i] > MW_H) {
299 MW_H = m_mw[i];
300 x_H = molefracs[i];
301 } else if (m_mw[i] < MW_L) {
302 MW_L = m_mw[i]; }
303
304 // Calculate reduced dipole moment for polar correction term:
305 double mu_ri = 52.46*100000*m_dipole(i,i)*m_dipole(i,i)
306 *Pcrit_i(i)/(Tc*Tc);
307 if (mu_ri < 0.022) {
308 FP_mix_o += molefracs[i];
309 } else if (mu_ri < 0.075) {
310 FP_mix_o += molefracs[i]*(1. + 30.55*pow(0.292 - Zc, 1.72));
311 } else { FP_mix_o += molefracs[i]*(1. + 30.55*pow(0.292 - Zc, 1.72)
312 *fabs(0.96 + 0.1*(Tr - 0.7)));
313 }
314
315 // Calculate contribution to quantum correction term.
316 // SCD Note: This assumes the species of interest (He, H2, and D2) have
317 // been named in this specific way. They are perhaps the most obvious
318 // names, but it would of course be preferred to have a more general
319 // approach, here.
320 vector<string> spnames = m_thermo->speciesNames();
321 if (spnames[i] == "He") {
322 FQ_mix_o += molefracs[i]*FQ_i(1.38,Tr,m_mw[i]);
323 } else if (spnames[i] == "H2") {
324 FQ_mix_o += molefracs[i]*(FQ_i(0.76,Tr,m_mw[i]));
325 } else if (spnames[i] == "D2") {
326 FQ_mix_o += molefracs[i]*(FQ_i(0.52,Tr,m_mw[i]));
327 } else {
328 FQ_mix_o += molefracs[i];
329 }
330 }
331
332 double Tr_mix = tKelvin/Tc_mix;
333 double Pc_mix = GasConstant*Tc_mix*Pc_mix_n/Pc_mix_d;
334 double Pr_mix = m_thermo->pressure()/Pc_mix;
335 double ratio = MW_H/MW_L;
336 double ksi = pow(GasConstant*Tc_mix*3.6277*pow(10.0,53.0)/(pow(MW_mix,3)
337 *pow(Pc_mix,4)),1.0/6.0);
338
339 if (ratio > 9 && x_H > 0.05 && x_H < 0.7) {
340 FQ_mix_o *= 1 - 0.01*pow(ratio,0.87);
341 }
342
343 // Calculate Z1m
344 double Z1m = (0.807*pow(Tr_mix,0.618) - 0.357*exp(-0.449*Tr_mix)
345 + 0.340*exp(-4.058*Tr_mix)+0.018)*FP_mix_o*FQ_mix_o;
346
347 // Calculate Z2m:
348 double Z2m;
349 if (Tr_mix <= 1.0) {
350 if (Pr_mix < Pvp_mix/Pc_mix) {
351 double alpha = 3.262 + 14.98*pow(Pr_mix,5.508);
352 double beta = 1.390 + 5.746*Pr_mix;
353 Z2m = 0.600 + 0.760*pow(Pr_mix,alpha) + (0.6990*pow(Pr_mix,beta) -
354 0.60)*(1- Tr_mix);
355 } else {
356 throw CanteraError("HighPressureGasTransport::viscosity",
357 "State is outside the limits of the Lucas model, Tr <= 1");
358 }
359 } else if ((Tr_mix > 1.0) && (Tr_mix < 40.0)) {
360 if ((Pr_mix > 0.0) && (Pr_mix <= 100.0)) {
361 double a_fac = 0.001245*exp(5.1726*pow(Tr_mix,-0.3286))/Tr_mix;
362 double b_fac = a_fac*(1.6553*Tr_mix - 1.2723);
363 double c_fac = 0.4489*exp(3.0578*pow(Tr_mix,-37.7332))/Tr_mix;
364 double d_fac = 1.7368*exp(2.2310*pow(Tr_mix,-7.6351))/Tr_mix;
365 double f_fac = 0.9425*exp(-0.1853*pow(Tr_mix,0.4489));
366
367 Z2m = Z1m*(1 + a_fac*pow(Pr_mix,1.3088)/(b_fac*pow(Pr_mix,f_fac)
368 + pow(1+c_fac*pow(Pr_mix,d_fac),-1)));
369 } else {
370 throw CanteraError("HighPressureGasTransport::viscosity",
371 "State is outside the limits of the Lucas model, 1.0 < Tr < 40");
372 }
373 } else {
374 throw CanteraError("HighPressureGasTransport::viscosity",
375 "State is outside the limits of the Lucas model, Tr > 40");
376 }
377
378 // Calculate Y:
379 double Y = Z2m/Z1m;
380
381 // Return the viscosity:
382 return Z2m*(1 + (FP_mix_o - 1)*pow(Y,-3))*(1 + (FQ_mix_o - 1)
383 *(1/Y - 0.007*pow(log(Y),4)))/(ksi*FP_mix_o*FQ_mix_o);
384}
385
386// Pure species critical properties - Tc, Pc, Vc, Zc:
387double HighPressureGasTransport::Tcrit_i(size_t i)
388{
389 // Store current molefracs and set temp molefrac of species i to 1.0:
390 vector<double> molefracs = store(i, m_thermo->nSpecies());
391
392 double tc = m_thermo->critTemperature();
393 // Restore actual molefracs:
394 m_thermo->setMoleFractions(&molefracs[0]);
395 return tc;
396}
397
398double HighPressureGasTransport::Pcrit_i(size_t i)
399{
400 // Store current molefracs and set temp molefrac of species i to 1.0:
401 vector<double> molefracs = store(i, m_thermo->nSpecies());
402
403 double pc = m_thermo->critPressure();
404 // Restore actual molefracs:
405 m_thermo->setMoleFractions(&molefracs[0]);
406 return pc;
407}
408
409double HighPressureGasTransport::Vcrit_i(size_t i)
410{
411 // Store current molefracs and set temp molefrac of species i to 1.0:
412 vector<double> molefracs = store(i, m_thermo->nSpecies());
413
414 double vc = m_thermo->critVolume();
415 // Restore actual molefracs:
416 m_thermo->setMoleFractions(&molefracs[0]);
417 return vc;
418}
419
420double HighPressureGasTransport::Zcrit_i(size_t i)
421{
422 // Store current molefracs and set temp molefrac of species i to 1.0:
423 vector<double> molefracs = store(i, m_thermo->nSpecies());
424
425 double zc = m_thermo->critCompressibility();
426 // Restore actual molefracs:
427 m_thermo->setMoleFractions(&molefracs[0]);
428 return zc;
429}
430
431vector<double> HighPressureGasTransport::store(size_t i, size_t nsp)
432{
433 vector<double> molefracs(nsp);
434 m_thermo->getMoleFractions(&molefracs[0]);
435 vector<double> mf_temp(nsp, 0.0);
436 mf_temp[i] = 1;
437 m_thermo->setMoleFractions(&mf_temp[0]);
438 return molefracs;
439}
440
441// Calculates quantum correction term for a species based on Tr and MW, used in
442// viscosity calculation:
443double HighPressureGasTransport::FQ_i(double Q, double Tr, double MW)
444{
445 return 1.22*pow(Q,0.15)*(1 + 0.00385*pow(pow(Tr - 12.,2.),1./MW)
446 *fabs(Tr-12)/(Tr-12));
447}
448
449// Set value of parameter values for Takahashi correlation, by interpolating
450// table of constants vs. Pr:
451double HighPressureGasTransport::setPcorr(double Pr, double Tr)
452{
453 const static double Pr_lookup[17] = {0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.8, 1.0,
454 1.2, 1.4, 1.6, 1.8, 2.0, 2.5, 3.0, 4.0, 5.0};
455 const static double DP_Rt_lookup[17] = {1.01, 1.01, 1.01, 1.01, 1.01, 1.01,
456 1.01, 1.02, 1.02, 1.02, 1.02, 1.03, 1.03, 1.04, 1.05, 1.06, 1.07};
457 const static double A_ij_lookup[17] = {0.038042, 0.067433, 0.098317,
458 0.137610, 0.175081, 0.216376, 0.314051, 0.385736, 0.514553, 0.599184,
459 0.557725, 0.593007, 0.696001, 0.790770, 0.502100, 0.837452, 0.890390};
460 const static double B_ij_lookup[17] = {1.52267, 2.16794, 2.42910, 2.77605,
461 2.98256, 3.11384, 3.50264, 3.07773, 3.54744, 3.61216, 3.41882, 3.18415,
462 3.37660, 3.27984, 3.39031, 3.23513, 3.13001};
463 const static double C_ij_lookup[17] = {0., 0., 0., 0., 0., 0., 0., 0.141211,
464 0.278407, 0.372683, 0.504894, 0.678469, 0.665702, 0., 0.602907, 0., 0.};
465 const static double E_ij_lookup[17] = {1., 1., 1., 1., 1., 1., 1., 13.45454,
466 14., 10.00900, 8.57519, 10.37483, 11.21674, 1., 6.19043, 1., 1.};
467
468 // Interpolate Pr vs. those used in Takahashi table:
469 int Pr_i = 0;
470 double frac = 0.;
471
472 if (Pr < 0.1) {
473 frac = (Pr - Pr_lookup[0])/(Pr_lookup[1] - Pr_lookup[0]);
474 } else {
475 for (int j = 1; j < 17; j++) {
476 if (Pr_lookup[j] > Pr) {
477 frac = (Pr - Pr_lookup[j-1])/(Pr_lookup[j] - Pr_lookup[j-1]);
478 break;
479 }
480 Pr_i++;
481 }
482 }
483 // If Pr is greater than the greatest value used by Takahashi (5.0), use the
484 // final table value. Should eventually add in an extrapolation:
485 if (Pr_i == 17) {
486 frac = 1.0;
487 }
488
489 double P_corr_1 = DP_Rt_lookup[Pr_i]*(1.0 - A_ij_lookup[Pr_i]
490 *pow(Tr,-B_ij_lookup[Pr_i]))*(1-C_ij_lookup[Pr_i]
491 *pow(Tr,-E_ij_lookup[Pr_i]));
492 double P_corr_2 = DP_Rt_lookup[Pr_i+1]*(1.0 - A_ij_lookup[Pr_i+1]
493 *pow(Tr,-B_ij_lookup[Pr_i+1]))*(1-C_ij_lookup[Pr_i+1]
494 *pow(Tr,-E_ij_lookup[Pr_i+1]));
495 return P_corr_1*(1.0-frac) + P_corr_2*frac;
496}
497
498}
Interface for class HighPressureGasTransport.
ThermoPhase object for the ideal gas equation of state - workhorse for Cantera (see Thermodynamic Pro...
Interface for class MultiTransport.
Header file defining class TransportFactory (see TransportFactory)
Base class for exceptions thrown by Cantera classes.
vector< double > m_mw
Local copy of the species molecular weights.
double m_temp
Current value of the temperature at which the properties in this object are calculated (Kelvin).
virtual void updateDiff_T()
Update the binary diffusion coefficients.
bool m_bindiff_ok
Update boolean for the binary diffusivities at unit pressure.
DenseMatrix m_bdiff
Matrix of binary diffusion coefficients at the reference pressure and the current temperature Size is...
DenseMatrix m_dipole
The effective dipole moment for (i,j) collisions.
vector< double > m_w_ac
Pitzer acentric factor.
void getBinaryDiffCoeffs(const size_t ld, double *const d) override
Returns the matrix of binary diffusion coefficients.
double thermalConductivity() override
Returns the mixture thermal conductivity in W/m/K.
HighPressureGasTransport(ThermoPhase *thermo=0)
default constructor
double viscosity() override
Viscosity of the mixture (kg /m /s)
void getThermalDiffCoeffs(double *const dt) override
Return the thermal diffusion coefficients (kg/m/s)
void getMultiDiffCoeffs(const size_t ld, double *const d) override
Return the Multicomponent diffusion coefficients. Units: [m^2/s].
Class MultiTransport implements multicomponent transport properties for ideal gas mixtures.
bool m_l0000_ok
Boolean indicating viscosity is up to date.
void update_T() override
Update basic temperature-dependent quantities if the temperature has changed.
void eval_L0000(const double *const x)
Evaluate the L0000 matrices.
void updateThermal_T()
Update the temperature-dependent terms needed to compute the thermal conductivity and thermal diffusi...
void update_C() override
Update basic concentration-dependent quantities if the concentrations have changed.
An error indicating that an unimplemented function has been called.
virtual void setMoleFractions(const double *const x)
Set the mole fractions to the specified values.
Definition Phase.cpp:304
size_t nSpecies() const
Returns the number of species in the phase.
Definition Phase.h:245
double temperature() const
Temperature (K).
Definition Phase.h:662
double meanMolecularWeight() const
The mean molecular weight. Units: (kg/kmol)
Definition Phase.h:760
void getMoleFractions(double *const x) const
Get the species mole fraction vector.
Definition Phase.cpp:540
const vector< string > & speciesNames() const
Return a const reference to the vector of species names.
Definition Phase.cpp:157
virtual double molarVolume() const
Molar volume (m^3/kmol).
Definition Phase.cpp:702
virtual double pressure() const
Return the thermodynamic pressure (Pa).
Definition Phase.h:680
Base class for a phase with thermodynamic properties.
virtual double critTemperature() const
Critical temperature (K).
virtual void getCp_R_ref(double *cprt) const
Returns the vector of nondimensional constant pressure heat capacities of the reference state at the ...
virtual double critPressure() const
Critical pressure (Pa).
virtual double critVolume() const
Critical volume (m3/kmol).
virtual double satPressure(double t)
Return the saturation pressure given the temperature.
virtual double critCompressibility() const
Critical compressibility (unitless).
ThermoPhase * m_thermo
pointer to the object representing the phase
Definition Transport.h:821
size_t m_nsp
Number of species.
Definition Transport.h:828
const double GasConstant
Universal Gas Constant [J/kmol/K].
Definition ct_defs.h:120
Namespace for the Cantera kernel.
Definition AnyMap.cpp:564
const double Tiny
Small number to compare differences of mole fractions against.
Definition ct_defs.h:173
int invert(DenseMatrix &A, size_t nn)
invert A. A is overwritten with A^-1.
const double BigNumber
largest number to compare to inf.
Definition ct_defs.h:160
Contains declarations for string manipulation functions within Cantera.
Various templated functions that carry out common vector and polynomial operations (see Templated Arr...