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
37
39
40#include <cassert>
41#include <istream>
42#include <ostream>
43#include <stdexcept>
44#include <string>
45
46namespace Opm {
47
53template <class TypeTag, bool enableBrineV = getPropValue<TypeTag, Properties::EnableBrine>()>
55{
68
69 using Toolbox = MathToolbox<Evaluation>;
70
71 using TabulatedFunction = typename BlackOilBrineParams<Scalar>::TabulatedFunction;
72
73 static constexpr unsigned saltConcentrationIdx = Indices::saltConcentrationIdx;
74 static constexpr unsigned contiBrineEqIdx = Indices::contiBrineEqIdx;
75 static constexpr unsigned waterPhaseIdx = FluidSystem::waterPhaseIdx;
76 static constexpr bool gasEnabled = Indices::gasEnabled;
77 static constexpr bool oilEnabled = Indices::oilEnabled;
78 static constexpr bool enableBrine = enableBrineV;
79 static constexpr bool enableSaltPrecipitation =
80 getPropValue<TypeTag, Properties::EnableSaltPrecipitation>();
81
82 static constexpr unsigned numEq = getPropValue<TypeTag, Properties::NumEq>();
83 static constexpr unsigned numPhases = FluidSystem::numPhases;
84
85public:
88 {
89 params_ = params;
90 }
91
95 static void registerParameters()
96 {}
97
98 static bool primaryVarApplies(unsigned pvIdx)
99 {
100 if constexpr (enableBrine) {
101 return pvIdx == saltConcentrationIdx;
102 }
103 else {
104 return false;
105 }
106 }
107
111 template <class FluidState>
112 static void assignPrimaryVars(PrimaryVariables& priVars,
113 const FluidState& fluidState)
114 {
115 if constexpr (enableBrine) {
116 priVars[saltConcentrationIdx] = fluidState.saltConcentration();
117 }
118 }
119
120 static std::string primaryVarName([[maybe_unused]] unsigned pvIdx)
121 {
122 assert(primaryVarApplies(pvIdx));
123
124 return "saltConcentration";
125 }
126
127 static Scalar primaryVarWeight([[maybe_unused]] unsigned pvIdx)
128 {
129 assert(primaryVarApplies(pvIdx));
130
131 // TODO: it may be beneficial to choose this differently.
132 return static_cast<Scalar>(1.0);
133 }
134
135 static bool eqApplies(unsigned eqIdx)
136 {
137 if constexpr (enableBrine) {
138 return eqIdx == contiBrineEqIdx;
139 }
140 else {
141 return false;
142 }
143 }
144
145 static std::string eqName([[maybe_unused]] unsigned eqIdx)
146 {
147 assert(eqApplies(eqIdx));
148
149 return "conti^brine";
150 }
151
152 static Scalar eqWeight([[maybe_unused]] unsigned eqIdx)
153 {
154 assert(eqApplies(eqIdx));
155
156 // TODO: it may be beneficial to choose this differently.
157 return static_cast<Scalar>(1.0);
158 }
159
160 // must be called after water storage is computed
161 template <class LhsEval>
162 static void addStorage(Dune::FieldVector<LhsEval, numEq>& storage,
163 const IntensiveQuantities& intQuants)
164 {
165 if constexpr (enableBrine) {
166 const auto& fs = intQuants.fluidState();
167
168 // avoid singular matrix if no water is present.
169 const LhsEval surfaceVolumeWater =
170 max(Toolbox::template decay<LhsEval>(fs.saturation(waterPhaseIdx)) *
171 Toolbox::template decay<LhsEval>(fs.invB(waterPhaseIdx)) *
172 Toolbox::template decay<LhsEval>(intQuants.porosity()),
173 1e-10);
174
175 // Brine in water phase
176 const LhsEval massBrine = surfaceVolumeWater *
177 Toolbox::template decay<LhsEval>(fs.saltConcentration());
178
179 if constexpr (enableSaltPrecipitation) {
180 const double saltDensity = intQuants.saltDensity(); // Solid salt density kg/m3
181 const LhsEval solidSalt =
182 Toolbox::template decay<LhsEval>(intQuants.porosity()) /
183 (1.0 - Toolbox::template decay<LhsEval>(intQuants.saltSaturation()) + 1.e-8) *
184 saltDensity *
185 Toolbox::template decay<LhsEval>(intQuants.saltSaturation());
186
187 storage[contiBrineEqIdx] += massBrine + solidSalt;
188 }
189 else {
190 storage[contiBrineEqIdx] += massBrine;
191 }
192 }
193 }
194
195 static void computeFlux([[maybe_unused]] RateVector& flux,
196 [[maybe_unused]] const ElementContext& elemCtx,
197 [[maybe_unused]] unsigned scvfIdx,
198 [[maybe_unused]] unsigned timeIdx)
199
200 {
201 if constexpr (enableBrine) {
202 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
203
204 const unsigned upIdx = extQuants.upstreamIndex(FluidSystem::waterPhaseIdx);
205 const unsigned inIdx = extQuants.interiorIndex();
206 const auto& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
207
208 if (upIdx == inIdx) {
209 flux[contiBrineEqIdx] =
210 extQuants.volumeFlux(waterPhaseIdx) *
211 up.fluidState().invB(waterPhaseIdx) *
212 up.fluidState().saltConcentration();
213 }
214 else {
215 flux[contiBrineEqIdx] =
216 extQuants.volumeFlux(waterPhaseIdx) *
217 decay<Scalar>(up.fluidState().invB(waterPhaseIdx)) *
218 decay<Scalar>(up.fluidState().saltConcentration());
219 }
220 }
221 }
222
226 static Scalar computeUpdateError(const PrimaryVariables&,
227 const EqVector&)
228 {
229 // do not consider the change of Brine primary variables for
230 // convergence
231 // TODO: maybe this should be changed
232 return static_cast<Scalar>(0.0);
233 }
234
235 template <class DofEntity>
236 static void serializeEntity(const Model& model, std::ostream& outstream, const DofEntity& dof)
237 {
238 if constexpr (enableBrine) {
239 const unsigned dofIdx = model.dofMapper().index(dof);
240 const PrimaryVariables& priVars = model.solution(/*timeIdx=*/0)[dofIdx];
241 outstream << priVars[saltConcentrationIdx];
242 }
243 }
244
245 template <class DofEntity>
246 static void deserializeEntity(Model& model, std::istream& instream, const DofEntity& dof)
247 {
248 if constexpr (enableBrine) {
249 const unsigned dofIdx = model.dofMapper().index(dof);
250 PrimaryVariables& priVars0 = model.solution(/*timeIdx=*/0)[dofIdx];
251 PrimaryVariables& priVars1 = model.solution(/*timeIdx=*/1)[dofIdx];
252
253 instream >> priVars0[saltConcentrationIdx];
254
255 // set the primary variables for the beginning of the current time step.
256 priVars1[saltConcentrationIdx] = priVars0[saltConcentrationIdx];
257 }
258 }
259
260 static const Scalar& referencePressure(const ElementContext& elemCtx,
261 unsigned scvIdx,
262 unsigned timeIdx)
263 {
264 const unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
265 return params_.referencePressure_[pvtnumRegionIdx];
266 }
267
268 static const TabulatedFunction& bdensityTable(const ElementContext& elemCtx,
269 unsigned scvIdx,
270 unsigned timeIdx)
271 {
272 const unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
273 return params_.bdensityTable_[pvtnumRegionIdx];
274 }
275
276 static const TabulatedFunction& pcfactTable(unsigned satnumRegionIdx)
277 { return params_.pcfactTable_[satnumRegionIdx]; }
278
279 static const TabulatedFunction& permfactTable(const ElementContext& elemCtx,
280 unsigned scvIdx,
281 unsigned timeIdx)
282 {
283 const unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
284 return params_.permfactTable_[satnumRegionIdx];
285 }
286
287 static const TabulatedFunction& permfactTable(unsigned satnumRegionIdx)
288 { return params_.permfactTable_[satnumRegionIdx]; }
289
290 static const Scalar saltsolTable(const ElementContext& elemCtx,
291 unsigned scvIdx,
292 unsigned timeIdx)
293 {
294 const unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
295 return params_.saltsolTable_[pvtnumRegionIdx];
296 }
297
298 static const Scalar saltdenTable(const ElementContext& elemCtx,
299 unsigned scvIdx,
300 unsigned timeIdx)
301 {
302 const unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
303 return params_.saltdenTable_[pvtnumRegionIdx];
304 }
305
306 static bool hasBDensityTables()
307 { return !params_.bdensityTable_.empty(); }
308
309 static bool hasSaltsolTables()
310 { return !params_.saltsolTable_.empty(); }
311
312 static bool hasPcfactTables()
313 {
314 if constexpr (enableSaltPrecipitation) {
315 return !params_.pcfactTable_.empty();
316 }
317 else {
318 return false;
319 }
320 }
321
322 static Scalar saltSol(unsigned regionIdx)
323 { return params_.saltsolTable_[regionIdx]; }
324
325private:
326 static BlackOilBrineParams<Scalar> params_;
327};
328
329template <class TypeTag, bool enableBrineV>
330BlackOilBrineParams<typename BlackOilBrineModule<TypeTag, enableBrineV>::Scalar>
331BlackOilBrineModule<TypeTag, enableBrineV>::params_;
332
333template <class TypeTag, bool enableBrineV>
335
343template <class TypeTag>
344class BlackOilBrineIntensiveQuantities<TypeTag, /*enableBrineV=*/true>
345{
347
355
357
358 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
359 static constexpr int saltConcentrationIdx = Indices::saltConcentrationIdx;
360 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
361 static constexpr int gasPhaseIdx = FluidSystem::gasPhaseIdx;
362 static constexpr int oilPhaseIdx = FluidSystem::oilPhaseIdx;
363 static constexpr bool enableBrine = true;
364 static constexpr bool enableSaltPrecipitation =
365 getPropValue<TypeTag, Properties::EnableSaltPrecipitation>();
366 static constexpr int contiBrineEqIdx = Indices::contiBrineEqIdx;
367
368public:
374 void updateSaltConcentration_(const ElementContext& elemCtx,
375 unsigned dofIdx,
376 unsigned timeIdx)
377 {
378 const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
379
380 auto& fs = asImp_().fluidState_;
381
382 if constexpr (enableSaltPrecipitation) {
383 const auto& saltsolTable = BrineModule::saltsolTable(elemCtx, dofIdx, timeIdx);
384 saltSolubility_ = saltsolTable;
385
386 const auto& saltdenTable = BrineModule::saltdenTable(elemCtx, dofIdx, timeIdx);
387 saltDensity_ = saltdenTable;
388
389 if (priVars.primaryVarsMeaningBrine() == PrimaryVariables::BrineMeaning::Sp) {
390 saltSaturation_ = priVars.makeEvaluation(saltConcentrationIdx, timeIdx);
391 fs.setSaltConcentration(saltSolubility_);
392 }
393 else {
394 saltConcentration_ = priVars.makeEvaluation(saltConcentrationIdx, timeIdx);
395 fs.setSaltConcentration(saltConcentration_);
396 saltSaturation_ = 0.0;
397 }
398 fs.setSaltSaturation(saltSaturation_);
399 }
400 else {
401 saltConcentration_ = priVars.makeEvaluation(saltConcentrationIdx, timeIdx);
402 fs.setSaltConcentration(saltConcentration_);
403 }
404 }
405
406 void saltPropertiesUpdate_([[maybe_unused]] const ElementContext& elemCtx,
407 [[maybe_unused]] unsigned dofIdx,
408 [[maybe_unused]] unsigned timeIdx)
409 {
410 if constexpr (enableSaltPrecipitation) {
411 const Evaluation porosityFactor = min(1.0 - saltSaturation(), 1.0); //phi/phi_0
412
413 const auto& permfactTable = BrineModule::permfactTable(elemCtx, dofIdx, timeIdx);
414
415 permFactor_ = permfactTable.eval(porosityFactor);
416 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
417 if (!FluidSystem::phaseIsActive(phaseIdx)) {
418 continue;
419 }
420
421 asImp_().mobility_[phaseIdx] *= permFactor_;
422 }
423 }
424 }
425
426 const Evaluation& saltConcentration() const
427 { return saltConcentration_; }
428
429 const Evaluation& brineRefDensity() const
430 { return refDensity_; }
431
432 const Evaluation& saltSaturation() const
433 { return saltSaturation_; }
434
435 Scalar saltSolubility() const
436 { return saltSolubility_; }
437
438 Scalar saltDensity() const
439 { return saltDensity_; }
440
441 const Evaluation& permFactor() const
442 { return permFactor_; }
443
444protected:
445 Implementation& asImp_()
446 { return *static_cast<Implementation*>(this); }
447
449 Evaluation refDensity_;
450 Evaluation saltSaturation_;
451 Evaluation permFactor_;
454};
455
456template <class TypeTag>
458{
462
463public:
464 void updateSaltConcentration_(const ElementContext&,
465 unsigned,
466 unsigned)
467 {}
468
469 void saltPropertiesUpdate_(const ElementContext&,
470 unsigned,
471 unsigned)
472 {}
473
474 const Evaluation& saltConcentration() const
475 { throw std::runtime_error("saltConcentration() called but brine are disabled"); }
476
477 const Evaluation& brineRefDensity() const
478 { throw std::runtime_error("brineRefDensity() called but brine are disabled"); }
479
480 const Evaluation& saltSaturation() const
481 { throw std::logic_error("saltSaturation() called but salt precipitation is disabled"); }
482
483 const Scalar saltSolubility() const
484 { throw std::logic_error("saltSolubility() called but salt precipitation is disabled"); }
485
486 const Scalar saltDensity() const
487 { throw std::logic_error("saltDensity() called but salt precipitation is disabled"); }
488
489 const Evaluation& permFactor() const
490 { throw std::logic_error("permFactor() called but salt precipitation is disabled"); }
491};
492
493} // namespace Opm
494
495#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:464
const Scalar saltDensity() const
Definition: blackoilbrinemodules.hh:486
const Evaluation & permFactor() const
Definition: blackoilbrinemodules.hh:489
const Evaluation & saltSaturation() const
Definition: blackoilbrinemodules.hh:480
void saltPropertiesUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilbrinemodules.hh:469
const Scalar saltSolubility() const
Definition: blackoilbrinemodules.hh:483
const Evaluation & brineRefDensity() const
Definition: blackoilbrinemodules.hh:477
const Evaluation & saltConcentration() const
Definition: blackoilbrinemodules.hh:474
void updateSaltConcentration_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Update the intensive properties needed to handle brine from the primary variables.
Definition: blackoilbrinemodules.hh:374
Evaluation refDensity_
Definition: blackoilbrinemodules.hh:449
const Evaluation & saltSaturation() const
Definition: blackoilbrinemodules.hh:432
const Evaluation & brineRefDensity() const
Definition: blackoilbrinemodules.hh:429
Evaluation permFactor_
Definition: blackoilbrinemodules.hh:451
Evaluation saltConcentration_
Definition: blackoilbrinemodules.hh:448
Scalar saltSolubility_
Definition: blackoilbrinemodules.hh:452
const Evaluation & permFactor() const
Definition: blackoilbrinemodules.hh:441
const Evaluation & saltConcentration() const
Definition: blackoilbrinemodules.hh:426
Scalar saltDensity_
Definition: blackoilbrinemodules.hh:453
Evaluation saltSaturation_
Definition: blackoilbrinemodules.hh:450
Scalar saltSolubility() const
Definition: blackoilbrinemodules.hh:435
Implementation & asImp_()
Definition: blackoilbrinemodules.hh:445
void saltPropertiesUpdate_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:406
Scalar saltDensity() const
Definition: blackoilbrinemodules.hh:438
Definition: blackoilbrinemodules.hh:334
Contains the high level supplements required to extend the black oil model by brine.
Definition: blackoilbrinemodules.hh:55
static bool hasBDensityTables()
Definition: blackoilbrinemodules.hh:306
static const TabulatedFunction & bdensityTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:268
static void serializeEntity(const Model &model, std::ostream &outstream, const DofEntity &dof)
Definition: blackoilbrinemodules.hh:236
static Scalar saltSol(unsigned regionIdx)
Definition: blackoilbrinemodules.hh:322
static std::string eqName(unsigned eqIdx)
Definition: blackoilbrinemodules.hh:145
static bool primaryVarApplies(unsigned pvIdx)
Definition: blackoilbrinemodules.hh:98
static const TabulatedFunction & pcfactTable(unsigned satnumRegionIdx)
Definition: blackoilbrinemodules.hh:276
static const Scalar saltdenTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:298
static const TabulatedFunction & permfactTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:279
static bool hasPcfactTables()
Definition: blackoilbrinemodules.hh:312
static std::string primaryVarName(unsigned pvIdx)
Definition: blackoilbrinemodules.hh:120
static void setParams(BlackOilBrineParams< Scalar > &&params)
Set parameters.
Definition: blackoilbrinemodules.hh:87
static void registerParameters()
Register all run-time parameters for the black-oil brine module.
Definition: blackoilbrinemodules.hh:95
static const TabulatedFunction & permfactTable(unsigned satnumRegionIdx)
Definition: blackoilbrinemodules.hh:287
static bool hasSaltsolTables()
Definition: blackoilbrinemodules.hh:309
static Scalar eqWeight(unsigned eqIdx)
Definition: blackoilbrinemodules.hh:152
static const Scalar saltsolTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:290
static Scalar primaryVarWeight(unsigned pvIdx)
Definition: blackoilbrinemodules.hh:127
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:195
static bool eqApplies(unsigned eqIdx)
Definition: blackoilbrinemodules.hh:135
static void assignPrimaryVars(PrimaryVariables &priVars, const FluidState &fluidState)
Assign the brine specific primary variables to a PrimaryVariables object.
Definition: blackoilbrinemodules.hh:112
static void deserializeEntity(Model &model, std::istream &instream, const DofEntity &dof)
Definition: blackoilbrinemodules.hh:246
static const Scalar & referencePressure(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:260
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilbrinemodules.hh:162
static Scalar computeUpdateError(const PrimaryVariables &, const EqVector &)
Return how much a Newton-Raphson update is considered an error.
Definition: blackoilbrinemodules.hh:226
Declare the properties used by the infrastructure code of the finite volume discretizations.
Definition: blackoilboundaryratevector.hh:39
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