blackoilbrinemodules.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
28#ifndef EWOMS_BLACK_OIL_BRINE_MODULE_HH
29#define EWOMS_BLACK_OIL_BRINE_MODULE_HH
30
31#include <dune/common/fvector.hh>
32
35
38
40
41#include <cassert>
42#include <istream>
43#include <ostream>
44#include <stdexcept>
45#include <string>
46
47namespace Opm {
48
54template <class TypeTag, bool enableBrineV = getPropValue<TypeTag, Properties::EnableBrine>()>
56{
69
70 using Toolbox = MathToolbox<Evaluation>;
71
72 using TabulatedFunction = typename BlackOilBrineParams<Scalar>::TabulatedFunction;
73
74 static constexpr unsigned saltConcentrationIdx = Indices::saltConcentrationIdx;
75 static constexpr unsigned contiBrineEqIdx = Indices::contiBrineEqIdx;
76 static constexpr unsigned waterPhaseIdx = FluidSystem::waterPhaseIdx;
77 static constexpr bool gasEnabled = Indices::gasEnabled;
78 static constexpr bool oilEnabled = Indices::oilEnabled;
79 static constexpr bool enableBrine = enableBrineV;
80 static constexpr bool enableSaltPrecipitation =
81 getPropValue<TypeTag, Properties::EnableSaltPrecipitation>();
82
83 static constexpr unsigned numEq = getPropValue<TypeTag, Properties::NumEq>();
84 static constexpr unsigned numPhases = FluidSystem::numPhases;
85
86public:
89 {
90 params_ = params;
91 }
92
96 static void registerParameters()
97 {}
98
99 static bool primaryVarApplies(unsigned pvIdx)
100 {
101 if constexpr (enableBrine) {
102 return pvIdx == saltConcentrationIdx;
103 }
104 else {
105 return false;
106 }
107 }
108
112 template <class FluidState>
113 static void assignPrimaryVars(PrimaryVariables& priVars,
114 const FluidState& fluidState)
115 {
116 if constexpr (enableBrine) {
117 priVars[saltConcentrationIdx] = fluidState.saltConcentration();
118 }
119 }
120
121 static std::string primaryVarName([[maybe_unused]] unsigned pvIdx)
122 {
123 assert(primaryVarApplies(pvIdx));
124
125 return "saltConcentration";
126 }
127
128 static Scalar primaryVarWeight([[maybe_unused]] unsigned pvIdx)
129 {
130 assert(primaryVarApplies(pvIdx));
131
132 // TODO: it may be beneficial to choose this differently.
133 return static_cast<Scalar>(1.0);
134 }
135
136 static bool eqApplies(unsigned eqIdx)
137 {
138 if constexpr (enableBrine) {
139 return eqIdx == contiBrineEqIdx;
140 }
141 else {
142 return false;
143 }
144 }
145
146 static std::string eqName([[maybe_unused]] unsigned eqIdx)
147 {
148 assert(eqApplies(eqIdx));
149
150 return "conti^brine";
151 }
152
153 static Scalar eqWeight([[maybe_unused]] unsigned eqIdx)
154 {
155 assert(eqApplies(eqIdx));
156
157 // TODO: it may be beneficial to choose this differently.
158 return static_cast<Scalar>(1.0);
159 }
160
161 // must be called after water storage is computed
162 template <class LhsEval>
163 static void addStorage(Dune::FieldVector<LhsEval, numEq>& storage,
164 const IntensiveQuantities& intQuants)
165 {
166 if constexpr (enableBrine) {
167 const auto& fs = intQuants.fluidState();
168
169 // avoid singular matrix if no water is present.
170 const LhsEval surfaceVolumeWater =
171 max(Toolbox::template decay<LhsEval>(fs.saturation(waterPhaseIdx)) *
172 Toolbox::template decay<LhsEval>(fs.invB(waterPhaseIdx)) *
173 Toolbox::template decay<LhsEval>(intQuants.porosity()),
174 1e-10);
175
176 // Brine in water phase
177 const LhsEval massBrine = surfaceVolumeWater *
178 Toolbox::template decay<LhsEval>(fs.saltConcentration());
179
180 if constexpr (enableSaltPrecipitation) {
181 const double saltDensity = intQuants.saltDensity(); // Solid salt density kg/m3
182 const LhsEval solidSalt =
183 Toolbox::template decay<LhsEval>(intQuants.porosity()) /
184 (1.0 - Toolbox::template decay<LhsEval>(fs.saltSaturation()) + 1.e-8) *
185 saltDensity *
186 Toolbox::template decay<LhsEval>(fs.saltSaturation());
187
188 storage[contiBrineEqIdx] += massBrine + solidSalt;
189 }
190 else {
191 storage[contiBrineEqIdx] += massBrine;
192 }
193 }
194 }
195
196 static void computeFlux([[maybe_unused]] RateVector& flux,
197 [[maybe_unused]] const ElementContext& elemCtx,
198 [[maybe_unused]] unsigned scvfIdx,
199 [[maybe_unused]] unsigned timeIdx)
200 {
201 if constexpr (enableBrine) {
202 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
203 unsigned focusIdx = elemCtx.focusDofIndex();
204 unsigned upIdx = extQuants.upstreamIndex(waterPhaseIdx);
205 flux[contiBrineEqIdx] = 0.0;
206 if (upIdx == focusIdx)
207 addBrineFluxes_<Evaluation>(flux, elemCtx, scvfIdx, timeIdx);
208 else
209 addBrineFluxes_<Scalar>(flux, elemCtx, scvfIdx, timeIdx);
210 }
211 }
212
213 template <class UpstreamEval>
214 static void addBrineFluxes_(RateVector& flux,
215 const ElementContext& elemCtx,
216 unsigned scvfIdx,
217 unsigned timeIdx)
218 {
219 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
220 unsigned upIdx = extQuants.upstreamIndex(waterPhaseIdx);
221 const auto& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
222 const auto& upFs = up.fluidState();
223 const auto& volFlux = extQuants.volumeFlux(waterPhaseIdx);
224 addBrineFluxes_<UpstreamEval>(flux, waterPhaseIdx, volFlux, upFs);
225 }
226
227 template <class UpEval, class FluidState>
228 static void addBrineFluxes_(RateVector& flux,
229 unsigned phaseIdx,
230 const Evaluation& volFlux,
231 const FluidState& upFs)
232 {
233 if constexpr (enableBrine) {
234 if (phaseIdx == waterPhaseIdx) {
235 flux[contiBrineEqIdx] =
236 decay<UpEval>(upFs.saltConcentration())
237 * decay<UpEval>(upFs.invB(waterPhaseIdx))
238 * volFlux;
239 }
240 }
241 }
242
246 static Scalar computeUpdateError(const PrimaryVariables&,
247 const EqVector&)
248 {
249 // do not consider the change of Brine primary variables for
250 // convergence
251 // TODO: maybe this should be changed
252 return static_cast<Scalar>(0.0);
253 }
254
255 template <class DofEntity>
256 static void serializeEntity(const Model& model, std::ostream& outstream, const DofEntity& dof)
257 {
258 if constexpr (enableBrine) {
259 const unsigned dofIdx = model.dofMapper().index(dof);
260 const PrimaryVariables& priVars = model.solution(/*timeIdx=*/0)[dofIdx];
261 outstream << priVars[saltConcentrationIdx];
262 }
263 }
264
265 template <class DofEntity>
266 static void deserializeEntity(Model& model, std::istream& instream, const DofEntity& dof)
267 {
268 if constexpr (enableBrine) {
269 const unsigned dofIdx = model.dofMapper().index(dof);
270 PrimaryVariables& priVars0 = model.solution(/*timeIdx=*/0)[dofIdx];
271 PrimaryVariables& priVars1 = model.solution(/*timeIdx=*/1)[dofIdx];
272
273 instream >> priVars0[saltConcentrationIdx];
274
275 // set the primary variables for the beginning of the current time step.
276 priVars1[saltConcentrationIdx] = priVars0[saltConcentrationIdx];
277 }
278 }
279
280 static Scalar referencePressure(const ElementContext& elemCtx,
281 unsigned scvIdx,
282 unsigned timeIdx)
283 {
284 const unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
285 return params_.referencePressure_[pvtnumRegionIdx];
286 }
287
288 static const TabulatedFunction& bdensityTable(const ElementContext& elemCtx,
289 unsigned scvIdx,
290 unsigned timeIdx)
291 {
292 const unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
293 return params_.bdensityTable_[pvtnumRegionIdx];
294 }
295
296 static const TabulatedFunction& pcfactTable(unsigned satnumRegionIdx)
297 { return params_.pcfactTable_[satnumRegionIdx]; }
298
299 static const TabulatedFunction& permfactTable(const ElementContext& elemCtx,
300 unsigned scvIdx,
301 unsigned timeIdx)
302 {
303 const unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
304 return params_.permfactTable_[satnumRegionIdx];
305 }
306
307 static const TabulatedFunction& permfactTable(unsigned satnumRegionIdx)
308 { return params_.permfactTable_[satnumRegionIdx]; }
309
310 static Scalar saltsolTable(const ElementContext& elemCtx,
311 unsigned scvIdx,
312 unsigned timeIdx)
313 {
314 const unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
315 return params_.saltsolTable_[pvtnumRegionIdx];
316 }
317
318 static Scalar saltsolTable(const unsigned pvtnumRegionIdx)
319 {
320 return params_.saltsolTable_[pvtnumRegionIdx];
321 }
322
323 static Scalar saltdenTable(const ElementContext& elemCtx,
324 unsigned scvIdx,
325 unsigned timeIdx)
326 {
327 const unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
328 return params_.saltdenTable_[pvtnumRegionIdx];
329 }
330
331 static Scalar saltdenTable(const unsigned pvtnumRegionIdx)
332 {
333 return params_.saltdenTable_[pvtnumRegionIdx];
334 }
335
336 static bool hasBDensityTables()
337 { return !params_.bdensityTable_.empty(); }
338
339 static bool hasSaltsolTables()
340 { return !params_.saltsolTable_.empty(); }
341
342 static bool hasPcfactTables()
343 {
344 if constexpr (enableSaltPrecipitation) {
345 return !params_.pcfactTable_.empty();
346 }
347 else {
348 return false;
349 }
350 }
351
352 static Scalar saltSol(unsigned regionIdx)
353 { return params_.saltsolTable_[regionIdx]; }
354
355private:
356 static BlackOilBrineParams<Scalar> params_;
357};
358
359template <class TypeTag, bool enableBrineV>
360BlackOilBrineParams<typename BlackOilBrineModule<TypeTag, enableBrineV>::Scalar>
361BlackOilBrineModule<TypeTag, enableBrineV>::params_;
362
363template <class TypeTag, bool enableBrineV>
365
373template <class TypeTag>
374class BlackOilBrineIntensiveQuantities<TypeTag, /*enableBrineV=*/true>
375{
377
385
387
388 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
389 static constexpr int saltConcentrationIdx = Indices::saltConcentrationIdx;
390 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
391 static constexpr int gasPhaseIdx = FluidSystem::gasPhaseIdx;
392 static constexpr int oilPhaseIdx = FluidSystem::oilPhaseIdx;
393 static constexpr bool enableBrine = true;
394 static constexpr bool enableSaltPrecipitation =
395 getPropValue<TypeTag, Properties::EnableSaltPrecipitation>();
396 static constexpr int contiBrineEqIdx = Indices::contiBrineEqIdx;
397
398public:
404 void updateSaltConcentration_(const ElementContext& elemCtx,
405 unsigned dofIdx,
406 unsigned timeIdx)
407 {
408 const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
409 const LinearizationType lintype = elemCtx.linearizationType();
410 updateSaltConcentration_(priVars, timeIdx, lintype);
411 }
412
413 void updateSaltConcentration_(const PrimaryVariables& priVars,
414 const unsigned timeIdx,
415 const LinearizationType lintype)
416 {
417 const unsigned pvtnumRegionIdx = priVars.pvtRegionIndex();
418 auto& fs = asImp_().fluidState_;
419
420 if constexpr (enableSaltPrecipitation) {
421 saltSolubility_ = BrineModule::saltsolTable(pvtnumRegionIdx);
422 saltDensity_ = BrineModule::saltdenTable(pvtnumRegionIdx);
423
424 if (priVars.primaryVarsMeaningBrine() == PrimaryVariables::BrineMeaning::Sp) {
425 saltSaturation_ = priVars.makeEvaluation(saltConcentrationIdx, timeIdx, lintype);
426 fs.setSaltConcentration(saltSolubility_);
427 }
428 else {
429 saltConcentration_ = priVars.makeEvaluation(saltConcentrationIdx, timeIdx, lintype);
430 fs.setSaltConcentration(saltConcentration_);
431 saltSaturation_ = 0.0;
432 }
433 fs.setSaltSaturation(saltSaturation_);
434 }
435 else {
436 saltConcentration_ = priVars.makeEvaluation(saltConcentrationIdx, timeIdx, lintype);
437 fs.setSaltConcentration(saltConcentration_);
438 }
439 }
440
441 void saltPropertiesUpdate_([[maybe_unused]] const ElementContext& elemCtx,
442 [[maybe_unused]] unsigned dofIdx,
443 [[maybe_unused]] unsigned timeIdx)
444 {
445 if constexpr (enableSaltPrecipitation) {
446 const Evaluation porosityFactor = min(1.0 - asImp_().fluidState_.saltSaturation(), 1.0); //phi/phi_0
447
448 const auto& permfactTable = BrineModule::permfactTable(elemCtx, dofIdx, timeIdx);
449
450 permFactor_ = permfactTable.eval(porosityFactor);
451 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
452 if (!FluidSystem::phaseIsActive(phaseIdx)) {
453 continue;
454 }
455
456 asImp_().mobility_[phaseIdx] *= permFactor_;
457 }
458 }
459 }
460
461 const Evaluation& brineRefDensity() const
462 { return refDensity_; }
463
464 Scalar saltSolubility() const
465 { return saltSolubility_; }
466
467 Scalar saltDensity() const
468 { return saltDensity_; }
469
470 const Evaluation& permFactor() const
471 { return permFactor_; }
472
473protected:
474 Implementation& asImp_()
475 { return *static_cast<Implementation*>(this); }
476
478 Evaluation refDensity_;
479 Evaluation saltSaturation_;
480 Evaluation permFactor_;
483};
484
485template <class TypeTag>
487{
491
492public:
493 void updateSaltConcentration_(const ElementContext&,
494 unsigned,
495 unsigned)
496 {}
497
498 void saltPropertiesUpdate_(const ElementContext&,
499 unsigned,
500 unsigned)
501 {}
502
503 const Evaluation& brineRefDensity() const
504 { throw std::runtime_error("brineRefDensity() called but brine are disabled"); }
505
506 const Scalar saltSolubility() const
507 { throw std::logic_error("saltSolubility() called but salt precipitation is disabled"); }
508
509 const Scalar saltDensity() const
510 { throw std::logic_error("saltDensity() called but salt precipitation is disabled"); }
511
512 const Evaluation& permFactor() const
513 { throw std::logic_error("permFactor() called but salt precipitation is disabled"); }
514};
515
516} // namespace Opm
517
518#endif
Defines a type tags and some fundamental properties all models.
Contains the parameters required to extend the black-oil model by brine.
Declares the properties required by the black oil model.
void updateSaltConcentration_(const ElementContext &, unsigned, unsigned)
Definition: blackoilbrinemodules.hh:493
const Scalar saltDensity() const
Definition: blackoilbrinemodules.hh:509
const Evaluation & permFactor() const
Definition: blackoilbrinemodules.hh:512
void saltPropertiesUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilbrinemodules.hh:498
const Scalar saltSolubility() const
Definition: blackoilbrinemodules.hh:506
const Evaluation & brineRefDensity() const
Definition: blackoilbrinemodules.hh:503
void updateSaltConcentration_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Update the intensive properties needed to handle brine from the primary variables.
Definition: blackoilbrinemodules.hh:404
Evaluation refDensity_
Definition: blackoilbrinemodules.hh:478
const Evaluation & brineRefDensity() const
Definition: blackoilbrinemodules.hh:461
Evaluation permFactor_
Definition: blackoilbrinemodules.hh:480
Evaluation saltConcentration_
Definition: blackoilbrinemodules.hh:477
Scalar saltSolubility_
Definition: blackoilbrinemodules.hh:481
const Evaluation & permFactor() const
Definition: blackoilbrinemodules.hh:470
Scalar saltDensity_
Definition: blackoilbrinemodules.hh:482
void updateSaltConcentration_(const PrimaryVariables &priVars, const unsigned timeIdx, const LinearizationType lintype)
Definition: blackoilbrinemodules.hh:413
Evaluation saltSaturation_
Definition: blackoilbrinemodules.hh:479
Scalar saltSolubility() const
Definition: blackoilbrinemodules.hh:464
Implementation & asImp_()
Definition: blackoilbrinemodules.hh:474
void saltPropertiesUpdate_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:441
Scalar saltDensity() const
Definition: blackoilbrinemodules.hh:467
Definition: blackoilbrinemodules.hh:364
Contains the high level supplements required to extend the black oil model by brine.
Definition: blackoilbrinemodules.hh:56
static bool hasBDensityTables()
Definition: blackoilbrinemodules.hh:336
static const TabulatedFunction & bdensityTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:288
static void serializeEntity(const Model &model, std::ostream &outstream, const DofEntity &dof)
Definition: blackoilbrinemodules.hh:256
static Scalar saltSol(unsigned regionIdx)
Definition: blackoilbrinemodules.hh:352
static std::string eqName(unsigned eqIdx)
Definition: blackoilbrinemodules.hh:146
static bool primaryVarApplies(unsigned pvIdx)
Definition: blackoilbrinemodules.hh:99
static const TabulatedFunction & pcfactTable(unsigned satnumRegionIdx)
Definition: blackoilbrinemodules.hh:296
static const TabulatedFunction & permfactTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:299
static bool hasPcfactTables()
Definition: blackoilbrinemodules.hh:342
static Scalar referencePressure(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:280
static Scalar saltsolTable(const unsigned pvtnumRegionIdx)
Definition: blackoilbrinemodules.hh:318
static std::string primaryVarName(unsigned pvIdx)
Definition: blackoilbrinemodules.hh:121
static void setParams(BlackOilBrineParams< Scalar > &&params)
Set parameters.
Definition: blackoilbrinemodules.hh:88
static void registerParameters()
Register all run-time parameters for the black-oil brine module.
Definition: blackoilbrinemodules.hh:96
static const TabulatedFunction & permfactTable(unsigned satnumRegionIdx)
Definition: blackoilbrinemodules.hh:307
static bool hasSaltsolTables()
Definition: blackoilbrinemodules.hh:339
static Scalar eqWeight(unsigned eqIdx)
Definition: blackoilbrinemodules.hh:153
static Scalar saltdenTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:323
static Scalar saltsolTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:310
static Scalar primaryVarWeight(unsigned pvIdx)
Definition: blackoilbrinemodules.hh:128
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:196
static bool eqApplies(unsigned eqIdx)
Definition: blackoilbrinemodules.hh:136
static Scalar saltdenTable(const unsigned pvtnumRegionIdx)
Definition: blackoilbrinemodules.hh:331
static void assignPrimaryVars(PrimaryVariables &priVars, const FluidState &fluidState)
Assign the brine specific primary variables to a PrimaryVariables object.
Definition: blackoilbrinemodules.hh:113
static void deserializeEntity(Model &model, std::istream &instream, const DofEntity &dof)
Definition: blackoilbrinemodules.hh:266
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilbrinemodules.hh:163
static void addBrineFluxes_(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:214
static void addBrineFluxes_(RateVector &flux, unsigned phaseIdx, const Evaluation &volFlux, const FluidState &upFs)
Definition: blackoilbrinemodules.hh:228
static Scalar computeUpdateError(const PrimaryVariables &, const EqVector &)
Return how much a Newton-Raphson update is considered an error.
Definition: blackoilbrinemodules.hh:246
Declare the properties used by the infrastructure code of the finite volume discretizations.
Definition: blackoilbioeffectsmodules.hh:43
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:233
Struct holding the parameters for the BlackoilBrineModule class.
Definition: blackoilbrineparams.hpp:44
Tabulated1DFunction< Scalar > TabulatedFunction
Definition: blackoilbrineparams.hpp:50
Definition: linearizationtype.hh:34