43 typedef ::Opm::IdealGas<Scalar>
IdealGas;
44 static const int liquidPhaseIdx = 0;
45 static const int gasPhaseIdx = 1;
55 template <
class Evaluation>
56 static Evaluation
gasDiffCoeff(
const Evaluation& temperature,
const Evaluation& pressure,
bool extrapolate =
false)
59 Scalar k = 1.3806504e-23;
61 Scalar R_h = 1.72e-10;
63 return k / (c * M_PI * R_h) * (temperature / mu);
72 template <
class Evaluation>
96 template <
class Evaluation>
99 const Evaluation& salinity,
100 const int knownPhaseIdx,
103 bool extrapolate =
false)
105 OPM_TIMEFUNCTION_LOCAL();
106 Evaluation A = computeA_(temperature, pg, extrapolate);
109 Evaluation x_NaCl = salinityToMolFrac_(salinity);
113 if (knownPhaseIdx < 0) {
114 Evaluation molalityNaCl = moleFracToMolality_(x_NaCl);
115 Evaluation m0_CO2 = molalityCO2inPureWater_(temperature, pg, extrapolate);
116 Evaluation gammaStar = activityCoefficient_(temperature, pg, molalityNaCl);
117 Evaluation m_CO2 = m0_CO2 / gammaStar;
118 xlCO2 = m_CO2 / (molalityNaCl + 55.508 + m_CO2);
119 ygH2O = A * (1 - xlCO2 - x_NaCl);
125 if (knownPhaseIdx == liquidPhaseIdx)
126 ygH2O = A * (1 - xlCO2 - x_NaCl);
131 if (knownPhaseIdx == gasPhaseIdx)
133 xlCO2 = 1 - x_NaCl - ygH2O / A;
139 template <
class Evaluation>
140 static Evaluation
henry(
const Evaluation& temperature,
bool extrapolate =
false)
151 template <
class Evaluation>
154 OPM_TIMEFUNCTION_LOCAL();
155 Valgrind::CheckDefined(temperature);
156 Valgrind::CheckDefined(pg);
159 Evaluation pg_bar = pg / 1.e5;
160 Evaluation a_CO2 = (7.54e7 - 4.13e4 * temperature);
165 lnPhiCO2 = log(V / (V - b_CO2));
166 lnPhiCO2 += b_CO2 / (V - b_CO2);
167 lnPhiCO2 -= 2 * a_CO2 / (R * pow(temperature, 1.5) * b_CO2) * log((V + b_CO2) / V);
171 * pow(temperature, 1.5)
174 * (log((V + b_CO2) / V)
175 - b_CO2 / (V + b_CO2));
176 lnPhiCO2 -= log(pg_bar * V / (R * temperature));
178 return exp(lnPhiCO2);
189 template <
class Evaluation>
192 OPM_TIMEFUNCTION_LOCAL();
194 const Evaluation& pg_bar = pg / 1.e5;
195 const Evaluation& a_CO2 = (7.54e7 - 4.13e4 * temperature);
196 Scalar a_CO2_H2O = 7.89e7;
198 Scalar b_H2O = 18.18;
204 + b_H2O/(V - b_CO2) - 2*a_CO2_H2O
205 / (R*pow(temperature, 1.5)*b_CO2)*log((V + b_CO2)/V)
206 + a_CO2*b_H2O/(R*pow(temperature, 1.5)*b_CO2*b_CO2)
207 *(log((V + b_CO2)/V) - b_CO2/(V + b_CO2))
208 - log(pg_bar*V/(R*temperature));
209 return exp(lnPhiH2O);
218 template <
class Evaluation>
219 static Evaluation salinityToMolFrac_(
const Evaluation& salinity) {
220 OPM_TIMEFUNCTION_LOCAL();
222 const Scalar Ms = 58.44e-3;
224 const Evaluation X_NaCl = salinity;
226 const Evaluation x_NaCl = -Mw * X_NaCl / ((Ms - Mw) * X_NaCl - Ms);
235 template <
class Evaluation>
236 static Evaluation moleFracToMolality_(
const Evaluation& x_NaCl)
239 return 55.508 * x_NaCl / (1 - x_NaCl);
249 template <
class Evaluation>
250 static Evaluation molalityCO2inPureWater_(
const Evaluation& temperature,
const Evaluation& pg,
bool extrapolate =
false)
252 OPM_TIMEFUNCTION_LOCAL();
253 const Evaluation& A = computeA_(temperature, pg, extrapolate);
254 const Evaluation& B = computeB_(temperature, pg, extrapolate);
255 const Evaluation& yH2OinGas = (1 - B) / (1. / A - B);
256 const Evaluation& xCO2inWater = B * (1 - yH2OinGas);
257 return (xCO2inWater * 55.508) / (1 - xCO2inWater);
269 template <
class Evaluation>
270 static Evaluation activityCoefficient_(
const Evaluation& temperature,
271 const Evaluation& pg,
272 const Evaluation& molalityNaCl)
274 OPM_TIMEFUNCTION_LOCAL();
275 const Evaluation& lambda = computeLambda_(temperature, pg);
276 const Evaluation& xi = computeXi_(temperature, pg);
277 const Evaluation& lnGammaStar =
278 2*molalityNaCl*lambda + xi*molalityNaCl*molalityNaCl;
279 return exp(lnGammaStar);
290 template <
class Evaluation>
291 static Evaluation computeA_(
const Evaluation& temperature,
const Evaluation& pg,
bool extrapolate =
false)
293 OPM_TIMEFUNCTION_LOCAL();
294 const Evaluation& deltaP = pg / 1e5 - 1;
295 Scalar v_av_H2O = 18.1;
297 const Evaluation& k0_H2O = equilibriumConstantH2O_(temperature);
299 const Evaluation& pg_bar = pg / 1.e5;
300 return k0_H2O/(phi_H2O*pg_bar)*exp(deltaP*v_av_H2O/(R*temperature));
311 template <
class Evaluation>
312 static Evaluation computeB_(
const Evaluation& temperature,
const Evaluation& pg,
bool extrapolate =
false)
314 OPM_TIMEFUNCTION_LOCAL();
315 const Evaluation& deltaP = pg / 1e5 - 1;
316 const Scalar v_av_CO2 = 32.6;
318 const Evaluation& k0_CO2 = equilibriumConstantCO2_(temperature);
320 const Evaluation& pg_bar = pg / 1.e5;
321 return phi_CO2*pg_bar/(55.508*k0_CO2)*exp(-(deltaP*v_av_CO2)/(R*temperature));
331 template <
class Evaluation>
332 static Evaluation computeLambda_(
const Evaluation& temperature,
const Evaluation& pg)
334 OPM_TIMEFUNCTION_LOCAL();
335 static const Scalar c[6] =
336 { -0.411370585, 6.07632013E-4, 97.5347708, -0.0237622469, 0.0170656236, 1.41335834E-5 };
338 Evaluation pg_bar = pg / 1.0E5;
343 + c[3]*pg_bar/temperature
344 + c[4]*pg_bar/(630.0 - temperature)
345 + c[5]*temperature*log(pg_bar);
355 template <
class Evaluation>
356 static Evaluation computeXi_(
const Evaluation& temperature,
const Evaluation& pg)
358 OPM_TIMEFUNCTION_LOCAL();
359 static const Scalar c[4] =
360 { 3.36389723E-4, -1.98298980E-5, 2.12220830E-3, -5.24873303E-3 };
362 Evaluation pg_bar = pg / 1.0E5;
363 return c[0] + c[1]*temperature + c[2]*pg_bar/temperature + c[3]*pg_bar/(630.0 - temperature);
372 template <
class Evaluation>
373 static Evaluation equilibriumConstantCO2_(
const Evaluation& temperature)
375 OPM_TIMEFUNCTION_LOCAL();
376 Evaluation temperatureCelcius = temperature - 273.15;
377 static const Scalar c[3] = { 1.189, 1.304e-2, -5.446e-5 };
378 Evaluation logk0_CO2 = c[0] + temperatureCelcius*(c[1] + temperatureCelcius*c[2]);
379 Evaluation k0_CO2 = pow(10.0, logk0_CO2);
389 template <
class Evaluation>
390 static Evaluation equilibriumConstantH2O_(
const Evaluation& temperature)
392 OPM_TIMEFUNCTION_LOCAL();
393 Evaluation temperatureCelcius = temperature - 273.15;
394 static const Scalar c[4] = { -2.209, 3.097e-2, -1.098e-4, 2.048e-7 };
395 Evaluation logk0_H2O =
396 c[0] + temperatureCelcius*(c[1] + temperatureCelcius*(c[2] + temperatureCelcius*c[3]));
397 return pow(10.0, logk0_H2O);
static void calculateMoleFractions(const Evaluation &temperature, const Evaluation &pg, const Evaluation &salinity, const int knownPhaseIdx, Evaluation &xlCO2, Evaluation &ygH2O, bool extrapolate=false)
Returns the mol (!) fraction of CO2 in the liquid phase and the mol_ (!) fraction of H2O in the gas p...
Definition Brine_CO2.hpp:97