Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #include "BrineFluidProperties.h"
11 :
12 : registerMooseObject("FluidPropertiesApp", BrineFluidProperties);
13 :
14 : InputParameters
15 98 : BrineFluidProperties::validParams()
16 : {
17 98 : InputParameters params = MultiComponentFluidProperties::validParams();
18 196 : params.addParam<UserObjectName>("water_fp",
19 : "The name of the FluidProperties UserObject for water");
20 98 : params.addClassDescription("Fluid properties for brine");
21 98 : return params;
22 0 : }
23 :
24 52 : BrineFluidProperties::BrineFluidProperties(const InputParameters & parameters)
25 52 : : MultiComponentFluidProperties(parameters), _water_fp_derivs(true)
26 : {
27 : // There are two possibilities to consider:
28 : // 1) No water_fp has been supplied (in which case one is constructed)
29 : // 2) A water_fp hase been supplied (in which case it is used)
30 : // In both cases, though, a Water97FluidProperties UserObject must be added
31 : // Note: this UserObject is only used to gain access to the Henry's constant
32 : // formulation. All property calculations are performed using _water_fp
33 52 : const std::string water_name = name() + ":water";
34 : {
35 52 : const std::string class_name = "Water97FluidProperties";
36 52 : InputParameters params = _app.getFactory().getValidParams(class_name);
37 52 : if (_tid == 0)
38 46 : _fe_problem.addUserObject(class_name, water_name, params);
39 52 : }
40 52 : _water97_fp = &_fe_problem.getUserObject<Water97FluidProperties>(water_name);
41 :
42 104 : if (parameters.isParamSetByUser("water_fp"))
43 : {
44 : // SinglePhaseFluidPropertiesPT UserObject for water
45 22 : _water_fp = &getUserObject<SinglePhaseFluidProperties>("water_fp");
46 :
47 : // Check that a water userobject has actually been supplied
48 44 : if (_water_fp->fluidName() != "water")
49 0 : paramError("water_fp", "A water FluidProperties UserObject must be supplied");
50 : }
51 : else
52 : {
53 : // Construct a SinglePhaseFluidProperties UserObject for water
54 30 : _water_fp = &_fe_problem.getUserObject<SinglePhaseFluidProperties>(water_name);
55 : }
56 :
57 : // SinglePhaseFluidProperties UserObject for NaCl
58 52 : const std::string nacl_name = name() + ":nacl";
59 : {
60 52 : const std::string class_name = "NaClFluidProperties";
61 52 : InputParameters params = _app.getFactory().getValidParams(class_name);
62 52 : if (_tid == 0)
63 46 : _fe_problem.addUserObject(class_name, nacl_name, params);
64 52 : }
65 52 : _nacl_fp = &_fe_problem.getUserObject<SinglePhaseFluidProperties>(nacl_name);
66 :
67 : // Molar mass of NaCl and H20
68 52 : _Mnacl = _nacl_fp->molarMass();
69 52 : _Mh2o = _water_fp->molarMass();
70 52 : }
71 :
72 52 : BrineFluidProperties::~BrineFluidProperties() {}
73 :
74 : const SinglePhaseFluidProperties &
75 10 : BrineFluidProperties::getComponent(unsigned int component) const
76 : {
77 10 : switch (component)
78 : {
79 9 : case WATER:
80 9 : return *_water_fp;
81 :
82 1 : case NACL:
83 1 : return *_nacl_fp;
84 :
85 0 : default:
86 0 : mooseError("BrineFluidProperties::getComponent has been provided an incorrect component");
87 : }
88 : }
89 :
90 : std::string
91 1 : BrineFluidProperties::fluidName() const
92 : {
93 1 : return "brine";
94 : }
95 :
96 : FPADReal
97 2679 : BrineFluidProperties::molarMass(const FPADReal & xnacl) const
98 : {
99 8037 : return 1.0 / (xnacl / _Mnacl + (1.0 - xnacl) / _Mh2o);
100 : }
101 :
102 : Real
103 436 : BrineFluidProperties::molarMass(Real xnacl) const
104 : {
105 436 : return 1.0 / (xnacl / _Mnacl + (1.0 - xnacl) / _Mh2o);
106 : }
107 :
108 : Real
109 2 : BrineFluidProperties::molarMassNaCl() const
110 : {
111 2 : return _Mnacl;
112 : }
113 :
114 : Real
115 2 : BrineFluidProperties::molarMassH2O() const
116 : {
117 2 : return _Mh2o;
118 : }
119 :
120 : FPADReal
121 895 : BrineFluidProperties::rho_from_p_T_X(const FPADReal & pressure,
122 : const FPADReal & temperature,
123 : const FPADReal & xnacl) const
124 : {
125 : using std::exp, std::sqrt;
126 : // The correlation requires the pressure in bar, not Pa.
127 895 : FPADReal pbar = pressure * 1.0e-5;
128 895 : FPADReal pbar2 = pbar * pbar;
129 895 : FPADReal pbar3 = pbar2 * pbar;
130 :
131 : // The correlation requires mole fraction
132 895 : const FPADReal Xnacl = massFractionToMoleFraction(xnacl);
133 :
134 895 : const FPADReal n11 = -54.2958 - 45.7623 * exp(-9.44785e-4 * pbar);
135 895 : const FPADReal n21 = -2.6142 - 2.39092e-4 * pbar;
136 895 : const FPADReal n22 = 0.0356828 + 4.37235e-6 * pbar + 2.0566e-9 * pbar2;
137 : const FPADReal n1x1 =
138 895 : 330.47 + 0.942876 * sqrt(pbar) + 0.0817193 * pbar - 2.47556e-8 * pbar2 + 3.45052e-10 * pbar3;
139 895 : const FPADReal n2x1 = -0.0370751 + 0.00237723 * sqrt(pbar) + 5.42049e-5 * pbar +
140 1790 : 5.84709e-9 * pbar2 - 5.99373e-13 * pbar3;
141 895 : const FPADReal n12 = -n1x1 - n11;
142 895 : const FPADReal n20 = 1.0 - n21 * sqrt(n22);
143 895 : const FPADReal n23 = n2x1 - n20 - n21 * sqrt(1.0 + n22);
144 :
145 : // The temperature Tv where the brine has the same molar volume as pure water
146 : // Note: correlation uses temperature in Celcius
147 895 : const FPADReal n1 = n1x1 + n11 * (1.0 - Xnacl) + n12 * (1.0 - Xnacl) * (1.0 - Xnacl);
148 895 : const FPADReal n2 = n20 + n21 * sqrt(Xnacl + n22) + n23 * Xnacl;
149 895 : const FPADReal Tv = n1 + n2 * (temperature - _T_c2k);
150 :
151 : // The density of water at temperature Tv
152 : // Note: convert Tv to Kelvin to calculate water density
153 : FPADReal water_density;
154 895 : if (_water_fp_derivs)
155 : {
156 : Real rho, drho_dp, drho_dT;
157 6 : _water_fp->rho_from_p_T(pressure.value(), Tv.value() + _T_c2k, rho, drho_dp, drho_dT);
158 : water_density = rho;
159 :
160 12 : water_density.derivatives() = pressure.derivatives() * drho_dp + Tv.derivatives() * drho_dT;
161 : }
162 : else
163 889 : water_density = _water_fp->rho_from_p_T(pressure.value(), Tv.value() + _T_c2k);
164 :
165 : // The brine density is given by the water density scaled by the ratio of
166 : // brine molar mass to pure water molar mass
167 1790 : return water_density * molarMass(xnacl) / _Mh2o;
168 : }
169 :
170 : Real
171 889 : BrineFluidProperties::rho_from_p_T_X(Real pressure, Real temperature, Real xnacl) const
172 : {
173 : // Initialise the AD value (no derivatives required)
174 : FPADReal p = pressure;
175 : FPADReal T = temperature;
176 : FPADReal x = xnacl;
177 :
178 889 : _water_fp_derivs = false;
179 889 : FPADReal ad_rho = this->rho_from_p_T_X(p, T, x);
180 :
181 889 : return ad_rho.value();
182 : }
183 :
184 : void
185 4 : BrineFluidProperties::rho_from_p_T_X(Real pressure,
186 : Real temperature,
187 : Real xnacl,
188 : Real & rho,
189 : Real & drho_dp,
190 : Real & drho_dT,
191 : Real & drho_dx) const
192 : {
193 : // Initialise the AD value and derivatives
194 : FPADReal p = pressure;
195 : Moose::derivInsert(p.derivatives(), 0, 1.0);
196 : FPADReal T = temperature;
197 : Moose::derivInsert(T.derivatives(), 1, 1.0);
198 : FPADReal x = xnacl;
199 : Moose::derivInsert(x.derivatives(), 2, 1.0);
200 :
201 4 : _water_fp_derivs = true;
202 4 : FPADReal ad_rho = this->rho_from_p_T_X(p, T, x);
203 :
204 4 : rho = ad_rho.value();
205 4 : drho_dp = ad_rho.derivatives()[0];
206 4 : drho_dT = ad_rho.derivatives()[1];
207 4 : drho_dx = ad_rho.derivatives()[2];
208 4 : }
209 :
210 : Real
211 13 : BrineFluidProperties::mu_from_p_T_X(Real pressure, Real temperature, Real xnacl) const
212 : {
213 : // Correlation requires molal concentration (mol/kg)
214 13 : const Real mol = massFractionToMolalConc(xnacl);
215 13 : const Real mol2 = mol * mol;
216 13 : const Real mol3 = mol2 * mol;
217 :
218 : // Correlation requires temperature in C
219 13 : const Real Tc = temperature - _T_c2k;
220 :
221 13 : const Real a = 1.0 + 0.0816 * mol + 0.0122 * mol2 + 0.128e-3 * mol3 +
222 13 : 0.629e-3 * Tc * (1.0 - std::exp(-0.7 * mol));
223 :
224 13 : const Real water_viscosity = _water_fp->mu_from_p_T(pressure, temperature);
225 :
226 13 : return a * water_viscosity;
227 : }
228 :
229 : void
230 4 : BrineFluidProperties::mu_from_p_T_X(Real pressure,
231 : Real temperature,
232 : Real xnacl,
233 : Real & mu,
234 : Real & dmu_dp,
235 : Real & dmu_dT,
236 : Real & dmu_dx) const
237 : {
238 : // Viscosity of water and derivatives wrt pressure and temperature
239 : Real muw, dmuw_dp, dmuw_dT;
240 4 : _water_fp->mu_from_p_T(pressure, temperature, muw, dmuw_dp, dmuw_dT);
241 :
242 : // Correlation requires molal concentration (mol/kg)
243 4 : Real mol = massFractionToMolalConc(xnacl);
244 4 : Real dmol_dx = 1.0 / ((1.0 - xnacl) * (1.0 - xnacl) * _Mnacl);
245 4 : Real mol2 = mol * mol;
246 4 : Real mol3 = mol2 * mol;
247 :
248 : // Correlation requires temperature in C
249 4 : Real Tc = temperature - _T_c2k;
250 :
251 4 : Real a = 1.0 + 0.0816 * mol + 0.0122 * mol2 + 0.128e-3 * mol3 +
252 4 : 0.629e-3 * Tc * (1.0 - std::exp(-0.7 * mol));
253 : Real da_dx =
254 4 : (0.0816 + 0.0244 * mol + 3.84e-4 * mol2 + 4.403e-4 * Tc * std::exp(-0.7 * mol)) * dmol_dx;
255 4 : Real da_dT = 0.629e-3 * (1.0 - std::exp(-0.7 * mol));
256 :
257 4 : mu = a * muw;
258 4 : dmu_dp = a * dmuw_dp;
259 4 : dmu_dx = da_dx * muw;
260 4 : dmu_dT = da_dT * muw + a * dmuw_dT;
261 4 : }
262 :
263 : FPADReal
264 889 : BrineFluidProperties::h_from_p_T_X(const FPADReal & pressure,
265 : const FPADReal & temperature,
266 : const FPADReal & xnacl) const
267 : {
268 : using std::sqrt;
269 :
270 : FPADReal q1, q2, q10, q11, q12, q20, q21, q22, q23, q1x1, q2x1, Th;
271 :
272 : // The correlation requires the pressure in bar, not Pa.
273 889 : const FPADReal pbar = pressure * 1.0e-5;
274 889 : const FPADReal pbar2 = pbar * pbar;
275 :
276 : // The correlation requires mole fraction
277 889 : const FPADReal Xnacl = massFractionToMoleFraction(xnacl);
278 :
279 889 : q11 = -32.1724 + 0.0621255 * pbar;
280 889 : q21 = -1.69513 - 4.52781e-4 * pbar - 6.04279e-8 * pbar2;
281 889 : q22 = 0.0612567 + 1.88082e-5 * pbar;
282 :
283 889 : q1x1 = 47.9048 - 9.36994e-3 * pbar + 6.51059e-6 * pbar2;
284 889 : q2x1 = 0.241022 + 3.45087e-5 * pbar - 4.28356e-9 * pbar2;
285 :
286 1778 : q12 = -q11 - q1x1;
287 : q10 = q1x1;
288 :
289 889 : q20 = 1.0 - q21 * sqrt(q22);
290 889 : q23 = q2x1 - q20 - q21 * sqrt(1.0 + q22);
291 :
292 889 : q1 = q10 + q11 * (1.0 - Xnacl) + q12 * (1.0 - Xnacl) * (1.0 - Xnacl);
293 1778 : q2 = q20 + q21 * sqrt(Xnacl + q22) + q23 * Xnacl;
294 : // The temperature Th where the brine has the same enthalpy as pure water
295 : // Note: correlation uses temperature in Celcius
296 1778 : Th = q1 + q2 * (temperature - _T_c2k);
297 :
298 : // The brine enthalpy is then given by the enthalpy of water at temperature Th
299 : // Note: water enthalpy requires temperature in Kelvin
300 : FPADReal enthalpy;
301 889 : if (_water_fp_derivs)
302 : {
303 : Real h, dh_dp, dh_dT;
304 4 : _water_fp->h_from_p_T(pressure.value(), Th.value() + _T_c2k, h, dh_dp, dh_dT);
305 : enthalpy = h;
306 :
307 8 : enthalpy.derivatives() = pressure.derivatives() * dh_dp + Th.derivatives() * dh_dT;
308 : }
309 : else
310 885 : enthalpy = _water_fp->h_from_p_T(pressure.value(), Th.value() + _T_c2k);
311 :
312 889 : return enthalpy;
313 : }
314 :
315 : Real
316 885 : BrineFluidProperties::h_from_p_T_X(Real pressure, Real temperature, Real xnacl) const
317 : {
318 : // Initialise the AD value (no derivatives required)
319 : FPADReal p = pressure;
320 : FPADReal T = temperature;
321 : FPADReal x = xnacl;
322 :
323 885 : _water_fp_derivs = false;
324 885 : return h_from_p_T_X(p, T, x).value();
325 : }
326 :
327 : void
328 2 : BrineFluidProperties::h_from_p_T_X(Real pressure,
329 : Real temperature,
330 : Real xnacl,
331 : Real & h,
332 : Real & dh_dp,
333 : Real & dh_dT,
334 : Real & dh_dx) const
335 : {
336 : // Initialise the AD value and derivatives
337 : FPADReal p = pressure;
338 : Moose::derivInsert(p.derivatives(), 0, 1.0);
339 : FPADReal T = temperature;
340 : Moose::derivInsert(T.derivatives(), 1, 1.0);
341 : FPADReal x = xnacl;
342 : Moose::derivInsert(x.derivatives(), 2, 1.0);
343 :
344 2 : _water_fp_derivs = true;
345 2 : FPADReal ad_h = h_from_p_T_X(p, T, x);
346 :
347 2 : h = ad_h.value();
348 2 : dh_dp = ad_h.derivatives()[0];
349 2 : dh_dT = ad_h.derivatives()[1];
350 2 : dh_dx = ad_h.derivatives()[2];
351 2 : }
352 :
353 : Real
354 435 : BrineFluidProperties::cp_from_p_T_X(Real pressure, Real temperature, Real xnacl) const
355 : {
356 : Real q1, q2, q10, q11, q12, q20, q21, q22, q23, q1x1, q2x1, Th;
357 :
358 : // The correlation requires the pressure in bar, not Pa.
359 435 : Real pbar = pressure * 1.0e-5;
360 435 : Real pbar2 = pbar * pbar;
361 :
362 : // The correlation requires mole fraction
363 435 : Real Xnacl = massFractionToMoleFraction(xnacl);
364 :
365 435 : q11 = -32.1724 + 0.0621255 * pbar;
366 435 : q21 = -1.69513 - 4.52781e-4 * pbar - 6.04279e-8 * pbar2;
367 435 : q22 = 0.0612567 + 1.88082e-5 * pbar;
368 :
369 435 : q1x1 = 47.9048 - 9.36994e-3 * pbar + 6.51059e-6 * pbar2;
370 435 : q2x1 = 0.241022 + 3.45087e-5 * pbar - 4.28356e-9 * pbar2;
371 :
372 435 : q12 = -q11 - q1x1;
373 : q10 = q1x1;
374 :
375 435 : q20 = 1.0 - q21 * std::sqrt(q22);
376 435 : q23 = q2x1 - q20 - q21 * std::sqrt(1.0 + q22);
377 :
378 435 : q1 = q10 + q11 * (1.0 - Xnacl) + q12 * (1.0 - Xnacl) * (1.0 - Xnacl);
379 435 : q2 = q20 + q21 * std::sqrt(Xnacl + q22) + q23 * Xnacl;
380 : // The temperature Th where the brine has the same isobaric heat capacity
381 : // as pure water. Note: correlation uses temperature in Celcius
382 435 : Th = q1 + q2 * (temperature - _T_c2k);
383 :
384 : // The brine isobaric heat capacity is then given by the isobaric heat
385 : // capacity of water at temperature Th multiplied by q2
386 : // Note: water isobaric heat capacity requires temperature in Kelvin
387 435 : return q2 * _water_fp->cp_from_p_T(pressure, Th + _T_c2k);
388 : }
389 :
390 : FPADReal
391 2 : BrineFluidProperties::e_from_p_T_X(const FPADReal & pressure,
392 : const FPADReal & temperature,
393 : const FPADReal & xnacl) const
394 : {
395 2 : FPADReal enthalpy = h_from_p_T_X(pressure, temperature, xnacl);
396 2 : FPADReal density = rho_from_p_T_X(pressure, temperature, xnacl);
397 :
398 2 : return enthalpy - pressure / density;
399 : }
400 :
401 : Real
402 444 : BrineFluidProperties::e_from_p_T_X(Real pressure, Real temperature, Real xnacl) const
403 : {
404 444 : Real enthalpy = h_from_p_T_X(pressure, temperature, xnacl);
405 444 : Real density = rho_from_p_T_X(pressure, temperature, xnacl);
406 :
407 444 : return enthalpy - pressure / density;
408 : }
409 :
410 : void
411 2 : BrineFluidProperties::e_from_p_T_X(Real pressure,
412 : Real temperature,
413 : Real xnacl,
414 : Real & e,
415 : Real & de_dp,
416 : Real & de_dT,
417 : Real & de_dx) const
418 : {
419 : // Initialise the AD value and derivatives
420 : FPADReal p = pressure;
421 : Moose::derivInsert(p.derivatives(), 0, 1.0);
422 : FPADReal T = temperature;
423 : Moose::derivInsert(T.derivatives(), 1, 1.0);
424 : FPADReal x = xnacl;
425 : Moose::derivInsert(x.derivatives(), 2, 1.0);
426 :
427 2 : _water_fp_derivs = true;
428 2 : FPADReal ad_e = e_from_p_T_X(p, T, x);
429 :
430 2 : e = ad_e.value();
431 2 : de_dp = ad_e.derivatives()[0];
432 2 : de_dT = ad_e.derivatives()[1];
433 2 : de_dx = ad_e.derivatives()[2];
434 2 : }
435 :
436 : Real
437 3 : BrineFluidProperties::k_from_p_T_X(Real pressure, Real temperature, Real xnacl) const
438 : {
439 : // Correlation requires molal concentration (mol/kg)
440 3 : Real mol = massFractionToMolalConc(xnacl);
441 : // Correlation requires temperature in C
442 3 : Real Tc = temperature - _T_c2k;
443 :
444 3 : Real S = 100.0 * _Mnacl * mol / (1.0 + _Mnacl * mol);
445 3 : Real lambdaw = _water_fp->k_from_p_T(pressure, temperature);
446 3 : Real lambda = 1.0 - (2.3434e-3 - 7.924e-6 * Tc + 3.924e-8 * Tc * Tc) * S +
447 3 : (1.06e-5 - 2.0e-8 * Tc - 1.2e-10 * Tc * Tc) * S * S;
448 :
449 3 : return lambda * lambdaw;
450 : }
451 :
452 : Real
453 3 : BrineFluidProperties::vaporPressure(Real temperature, Real xnacl) const
454 : {
455 : // Correlation requires molal concentration (mol/kg)
456 3 : Real mol = massFractionToMolalConc(xnacl);
457 3 : Real mol2 = mol * mol;
458 3 : Real mol3 = mol2 * mol;
459 :
460 3 : Real a = 1.0 + 5.93582e-6 * mol - 5.19386e-5 * mol2 + 1.23156e-5 * mol3;
461 3 : Real b = 1.1542e-6 * mol + 1.41254e-7 * mol2 - 1.92476e-8 * mol3 - 1.70717e-9 * mol * mol3 +
462 3 : 1.0539e-10 * mol2 * mol3;
463 :
464 : // The temperature of pure water at the same pressure as the brine is given by
465 3 : Real th20 = std::exp(std::log(temperature) / (a + b * temperature));
466 :
467 : // The brine vapour pressure is then found by evaluating the saturation pressure for pure water
468 : // using this effective temperature
469 3 : return _water_fp->vaporPressure(th20);
470 : }
471 :
472 : Real
473 3 : BrineFluidProperties::haliteSolubility(Real temperature) const
474 : {
475 : // This correlation requires temperature in Celcius
476 3 : Real Tc = temperature - _T_c2k;
477 :
478 3 : return (26.18 + 7.2e-3 * Tc + 1.06e-4 * Tc * Tc) / 100.0;
479 : }
480 :
481 : Real
482 23 : BrineFluidProperties::massFractionToMolalConc(Real xnacl) const
483 : {
484 23 : return xnacl / ((1.0 - xnacl) * _Mnacl);
485 : }
486 :
487 : Real
488 435 : BrineFluidProperties::massFractionToMoleFraction(Real xnacl) const
489 : {
490 : // The average molar mass of brine from the mass fraction
491 435 : Real Mbrine = molarMass(xnacl);
492 : // The mole fraction is then
493 435 : return xnacl * Mbrine / _Mnacl;
494 : }
495 :
496 : FPADReal
497 1784 : BrineFluidProperties::massFractionToMoleFraction(const FPADReal & xnacl) const
498 : {
499 : // The average molar mass of brine from the mass fraction
500 1784 : FPADReal Mbrine = molarMass(xnacl);
501 : // The mole fraction is then
502 1784 : return xnacl * Mbrine / _Mnacl;
503 : }
504 :
505 : Real
506 0 : BrineFluidProperties::henryConstant(Real temperature, const std::vector<Real> & coeffs) const
507 : {
508 0 : return _water97_fp->henryConstant(temperature, coeffs);
509 : }
510 :
511 : void
512 0 : BrineFluidProperties::henryConstant(Real temperature,
513 : const std::vector<Real> & coeffs,
514 : Real & Kh,
515 : Real & dKh_dT) const
516 : {
517 0 : _water97_fp->henryConstant(temperature, coeffs, Kh, dKh_dT);
518 0 : }
519 :
520 : ADReal
521 0 : BrineFluidProperties::henryConstant(const ADReal & temperature,
522 : const std::vector<Real> & coeffs) const
523 : {
524 : Real Kh, dKh_dT;
525 0 : henryConstant(temperature.value(), coeffs, Kh, dKh_dT);
526 :
527 0 : ADReal henry = Kh;
528 0 : henry.derivatives() = temperature.derivatives() * dKh_dT;
529 :
530 0 : return henry;
531 : }
|