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 "blackoilproperties.hh"
32
34
35#if HAVE_ECL_INPUT
36#include <opm/input/eclipse/EclipseState/EclipseState.hpp>
37#include <opm/input/eclipse/EclipseState/Tables/PvtwsaltTable.hpp>
38#include <opm/input/eclipse/EclipseState/Tables/PcfactTable.hpp>
39#include <opm/input/eclipse/EclipseState/Tables/PermfactTable.hpp>
40#include <opm/input/eclipse/EclipseState/Tables/SaltSolubilityTable.hpp>
41#include <opm/input/eclipse/EclipseState/Tables/TableManager.hpp>
42#include <opm/input/eclipse/EclipseState/Tables/SimpleTable.hpp>
43#endif
44
45#include <dune/common/fvector.hh>
46
47#include <string>
48#include <math.h>
49
50namespace Opm {
56template <class TypeTag, bool enableBrineV = getPropValue<TypeTag, Properties::EnableBrine>()>
58{
71
72 using Toolbox = MathToolbox<Evaluation>;
73
74 using TabulatedFunction = typename BlackOilBrineParams<Scalar>::TabulatedFunction;
75
76 static constexpr unsigned saltConcentrationIdx = Indices::saltConcentrationIdx;
77 static constexpr unsigned contiBrineEqIdx = Indices::contiBrineEqIdx;
78 static constexpr unsigned waterPhaseIdx = FluidSystem::waterPhaseIdx;
79 static constexpr bool gasEnabled = Indices::gasEnabled;
80 static constexpr bool oilEnabled = Indices::oilEnabled;
81 static constexpr unsigned enableBrine = enableBrineV;
82 static constexpr unsigned enableSaltPrecipitation = getPropValue<TypeTag, Properties::EnableSaltPrecipitation>();
83
84 static constexpr unsigned numEq = getPropValue<TypeTag, Properties::NumEq>();
85 static constexpr unsigned numPhases = FluidSystem::numPhases;
86
87public:
88
89#if HAVE_ECL_INPUT
93 static void initFromState(const EclipseState& eclState)
94 {
95 // some sanity checks: if brine are enabled, the BRINE keyword must be
96 // present, if brine are disabled the keyword must not be present.
97 if (enableBrine && !eclState.runspec().phases().active(Phase::BRINE)) {
98 throw std::runtime_error("Non-trivial brine treatment requested at compile time, but "
99 "the deck does not contain the BRINE keyword");
100 }
101 else if (!enableBrine && eclState.runspec().phases().active(Phase::BRINE)) {
102 throw std::runtime_error("Brine treatment disabled at compile time, but the deck "
103 "contains the BRINE keyword");
104 }
105
106 if (!eclState.runspec().phases().active(Phase::BRINE))
107 return; // brine treatment is supposed to be disabled
108
109 const auto& tableManager = eclState.getTableManager();
110
111 unsigned numPvtRegions = tableManager.getTabdims().getNumPVTTables();
112 params_.referencePressure_.resize(numPvtRegions);
113
114 const auto& pvtwsaltTables = tableManager.getPvtwSaltTables();
115
116 // initialize the objects which deal with the BDENSITY keyword
117 const auto& bdensityTables = tableManager.getBrineDensityTables();
118 if (!bdensityTables.empty()) {
119 params_.bdensityTable_.resize(numPvtRegions);
120 assert(numPvtRegions == bdensityTables.size());
121 for (unsigned pvtRegionIdx = 0; pvtRegionIdx < numPvtRegions; ++ pvtRegionIdx) {
122 const auto& bdensityTable = bdensityTables[pvtRegionIdx];
123 const auto& pvtwsaltTable = pvtwsaltTables[pvtRegionIdx];
124 const auto& c = pvtwsaltTable.getSaltConcentrationColumn();
125 params_.bdensityTable_[pvtRegionIdx].setXYContainers(c, bdensityTable);
126 }
127 }
128
129 if constexpr (enableSaltPrecipitation) {
130 const TableContainer& permfactTables = tableManager.getPermfactTables();
131 params_.permfactTable_.resize(numPvtRegions);
132 for (size_t i = 0; i < permfactTables.size(); ++i) {
133 const PermfactTable& permfactTable = permfactTables.getTable<PermfactTable>(i);
134 params_.permfactTable_[i].setXYContainers(permfactTable.getPorosityChangeColumn(), permfactTable.getPermeabilityMultiplierColumn());
135 }
136
137 const TableContainer& saltsolTables = tableManager.getSaltsolTables();
138 if (!saltsolTables.empty()) {
139 params_.saltsolTable_.resize(numPvtRegions);
140 params_.saltdenTable_.resize(numPvtRegions);
141 assert(numPvtRegions == saltsolTables.size());
142 for (unsigned pvtRegionIdx = 0; pvtRegionIdx < numPvtRegions; ++ pvtRegionIdx) {
143 const SaltsolTable& saltsolTable = saltsolTables.getTable<SaltsolTable>(pvtRegionIdx );
144 params_.saltsolTable_[pvtRegionIdx] = saltsolTable.getSaltsolColumn().front();
145 params_.saltdenTable_[pvtRegionIdx] = saltsolTable.getSaltdenColumn().front();
146 }
147 }
148
149 const TableContainer& pcfactTables = tableManager.getPcfactTables();
150 if (!pcfactTables.empty()) {
151 unsigned numSatRegions = tableManager.getTabdims().getNumSatTables();
152 params_.pcfactTable_.resize(numSatRegions);
153 for (size_t i = 0; i < pcfactTables.size(); ++i) {
154 const PcfactTable& pcfactTable = pcfactTables.getTable<PcfactTable>(i);
155 params_.pcfactTable_[i].setXYContainers(pcfactTable.getPorosityChangeColumn(), pcfactTable.getPcMultiplierColumn());
156 }
157 }
158 }
159 }
160#endif
161
165 static void registerParameters()
166 {
167 }
168
169 static bool primaryVarApplies(unsigned pvIdx)
170 {
171 if constexpr (enableBrine)
172 return pvIdx == saltConcentrationIdx;
173 else
174 return false;
175 }
176
180 template <class FluidState>
181 static void assignPrimaryVars(PrimaryVariables& priVars,
182 const FluidState& fluidState)
183 {
184 if constexpr (enableBrine)
185 priVars[saltConcentrationIdx] = fluidState.saltConcentration();
186 }
187
188 static std::string primaryVarName([[maybe_unused]] unsigned pvIdx)
189 {
190 assert(primaryVarApplies(pvIdx));
191
192 return "saltConcentration";
193 }
194
195 static Scalar primaryVarWeight([[maybe_unused]] unsigned pvIdx)
196 {
197 assert(primaryVarApplies(pvIdx));
198
199 // TODO: it may be beneficial to chose this differently.
200 return static_cast<Scalar>(1.0);
201 }
202
203 static bool eqApplies(unsigned eqIdx)
204 {
205 if constexpr (enableBrine)
206 return eqIdx == contiBrineEqIdx;
207 else
208 return false;
209 }
210
211 static std::string eqName([[maybe_unused]] unsigned eqIdx)
212 {
213 assert(eqApplies(eqIdx));
214
215 return "conti^brine";
216 }
217
218 static Scalar eqWeight([[maybe_unused]] unsigned eqIdx)
219 {
220 assert(eqApplies(eqIdx));
221
222 // TODO: it may be beneficial to chose this differently.
223 return static_cast<Scalar>(1.0);
224 }
225
226 // must be called after water storage is computed
227 template <class LhsEval>
228 static void addStorage(Dune::FieldVector<LhsEval, numEq>& storage,
229 const IntensiveQuantities& intQuants)
230 {
231 if constexpr (enableBrine) {
232 const auto& fs = intQuants.fluidState();
233
234 LhsEval surfaceVolumeWater =
235 Toolbox::template decay<LhsEval>(fs.saturation(waterPhaseIdx))
236 * Toolbox::template decay<LhsEval>(fs.invB(waterPhaseIdx))
237 * Toolbox::template decay<LhsEval>(intQuants.porosity());
238
239 // avoid singular matrix if no water is present.
240 surfaceVolumeWater = max(surfaceVolumeWater, 1e-10);
241
242 // Brine in water phase
243 const LhsEval massBrine = surfaceVolumeWater
244 * Toolbox::template decay<LhsEval>(fs.saltConcentration());
245
246 if (enableSaltPrecipitation){
247 double saltDensity = intQuants.saltDensity(); // Solid salt density kg/m3
248 const LhsEval solidSalt =
249 Toolbox::template decay<LhsEval>(intQuants.porosity())
250 / (1.0 - Toolbox::template decay<LhsEval>(intQuants.saltSaturation()) + 1.e-8)
251 * saltDensity
252 * Toolbox::template decay<LhsEval>(intQuants.saltSaturation());
253
254 storage[contiBrineEqIdx] += massBrine + solidSalt;
255 }
256 else { storage[contiBrineEqIdx] += massBrine;}
257 }
258 }
259
260 static void computeFlux([[maybe_unused]] RateVector& flux,
261 [[maybe_unused]] const ElementContext& elemCtx,
262 [[maybe_unused]] unsigned scvfIdx,
263 [[maybe_unused]] unsigned timeIdx)
264
265 {
266 if constexpr (enableBrine) {
267 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
268
269 const unsigned upIdx = extQuants.upstreamIndex(FluidSystem::waterPhaseIdx);
270 const unsigned inIdx = extQuants.interiorIndex();
271 const auto& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
272
273 if (upIdx == inIdx) {
274 flux[contiBrineEqIdx] =
275 extQuants.volumeFlux(waterPhaseIdx)
276 *up.fluidState().invB(waterPhaseIdx)
277 *up.fluidState().saltConcentration();
278 }
279 else {
280 flux[contiBrineEqIdx] =
281 extQuants.volumeFlux(waterPhaseIdx)
282 *decay<Scalar>(up.fluidState().invB(waterPhaseIdx))
283 *decay<Scalar>(up.fluidState().saltConcentration());
284 }
285 }
286 }
287
291 static Scalar computeUpdateError(const PrimaryVariables&,
292 const EqVector&)
293 {
294 // do not consider consider the change of Brine primary variables for
295 // convergence
296 // TODO: maybe this should be changed
297 return static_cast<Scalar>(0.0);
298 }
299
300 template <class DofEntity>
301 static void serializeEntity(const Model& model, std::ostream& outstream, const DofEntity& dof)
302 {
303 if constexpr (enableBrine) {
304 unsigned dofIdx = model.dofMapper().index(dof);
305 const PrimaryVariables& priVars = model.solution(/*timeIdx=*/0)[dofIdx];
306 outstream << priVars[saltConcentrationIdx];
307 }
308 }
309
310 template <class DofEntity>
311 static void deserializeEntity(Model& model, std::istream& instream, const DofEntity& dof)
312 {
313 if constexpr (enableBrine) {
314 unsigned dofIdx = model.dofMapper().index(dof);
315 PrimaryVariables& priVars0 = model.solution(/*timeIdx=*/0)[dofIdx];
316 PrimaryVariables& priVars1 = model.solution(/*timeIdx=*/1)[dofIdx];
317
318 instream >> priVars0[saltConcentrationIdx];
319
320 // set the primary variables for the beginning of the current time step.
321 priVars1[saltConcentrationIdx] = priVars0[saltConcentrationIdx];
322 }
323 }
324
325 static const Scalar& referencePressure(const ElementContext& elemCtx,
326 unsigned scvIdx,
327 unsigned timeIdx)
328 {
329 unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
330 return params_.referencePressure_[pvtnumRegionIdx];
331 }
332
333
334 static const TabulatedFunction& bdensityTable(const ElementContext& elemCtx,
335 unsigned scvIdx,
336 unsigned timeIdx)
337 {
338 unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
339 return params_.bdensityTable_[pvtnumRegionIdx];
340 }
341
342 static const TabulatedFunction& pcfactTable(unsigned satnumRegionIdx)
343 {
344 return params_.pcfactTable_[satnumRegionIdx];
345 }
346
347 static const TabulatedFunction& permfactTable(const ElementContext& elemCtx,
348 unsigned scvIdx,
349 unsigned timeIdx)
350 {
351 unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
352 return params_.permfactTable_[pvtnumRegionIdx];
353 }
354
355 static const TabulatedFunction& permfactTable(unsigned pvtnumRegionIdx)
356 {
357 return params_.permfactTable_[pvtnumRegionIdx];
358 }
359
360 static const Scalar saltsolTable(const ElementContext& elemCtx,
361 unsigned scvIdx,
362 unsigned timeIdx)
363 {
364 unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
365 return params_.saltsolTable_[pvtnumRegionIdx];
366 }
367
368 static const Scalar saltdenTable(const ElementContext& elemCtx,
369 unsigned scvIdx,
370 unsigned timeIdx)
371 {
372 unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
373 return params_.saltdenTable_[pvtnumRegionIdx];
374 }
375
376 static bool hasBDensityTables()
377 {
378 return !params_.bdensityTable_.empty();
379 }
380
381 static bool hasSaltsolTables()
382 {
383 return !params_.saltsolTable_.empty();
384 }
385
386 static bool hasPcfactTables()
387 {
388 if constexpr (enableSaltPrecipitation)
389 return !params_.pcfactTable_.empty();
390 else
391 return false;
392 }
393
394 static Scalar saltSol(unsigned regionIdx) {
395 return params_.saltsolTable_[regionIdx];
396 }
397
398private:
399 static BlackOilBrineParams<Scalar> params_;
400};
401
402
403template <class TypeTag, bool enableBrineV>
404BlackOilBrineParams<typename BlackOilBrineModule<TypeTag, enableBrineV>::Scalar>
405BlackOilBrineModule<TypeTag, enableBrineV>::params_;
406
414template <class TypeTag, bool enableBrineV = getPropValue<TypeTag, Properties::EnableBrine>()>
416{
418
426
428
429 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
430 static constexpr int saltConcentrationIdx = Indices::saltConcentrationIdx;
431 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
432 static constexpr int gasPhaseIdx = FluidSystem::gasPhaseIdx;
433 static constexpr int oilPhaseIdx = FluidSystem::oilPhaseIdx;
434 static constexpr unsigned enableBrine = enableBrineV;
435 static constexpr unsigned enableSaltPrecipitation = getPropValue<TypeTag, Properties::EnableSaltPrecipitation>();
436 static constexpr int contiBrineEqIdx = Indices::contiBrineEqIdx;
437
438public:
439
445 void updateSaltConcentration_(const ElementContext& elemCtx,
446 unsigned dofIdx,
447 unsigned timeIdx)
448 {
449 const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
450
451 auto& fs = asImp_().fluidState_;
452
453 if constexpr (enableSaltPrecipitation) {
454 const auto& saltsolTable = BrineModule::saltsolTable(elemCtx, dofIdx, timeIdx);
455 saltSolubility_ = saltsolTable;
456
457 const auto& saltdenTable = BrineModule::saltdenTable(elemCtx, dofIdx, timeIdx);
458 saltDensity_ = saltdenTable;
459
460 if (priVars.primaryVarsMeaningBrine() == PrimaryVariables::BrineMeaning::Sp) {
461 saltSaturation_ = priVars.makeEvaluation(saltConcentrationIdx, timeIdx);
462 fs.setSaltConcentration(saltSolubility_);
463 }
464 else {
465 saltConcentration_ = priVars.makeEvaluation(saltConcentrationIdx, timeIdx);
466 fs.setSaltConcentration(saltConcentration_);
467 saltSaturation_ = 0.0;
468 }
469 fs.setSaltSaturation(saltSaturation_);
470 }
471 else {
472 saltConcentration_ = priVars.makeEvaluation(saltConcentrationIdx, timeIdx);
473 fs.setSaltConcentration(saltConcentration_);
474 }
475 }
476
477 void saltPropertiesUpdate_([[maybe_unused]] const ElementContext& elemCtx,
478 [[maybe_unused]] unsigned dofIdx,
479 [[maybe_unused]] unsigned timeIdx)
480 {
481 if constexpr (enableSaltPrecipitation) {
482 const Evaluation porosityFactor = min(1.0 - saltSaturation(), 1.0); //phi/phi_0
483
484 const auto& permfactTable = BrineModule::permfactTable(elemCtx, dofIdx, timeIdx);
485
486 permFactor_ = permfactTable.eval(porosityFactor);
487 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
488 if (!FluidSystem::phaseIsActive(phaseIdx))
489 continue;
490
491 asImp_().mobility_[phaseIdx] *= permFactor_;
492 }
493 }
494 }
495
496 const Evaluation& saltConcentration() const
497 { return saltConcentration_; }
498
499 const Evaluation& brineRefDensity() const
500 { return refDensity_; }
501
502 const Evaluation& saltSaturation() const
503 { return saltSaturation_; }
504
505 Scalar saltSolubility() const
506 { return saltSolubility_; }
507
508 Scalar saltDensity() const
509 { return saltDensity_; }
510
511 const Evaluation& permFactor() const
512 { return permFactor_; }
513
514protected:
515 Implementation& asImp_()
516 { return *static_cast<Implementation*>(this); }
517
519 Evaluation refDensity_;
520 Evaluation saltSaturation_;
521 Evaluation permFactor_;
524
525};
526
527template <class TypeTag>
529{
533
534public:
535 void updateSaltConcentration_(const ElementContext&,
536 unsigned,
537 unsigned)
538 { }
539
540 void saltPropertiesUpdate_(const ElementContext&,
541 unsigned,
542 unsigned)
543 { }
544
545 const Evaluation& saltConcentration() const
546 { throw std::runtime_error("saltConcentration() called but brine are disabled"); }
547
548 const Evaluation& brineRefDensity() const
549 { throw std::runtime_error("brineRefDensity() called but brine are disabled"); }
550
551 const Evaluation& saltSaturation() const
552 { throw std::logic_error("saltSaturation() called but salt precipitation is disabled"); }
553
554 const Scalar saltSolubility() const
555 { throw std::logic_error("saltSolubility() called but salt precipitation is disabled"); }
556
557 const Scalar saltDensity() const
558 { throw std::logic_error("saltDensity() called but salt precipitation is disabled"); }
559
560 const Evaluation& permFactor() const
561 { throw std::logic_error("permFactor() called but salt precipitation is disabled"); }
562
563};
564
565} // namespace Ewoms
566
567#endif
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:535
const Scalar saltDensity() const
Definition: blackoilbrinemodules.hh:557
const Evaluation & permFactor() const
Definition: blackoilbrinemodules.hh:560
const Evaluation & saltSaturation() const
Definition: blackoilbrinemodules.hh:551
void saltPropertiesUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilbrinemodules.hh:540
const Scalar saltSolubility() const
Definition: blackoilbrinemodules.hh:554
const Evaluation & brineRefDensity() const
Definition: blackoilbrinemodules.hh:548
const Evaluation & saltConcentration() const
Definition: blackoilbrinemodules.hh:545
Definition: blackoilbrinemodules.hh:416
void updateSaltConcentration_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Update the intensive properties needed to handle brine from the primary variables.
Definition: blackoilbrinemodules.hh:445
const Evaluation & brineRefDensity() const
Definition: blackoilbrinemodules.hh:499
const Evaluation & saltSaturation() const
Definition: blackoilbrinemodules.hh:502
Evaluation permFactor_
Definition: blackoilbrinemodules.hh:521
Evaluation refDensity_
Definition: blackoilbrinemodules.hh:519
Scalar saltDensity_
Definition: blackoilbrinemodules.hh:523
const Evaluation & saltConcentration() const
Definition: blackoilbrinemodules.hh:496
Scalar saltSolubility_
Definition: blackoilbrinemodules.hh:522
Implementation & asImp_()
Definition: blackoilbrinemodules.hh:515
Evaluation saltSaturation_
Definition: blackoilbrinemodules.hh:520
Scalar saltDensity() const
Definition: blackoilbrinemodules.hh:508
Evaluation saltConcentration_
Definition: blackoilbrinemodules.hh:518
const Evaluation & permFactor() const
Definition: blackoilbrinemodules.hh:511
Scalar saltSolubility() const
Definition: blackoilbrinemodules.hh:505
void saltPropertiesUpdate_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:477
Contains the high level supplements required to extend the black oil model by brine.
Definition: blackoilbrinemodules.hh:58
static bool hasBDensityTables()
Definition: blackoilbrinemodules.hh:376
static const TabulatedFunction & bdensityTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:334
static void serializeEntity(const Model &model, std::ostream &outstream, const DofEntity &dof)
Definition: blackoilbrinemodules.hh:301
static const TabulatedFunction & permfactTable(unsigned pvtnumRegionIdx)
Definition: blackoilbrinemodules.hh:355
static Scalar saltSol(unsigned regionIdx)
Definition: blackoilbrinemodules.hh:394
static std::string eqName(unsigned eqIdx)
Definition: blackoilbrinemodules.hh:211
static bool primaryVarApplies(unsigned pvIdx)
Definition: blackoilbrinemodules.hh:169
static const TabulatedFunction & pcfactTable(unsigned satnumRegionIdx)
Definition: blackoilbrinemodules.hh:342
static const Scalar saltdenTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:368
static const TabulatedFunction & permfactTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:347
static bool hasPcfactTables()
Definition: blackoilbrinemodules.hh:386
static std::string primaryVarName(unsigned pvIdx)
Definition: blackoilbrinemodules.hh:188
static void registerParameters()
Register all run-time parameters for the black-oil brine module.
Definition: blackoilbrinemodules.hh:165
static bool hasSaltsolTables()
Definition: blackoilbrinemodules.hh:381
static Scalar eqWeight(unsigned eqIdx)
Definition: blackoilbrinemodules.hh:218
static const Scalar saltsolTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:360
static Scalar primaryVarWeight(unsigned pvIdx)
Definition: blackoilbrinemodules.hh:195
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:260
static bool eqApplies(unsigned eqIdx)
Definition: blackoilbrinemodules.hh:203
static void assignPrimaryVars(PrimaryVariables &priVars, const FluidState &fluidState)
Assign the brine specific primary variables to a PrimaryVariables object.
Definition: blackoilbrinemodules.hh:181
static void deserializeEntity(Model &model, std::istream &instream, const DofEntity &dof)
Definition: blackoilbrinemodules.hh:311
static const Scalar & referencePressure(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:325
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
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:291
Definition: blackoilboundaryratevector.hh:37
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:235
Struct holding the parameters for the BlackoilBrineModule class.
Definition: blackoilbrineparams.hh:39
Tabulated1DFunction< Scalar > TabulatedFunction
Definition: blackoilbrineparams.hh:40