blackoildiffusionmodule.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 OPM_BLACKOIL_DIFFUSION_MODULE_HH
29#define OPM_BLACKOIL_DIFFUSION_MODULE_HH
30
31#include <dune/common/fvector.hh>
32
33#include <opm/material/common/Valgrind.hpp>
34
37
38#include <array>
39#include <stdexcept>
40
41namespace Opm {
42
49template <class TypeTag, bool enableDiffusion>
51
52template <class TypeTag, bool enableDiffusion>
54
58template <class TypeTag>
59class BlackOilDiffusionModule<TypeTag, /*enableDiffusion=*/false>
60{
62
63public:
64 #if HAVE_ECL_INPUT
68 static void initFromState(const EclipseState&)
69 {}
70 #endif
71
75 static void registerParameters()
76 {}
77
82 template <class Context>
83 static void addDiffusiveFlux(RateVector&,
84 const Context&,
85 unsigned,
86 unsigned)
87 {}
88};
89
93template <class TypeTag>
94class BlackOilDiffusionModule<TypeTag, /*enableDiffusion=*/true>
95{
104
105 enum { numPhases = FluidSystem::numPhases };
106 enum { conti0EqIdx = Indices::conti0EqIdx };
107
108 enum { enableBioeffects = getPropValue<TypeTag, Properties::EnableBioeffects>() };
109 enum { enableMICP = Indices::enableMICP };
110
111 static constexpr unsigned contiMicrobialEqIdx = Indices::contiMicrobialEqIdx;
112 static constexpr unsigned contiOxygenEqIdx = Indices::contiOxygenEqIdx;
113 static constexpr unsigned waterPhaseIdx = FluidSystem::waterPhaseIdx;
114 static constexpr unsigned contiUreaEqIdx = Indices::contiUreaEqIdx;
115
116 using Toolbox = MathToolbox<Evaluation>;
117
118public:
120
121 #if HAVE_ECL_INPUT
125 static void initFromState(const EclipseState& eclState)
126 {
127 use_mole_fraction_ = eclState.getTableManager().diffMoleFraction();
128 }
129 #endif
130
134 static void registerParameters()
135 {}
136
159 template <class Context>
160 static void addDiffusiveFlux(RateVector& flux, const Context& context,
161 unsigned spaceIdx, unsigned timeIdx)
162 {
163 // Only work if diffusion is enabled run-time by DIFFUSE in the deck
164 if (!FluidSystem::enableDiffusion()) {
165 return;
166 }
167
168 const auto& extQuants = context.extensiveQuantities(spaceIdx, timeIdx);
169 const auto& inIq = context.intensiveQuantities(extQuants.interiorIndex(), timeIdx);
170 const auto& exIq = context.intensiveQuantities(extQuants.exteriorIndex(), timeIdx);
171 const auto& diffusivity = extQuants.diffusivity();
172 const auto& effectiveDiffusionCoefficient = extQuants.effectiveDiffusionCoefficient();
173 addDiffusiveFlux(flux, inIq, exIq, diffusivity, effectiveDiffusionCoefficient);
174 if constexpr(enableBioeffects) {
175 const auto& effectiveBioDiffCoefficient = extQuants.effectiveBioDiffCoefficient();
176 addBioDiffFlux(flux, inIq, exIq, diffusivity, effectiveBioDiffCoefficient);
177 }
178 }
179
180 template<class IntensiveQuantities,class EvaluationArray>
181 static void addDiffusiveFlux(RateVector& flux,
182 const IntensiveQuantities& inIq,
183 const IntensiveQuantities& exIq,
184 const Evaluation& diffusivity,
185 const EvaluationArray& effectiveDiffusionCoefficient)
186 {
187 const auto& inFs = inIq.fluidState();
188 const auto& exFs = exIq.fluidState();
189 unsigned pvtRegionIndex = inFs.pvtRegionIndex();
190 Evaluation diffR = 0.0;
191
192 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
193 if (!FluidSystem::phaseIsActive(phaseIdx)) {
194 continue;
195 }
196
197 // no diffusion in water for blackoil models
198 if (!FluidSystem::enableDissolvedGasInWater() && FluidSystem::waterPhaseIdx == phaseIdx) {
199 continue;
200 }
201
202 // no diffusion in gas phase in water + gas system.
203 if (FluidSystem::gasPhaseIdx == phaseIdx && !FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
204 continue;
205 }
206
207 // arithmetic mean of the phase's b factor weighed by saturation
208 Evaluation bSAvg = inFs.saturation(phaseIdx) * inFs.invB(phaseIdx);
209 bSAvg += Toolbox::value(exFs.saturation(phaseIdx)) * Toolbox::value(exFs.invB(phaseIdx));
210 bSAvg /= 2;
211
212 // phase not present, skip
213 if (bSAvg < 1.0e-6) {
214 continue;
215 }
216 Evaluation convFactor = 1.0;
217 if (FluidSystem::enableDissolvedGas() &&
218 FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) &&
219 phaseIdx == FluidSystem::oilPhaseIdx)
220 {
221 const Evaluation rsAvg = (inFs.Rs() + Toolbox::value(exFs.Rs())) / 2;
222 convFactor = 1.0 / (toFractionGasOil(pvtRegionIndex) + rsAvg);
223 diffR = inFs.Rs() - Toolbox::value(exFs.Rs());
224 }
225 if (FluidSystem::enableVaporizedOil() &&
226 FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) &&
227 phaseIdx == FluidSystem::gasPhaseIdx)
228 {
229 const Evaluation rvAvg = (inFs.Rv() + Toolbox::value(exFs.Rv())) / 2;
230 convFactor = toFractionGasOil(pvtRegionIndex) /
231 (1.0 + rvAvg * toFractionGasOil(pvtRegionIndex));
232 diffR = inFs.Rv() - Toolbox::value(exFs.Rv());
233 }
234 if (FluidSystem::enableDissolvedGasInWater() && phaseIdx == FluidSystem::waterPhaseIdx) {
235 const Evaluation rsAvg = (inFs.Rsw() + Toolbox::value(exFs.Rsw())) / 2;
236 convFactor = 1.0 / (toFractionGasWater(pvtRegionIndex) + rsAvg);
237 diffR = inFs.Rsw() - Toolbox::value(exFs.Rsw());
238 }
239 if (FluidSystem::enableVaporizedWater() && phaseIdx == FluidSystem::gasPhaseIdx) {
240 const Evaluation rvAvg = (inFs.Rvw() + Toolbox::value(exFs.Rvw())) / 2;
241 convFactor = toFractionGasWater(pvtRegionIndex) /
242 (1.0 + rvAvg * toFractionGasWater(pvtRegionIndex));
243 diffR = inFs.Rvw() - Toolbox::value(exFs.Rvw());
244 }
245
246 // mass flux of solvent component (oil in oil or gas in gas)
247 const unsigned solventCompIdx = FluidSystem::solventComponentIndex(phaseIdx);
248 const unsigned activeSolventCompIdx = FluidSystem::canonicalToActiveCompIdx(solventCompIdx);
249 flux[conti0EqIdx + activeSolventCompIdx] +=
250 -bSAvg *
251 convFactor *
252 diffR *
253 diffusivity *
254 effectiveDiffusionCoefficient[phaseIdx][solventCompIdx];
255
256 // mass flux of solute component (gas in oil or oil in gas)
257 const unsigned soluteCompIdx = FluidSystem::soluteComponentIndex(phaseIdx);
258 const unsigned activeSoluteCompIdx = FluidSystem::canonicalToActiveCompIdx(soluteCompIdx);
259 flux[conti0EqIdx + activeSoluteCompIdx] +=
260 bSAvg *
261 diffR *
262 convFactor *
263 diffusivity *
264 effectiveDiffusionCoefficient[phaseIdx][soluteCompIdx];
265 }
266 }
267
268 template<class IntensiveQuantities,class EvaluationArray>
269 static void addBioDiffFlux(RateVector& flux,
270 const IntensiveQuantities& inIq,
271 const IntensiveQuantities& exIq,
272 const Evaluation& diffusivity,
273 const EvaluationArray& effectiveBioDiffCoefficient)
274 {
275 const auto& inFs = inIq.fluidState();
276 const auto& exFs = exIq.fluidState();
277 Evaluation diffR = 0.0;
278
279 // The diffusion coefficients are given for mass concentrations
280 Evaluation bAvg = (inFs.saturation(waterPhaseIdx) * inFs.invB(waterPhaseIdx) +
281 Toolbox::value(exFs.saturation(waterPhaseIdx)) * Toolbox::value(exFs.invB(waterPhaseIdx))) / 2;
282 diffR = inIq.microbialConcentration() - Toolbox::value(exIq.microbialConcentration());
283 flux[contiMicrobialEqIdx] +=
284 bAvg *
285 diffR *
286 diffusivity *
287 effectiveBioDiffCoefficient[BioeffectsParams::micrDiffIdx];
288 if constexpr(enableMICP) {
289 diffR = inIq.oxygenConcentration() - Toolbox::value(exIq.oxygenConcentration());
290 flux[contiOxygenEqIdx] +=
291 bAvg *
292 diffR *
293 diffusivity *
294 effectiveBioDiffCoefficient[BioeffectsParams::oxygDiffIdx];
295 diffR = inIq.ureaConcentration() - Toolbox::value(exIq.ureaConcentration());
296 flux[contiUreaEqIdx] +=
297 bAvg *
298 diffR *
299 diffusivity *
300 effectiveBioDiffCoefficient[BioeffectsParams::ureaDiffIdx];
301 }
302 }
303
304private:
305 static Scalar toFractionGasOil (unsigned regionIdx)
306 {
307 const Scalar mMOil = use_mole_fraction_ ? FluidSystem::molarMass(FluidSystem::oilCompIdx, regionIdx) : 1;
308 const Scalar rhoO = FluidSystem::referenceDensity(FluidSystem::oilPhaseIdx, regionIdx);
309 const Scalar mMGas = use_mole_fraction_ ? FluidSystem::molarMass(FluidSystem::gasCompIdx, regionIdx) : 1;
310 const Scalar rhoG = FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx, regionIdx);
311 return rhoO * mMGas / (rhoG * mMOil);
312 }
313
314 static Scalar toFractionGasWater (unsigned regionIdx)
315 {
316 const Scalar mMWater = use_mole_fraction_ ? FluidSystem::molarMass(FluidSystem::waterCompIdx, regionIdx) : 1;
317 const Scalar rhoW = FluidSystem::referenceDensity(FluidSystem::waterPhaseIdx, regionIdx);
318 const Scalar mMGas = use_mole_fraction_ ? FluidSystem::molarMass(FluidSystem::gasCompIdx, regionIdx) : 1;
319 const Scalar rhoG = FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx, regionIdx);
320 return rhoW * mMGas / (rhoG * mMWater);
321 }
322
323 static bool use_mole_fraction_;
324};
325
326template <typename TypeTag>
327bool
328BlackOilDiffusionModule<TypeTag, true>::use_mole_fraction_;
329
337template <class TypeTag, bool enableDiffusion>
339
343template <class TypeTag>
344class BlackOilDiffusionIntensiveQuantities<TypeTag, /*enableDiffusion=*/false>
345{
348
349public:
354 Scalar tortuosity(unsigned) const
355 {
356 throw std::logic_error("Method tortuosity() does not make sense "
357 "if diffusion is disabled");
358 }
359
364 Scalar diffusionCoefficient(unsigned, unsigned) const
365 {
366 throw std::logic_error("Method diffusionCoefficient() does not "
367 "make sense if diffusion is disabled");
368 }
369
374 Scalar effectiveDiffusionCoefficient(unsigned, unsigned) const
375 {
376 throw std::logic_error("Method effectiveDiffusionCoefficient() "
377 "does not make sense if diffusion is disabled");
378 }
379
380protected:
385 template <class FluidState>
386 void update_(FluidState&,
387 const unsigned,
388 const ElementContext&,
389 unsigned,
390 unsigned)
391 {}
392};
393
397template <class TypeTag>
398class BlackOilDiffusionIntensiveQuantities<TypeTag, /*enableDiffusion=*/true>
399{
407
408 enum { numPhases = FluidSystem::numPhases };
409 enum { numComponents = FluidSystem::numComponents };
410 enum { enableBioeffects = getPropValue<TypeTag, Properties::EnableBioeffects>() };
411 enum { enableMICP = Indices::enableMICP };
412 enum { numBioInWat = Indices::numBioInWat };
413
414 static constexpr unsigned waterPhaseIdx = FluidSystem::waterPhaseIdx;
415
416public:
420
422
425 {
426 if (this == &rhs) {
427 return *this;
428 }
429
430 if (FluidSystem::enableDiffusion()) {
431 tortuosity_ = rhs.tortuosity_;
432 diffusionCoefficient_ = rhs.diffusionCoefficient_;
433 }
434 return *this;
435 }
436
441 Evaluation diffusionCoefficient(unsigned phaseIdx, unsigned compIdx) const
442 { return diffusionCoefficient_[phaseIdx][compIdx]; }
443
448 Evaluation tortuosity(unsigned phaseIdx) const
449 { return tortuosity_[phaseIdx]; }
450
455 Evaluation effectiveDiffusionCoefficient(unsigned phaseIdx, unsigned compIdx) const
456 {
457 // For the blackoil model tortuosity is disabled.
458 // TODO add a run-time parameter to enable tortuosity
459 static constexpr bool enableTortuosity = false;
460 if constexpr (enableTortuosity) {
461 return tortuosity_[phaseIdx] * diffusionCoefficient_[phaseIdx][compIdx];
462 } else {
463 return diffusionCoefficient_[phaseIdx][compIdx];
464 }
465 }
466
471 Evaluation bioDiffCoefficient(unsigned compIdx) const
472 { return bioDiffCoefficient_[compIdx]; }
473
478 Evaluation effectiveBioDiffCoefficient(unsigned compIdx) const
479 {
480 // TODO add a run-time parameter to enable tortuosity
481 static bool enableTortuosity = false;
482 if (enableTortuosity)
483 return tortuosity_[waterPhaseIdx] * bioDiffCoefficient_[compIdx];
484
485 return bioDiffCoefficient_[compIdx];
486 }
487
488protected:
493 template <class FluidState>
494 void update_(FluidState& fluidState,
495 const unsigned regionIdx,
496 const ElementContext& elemCtx,
497 unsigned dofIdx,
498 unsigned timeIdx)
499 {
500 // Only work if diffusion is enabled run-time by DIFFUSE in the deck
501 if (!FluidSystem::enableDiffusion()) {
502 return;
503 }
504
505 const auto& intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
506 update_(fluidState, regionIdx, intQuants);
507 }
508
509 template<class FluidState>
510 void update_(FluidState& fluidState,
511 const unsigned regionIdx,
512 const IntensiveQuantities& intQuants)
513 {
514 using Toolbox = MathToolbox<Evaluation>;
515
516 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
517 if (!FluidSystem::phaseIsActive(phaseIdx)) {
518 continue;
519 }
520
521 // no diffusion in water for blackoil models
522 if (!FluidSystem::enableDissolvedGasInWater() && FluidSystem::waterPhaseIdx == phaseIdx) {
523 continue;
524 }
525
526 // Based on Millington, R. J., & Quirk, J. P. (1961).
527 // \Note: it is possible to use NumericalConstants later
528 // constexpr auto& numconst = GetPropValue<TypeTag, Properties::NumericalConstants>;
529 constexpr double myeps = 0.0001; //numconst.blackoildiffusionmoduleeps;
530 const Evaluation& base =
531 Toolbox::max(myeps, //0.0001,
532 intQuants.porosity() *
533 intQuants.fluidState().saturation(phaseIdx));
534 tortuosity_[phaseIdx] =
535 1.0 / (intQuants.porosity() * intQuants.porosity()) *
536 Toolbox::pow(base, 10.0 / 3.0);
537
538 using PCache = typename FluidSystem::template ParameterCache<Scalar>;
539 PCache pcache(regionIdx);
540 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
541 diffusionCoefficient_[phaseIdx][compIdx] =
542 FluidSystem::diffusionCoefficient(fluidState,
543 pcache,
544 phaseIdx,
545 compIdx);
546 }
547 }
548
549 if constexpr(enableBioeffects) {
550 unsigned pvtRegionIndex = intQuants.pvtRegionIndex();
551 for (unsigned compIdx = 0; compIdx < numBioInWat; ++compIdx) {
552 bioDiffCoefficient_[compIdx] =
553 BioeffectsModule::bioDiffCoefficient(pvtRegionIndex, compIdx);
554 }
555 }
556 }
557
558private:
559 std::array<Evaluation, numPhases> tortuosity_{};
560 std::array<std::array<Evaluation, numComponents>, numPhases> diffusionCoefficient_{};
561 std::array<Evaluation, numBioInWat> bioDiffCoefficient_{};
562};
563
570template <class TypeTag, bool enableDiffusion>
571class BlackOilDiffusionExtensiveQuantities;
572
576template <class TypeTag>
577class BlackOilDiffusionExtensiveQuantities<TypeTag, /*enableDiffusion=*/false>
578{
582
583protected:
588 void update_(const ElementContext&,
589 unsigned,
590 unsigned)
591 {}
592
593 template <class Context, class FluidState>
594 void updateBoundary_(const Context&,
595 unsigned,
596 unsigned,
597 const FluidState&)
598 {}
599
600public:
605 Scalar diffusivity() const
606 {
607 throw std::logic_error("The method diffusivity() does not "
608 "make sense if diffusion is disabled.");
609 }
610
618 const Evaluation& effectiveDiffusionCoefficient(unsigned,
619 unsigned) const
620 {
621 throw std::logic_error("The method effectiveDiffusionCoefficient() "
622 "does not make sense if diffusion is disabled.");
623 }
624};
625
629template <class TypeTag>
630class BlackOilDiffusionExtensiveQuantities<TypeTag, /*enableDiffusion=*/true>
631{
638
639 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
640 enum { numComponents = getPropValue<TypeTag, Properties::NumComponents>() };
641 enum { enableBioeffects = getPropValue<TypeTag, Properties::EnableBioeffects>() };
642 enum { enableMICP = Indices::enableMICP };
643 enum { numBioInWat = Indices::numBioInWat };
644
645public:
646 using EvaluationArray = std::array<std::array<Evaluation, numComponents>, numPhases>;
647
648protected:
653 void update_(const ElementContext& elemCtx, unsigned faceIdx, unsigned timeIdx)
654 {
655 // Only work if diffusion is enabled run-time by DIFFUSE in the deck
656 if (!FluidSystem::enableDiffusion()) {
657 return;
658 }
659
660 const auto& stencil = elemCtx.stencil(timeIdx);
661 const auto& face = stencil.interiorFace(faceIdx);
662 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, timeIdx);
663 const auto& intQuantsInside = elemCtx.intensiveQuantities(extQuants.interiorIndex(), timeIdx);
664 const auto& intQuantsOutside = elemCtx.intensiveQuantities(extQuants.exteriorIndex(), timeIdx);
665
666 const Scalar diff = elemCtx.problem().diffusivity(elemCtx, face.interiorIndex(), face.exteriorIndex());
667 const Scalar faceArea = face.area();
668 diffusivity_ = diff / faceArea;
669 update(effectiveDiffusionCoefficient_, intQuantsInside, intQuantsOutside);
670 Valgrind::CheckDefined(diffusivity_);
671 if constexpr(enableBioeffects) {
672 updateBio(effectiveBioDiffCoefficient_, intQuantsInside, intQuantsOutside);
673 }
674 }
675
676public:
677 static void update(EvaluationArray& effectiveDiffusionCoefficient,
678 const IntensiveQuantities& intQuantsInside,
679 const IntensiveQuantities& intQuantsOutside)
680 {
681 // opm-models expects per area flux
682 // use the arithmetic average for the effective
683 // diffusion coefficients.
684 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
685 if (!FluidSystem::phaseIsActive(phaseIdx)) {
686 continue;
687 }
688 // no diffusion in water for blackoil models
689 if (!FluidSystem::enableDissolvedGasInWater() && FluidSystem::waterPhaseIdx == phaseIdx) {
690 continue;
691 }
692 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
693 effectiveDiffusionCoefficient[phaseIdx][compIdx] =
694 0.5 * (intQuantsInside.effectiveDiffusionCoefficient(phaseIdx, compIdx) +
695 intQuantsOutside.effectiveDiffusionCoefficient(phaseIdx, compIdx));
696 Valgrind::CheckDefined(effectiveDiffusionCoefficient[phaseIdx][compIdx]);
697 }
698 }
699 }
700
701 static void updateBio(std::array<Evaluation, numBioInWat>& effectiveBioDiffCoefficient,
702 const IntensiveQuantities& intQuantsInside,
703 const IntensiveQuantities& intQuantsOutside) {
704
705 for (unsigned compIdx = 0; compIdx < numBioInWat; ++compIdx) {
706 effectiveBioDiffCoefficient[compIdx] =
707 0.5 * (intQuantsInside.effectiveBioDiffCoefficient(compIdx) +
708 intQuantsOutside.effectiveBioDiffCoefficient(compIdx));
709 Valgrind::CheckDefined(effectiveBioDiffCoefficient[compIdx]);
710 }
711 }
712
713protected:
714 template <class Context, class FluidState>
715 void updateBoundary_(const Context&,
716 unsigned,
717 unsigned,
718 const FluidState&)
719 {
720 throw std::runtime_error("Not implemented: Diffusion across boundary not implemented for blackoil");
721 }
722
723public:
730 Scalar diffusivity() const
731 { return diffusivity_; }
732
740 const Evaluation& effectiveDiffusionCoefficient(unsigned phaseIdx, unsigned compIdx) const
741 { return effectiveDiffusionCoefficient_[phaseIdx][compIdx]; }
742
744 { return effectiveDiffusionCoefficient_; }
745
746 const Evaluation& effectiveBioDiffCoefficient(unsigned compIdx) const
747 { return effectiveBioDiffCoefficient_[compIdx]; }
748
749 const auto& effectiveBioDiffCoefficient() const{
750 return effectiveBioDiffCoefficient_;
751 }
752
753private:
754 Scalar diffusivity_;
755 EvaluationArray effectiveDiffusionCoefficient_;
756 std::array<Evaluation, numBioInWat> effectiveBioDiffCoefficient_;
757};
758
759} // namespace Opm
760
761#endif
Contains the classes required to extend the black-oil model by bioeffects.
Contains the high level supplements required to extend the black oil model by bioeffects.
Definition: blackoilbioeffectsmodules.hh:93
Scalar diffusivity() const
The diffusivity the face.
Definition: blackoildiffusionmodule.hh:605
const Evaluation & effectiveDiffusionCoefficient(unsigned, unsigned) const
The effective diffusion coeffcient of a component in a fluid phase at the face's integration point.
Definition: blackoildiffusionmodule.hh:618
void update_(const ElementContext &, unsigned, unsigned)
Update the quantities required to calculate the diffusive mass fluxes.
Definition: blackoildiffusionmodule.hh:588
void updateBoundary_(const Context &, unsigned, unsigned, const FluidState &)
Definition: blackoildiffusionmodule.hh:594
Provides the quantities required to calculate diffusive mass fluxes.
Definition: blackoildiffusionmodule.hh:631
Scalar diffusivity() const
The diffusivity of the face.
Definition: blackoildiffusionmodule.hh:730
static void update(EvaluationArray &effectiveDiffusionCoefficient, const IntensiveQuantities &intQuantsInside, const IntensiveQuantities &intQuantsOutside)
Definition: blackoildiffusionmodule.hh:677
const Evaluation & effectiveDiffusionCoefficient(unsigned phaseIdx, unsigned compIdx) const
The effective diffusion coefficient of a component in a fluid phase at the face's integration point.
Definition: blackoildiffusionmodule.hh:740
const Evaluation & effectiveBioDiffCoefficient(unsigned compIdx) const
Definition: blackoildiffusionmodule.hh:746
static void updateBio(std::array< Evaluation, numBioInWat > &effectiveBioDiffCoefficient, const IntensiveQuantities &intQuantsInside, const IntensiveQuantities &intQuantsOutside)
Definition: blackoildiffusionmodule.hh:701
std::array< std::array< Evaluation, numComponents >, numPhases > EvaluationArray
Definition: blackoildiffusionmodule.hh:646
const auto & effectiveBioDiffCoefficient() const
Definition: blackoildiffusionmodule.hh:749
const auto & effectiveDiffusionCoefficient() const
Definition: blackoildiffusionmodule.hh:743
void updateBoundary_(const Context &, unsigned, unsigned, const FluidState &)
Definition: blackoildiffusionmodule.hh:715
void update_(const ElementContext &elemCtx, unsigned faceIdx, unsigned timeIdx)
Update the quantities required to calculate the diffusive mass fluxes.
Definition: blackoildiffusionmodule.hh:653
Provides the quantities required to calculate diffusive mass fluxes.
Definition: blackoildiffusionmodule.hh:53
void update_(FluidState &, const unsigned, const ElementContext &, unsigned, unsigned)
Update the quantities required to calculate diffusive mass fluxes.
Definition: blackoildiffusionmodule.hh:386
Scalar effectiveDiffusionCoefficient(unsigned, unsigned) const
Returns the effective molecular diffusion coefficient of the porous medium for a component in a phase...
Definition: blackoildiffusionmodule.hh:374
Scalar diffusionCoefficient(unsigned, unsigned) const
Returns the molecular diffusion coefficient for a component in a phase.
Definition: blackoildiffusionmodule.hh:364
Scalar tortuosity(unsigned) const
Returns the tortuousity of the sub-domain of a fluid phase in the porous medium.
Definition: blackoildiffusionmodule.hh:354
Evaluation bioDiffCoefficient(unsigned compIdx) const
Returns the molecular diffusion coefficient for a biocomponent in the water phase.
Definition: blackoildiffusionmodule.hh:471
Evaluation effectiveBioDiffCoefficient(unsigned compIdx) const
Returns the effective molecular diffusion coefficient of the porous medium for a biocomponent in the ...
Definition: blackoildiffusionmodule.hh:478
Evaluation diffusionCoefficient(unsigned phaseIdx, unsigned compIdx) const
Returns the molecular diffusion coefficient for a component in a phase.
Definition: blackoildiffusionmodule.hh:441
Evaluation tortuosity(unsigned phaseIdx) const
Returns the tortuousity of the sub-domain of a fluid phase in the porous medium.
Definition: blackoildiffusionmodule.hh:448
BlackOilDiffusionIntensiveQuantities(BlackOilDiffusionIntensiveQuantities &&) noexcept=default
Evaluation effectiveDiffusionCoefficient(unsigned phaseIdx, unsigned compIdx) const
Returns the effective molecular diffusion coefficient of the porous medium for a component in a phase...
Definition: blackoildiffusionmodule.hh:455
void update_(FluidState &fluidState, const unsigned regionIdx, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Update the quantities required to calculate diffusive mass fluxes.
Definition: blackoildiffusionmodule.hh:494
void update_(FluidState &fluidState, const unsigned regionIdx, const IntensiveQuantities &intQuants)
Definition: blackoildiffusionmodule.hh:510
Provides the volumetric quantities required for the calculation of molecular diffusive fluxes.
Definition: blackoildiffusionmodule.hh:338
static void registerParameters()
Register all run-time parameters for the diffusion module.
Definition: blackoildiffusionmodule.hh:75
static void addDiffusiveFlux(RateVector &, const Context &, unsigned, unsigned)
Adds the diffusive mass flux flux to the flux vector over a flux integration point.
Definition: blackoildiffusionmodule.hh:83
static void addDiffusiveFlux(RateVector &flux, const Context &context, unsigned spaceIdx, unsigned timeIdx)
Adds the mass flux due to molecular diffusion to the flux vector over the integration point....
Definition: blackoildiffusionmodule.hh:160
static void addDiffusiveFlux(RateVector &flux, const IntensiveQuantities &inIq, const IntensiveQuantities &exIq, const Evaluation &diffusivity, const EvaluationArray &effectiveDiffusionCoefficient)
Definition: blackoildiffusionmodule.hh:181
static void addBioDiffFlux(RateVector &flux, const IntensiveQuantities &inIq, const IntensiveQuantities &exIq, const Evaluation &diffusivity, const EvaluationArray &effectiveBioDiffCoefficient)
Definition: blackoildiffusionmodule.hh:269
static void registerParameters()
Register all run-time parameters for the diffusion module.
Definition: blackoildiffusionmodule.hh:134
Provides the auxiliary methods required for consideration of the diffusion equation.
Definition: blackoildiffusionmodule.hh:50
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 BlackOilBioeffectsModule class.
Definition: blackoilbioeffectsparams.hpp:44