28#ifndef OPM_BLACKOIL_DIFFUSION_MODULE_HH
29#define OPM_BLACKOIL_DIFFUSION_MODULE_HH
34#include <opm/material/common/Valgrind.hpp>
36#include <dune/common/fvector.hh>
48template <
class TypeTag,
bool enableDiffusion>
51template <
class TypeTag,
bool enableDiffusion>
57template <
class TypeTag>
85 template <
class Context>
96template <
class TypeTag>
106 enum { numPhases = FluidSystem::numPhases };
107 enum { numComponents = FluidSystem::numComponents };
108 enum { conti0EqIdx = Indices::conti0EqIdx };
112 static constexpr unsigned contiMicrobialEqIdx = Indices::contiMicrobialEqIdx;
113 static constexpr unsigned contiOxygenEqIdx = Indices::contiOxygenEqIdx;
114 static constexpr unsigned waterPhaseIdx = FluidSystem::waterPhaseIdx;
115 static constexpr unsigned contiUreaEqIdx = Indices::contiUreaEqIdx;
129 use_mole_fraction_ = eclState.getTableManager().diffMoleFraction();
161 template <
class Context>
166 if(!FluidSystem::enableDiffusion())
171 const auto& diffusivity =
extQuants.diffusivity();
172 const auto& effectiveDiffusionCoefficient =
extQuants.effectiveDiffusionCoefficient();
173 addDiffusiveFlux(
flux,
inIq,
exIq, diffusivity, effectiveDiffusionCoefficient);
176 template<
class IntensiveQuantities,
class EvaluationArray>
177 static void addDiffusiveFlux(RateVector&
flux,
178 const IntensiveQuantities&
inIq,
179 const IntensiveQuantities&
exIq,
180 const Evaluation& diffusivity,
181 const EvaluationArray& effectiveDiffusionCoefficient)
183 const auto&
inFs =
inIq.fluidState();
184 const auto&
exFs =
exIq.fluidState();
185 Evaluation
diffR = 0.0;
186 if constexpr(enableMICP) {
188 Evaluation
bAvg = (
inFs.invB(waterPhaseIdx) + Toolbox::value(
exFs.invB(waterPhaseIdx))) / 2;
189 diffR =
inIq.microbialConcentration() - Toolbox::value(
exIq.microbialConcentration());
190 flux[contiMicrobialEqIdx] +=
194 * effectiveDiffusionCoefficient[waterPhaseIdx][contiMicrobialEqIdx - 1];
195 diffR =
inIq.oxygenConcentration() - Toolbox::value(
exIq.oxygenConcentration());
196 flux[contiOxygenEqIdx] +=
200 * effectiveDiffusionCoefficient[waterPhaseIdx][contiOxygenEqIdx - 1];
201 diffR =
inIq.ureaConcentration() - Toolbox::value(
exIq.ureaConcentration());
202 flux[contiUreaEqIdx] +=
206 * effectiveDiffusionCoefficient[waterPhaseIdx][contiUreaEqIdx - 1];
210 unsigned pvtRegionIndex =
inFs.pvtRegionIndex();
212 if (!FluidSystem::phaseIsActive(
phaseIdx)) {
217 if (!FluidSystem::enableDissolvedGasInWater() && FluidSystem::waterPhaseIdx ==
phaseIdx) {
222 if (FluidSystem::gasPhaseIdx ==
phaseIdx && !FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
235 if (FluidSystem::enableDissolvedGas() && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) &&
phaseIdx == FluidSystem::oilPhaseIdx) {
236 Evaluation
rsAvg = (
inFs.Rs() + Toolbox::value(
exFs.Rs())) / 2;
240 if (FluidSystem::enableVaporizedOil() && FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) &&
phaseIdx == FluidSystem::gasPhaseIdx) {
241 Evaluation
rvAvg = (
inFs.Rv() + Toolbox::value(
exFs.Rv())) / 2;
242 convFactor = toFractionGasOil(pvtRegionIndex) / (1.0 +
rvAvg*toFractionGasOil(pvtRegionIndex));
245 if (FluidSystem::enableDissolvedGasInWater() &&
phaseIdx == FluidSystem::waterPhaseIdx) {
246 Evaluation
rsAvg = (
inFs.Rsw() + Toolbox::value(
exFs.Rsw())) / 2;
250 if (FluidSystem::enableVaporizedWater() &&
phaseIdx == FluidSystem::gasPhaseIdx) {
251 Evaluation
rvAvg = (
inFs.Rvw() + Toolbox::value(
exFs.Rvw())) / 2;
252 convFactor = toFractionGasWater(pvtRegionIndex)/ (1.0 +
rvAvg*toFractionGasWater(pvtRegionIndex));
279 static Scalar toFractionGasOil (
unsigned regionIdx) {
280 Scalar
mMOil = use_mole_fraction_? FluidSystem::molarMass(FluidSystem::oilCompIdx,
regionIdx) : 1;
281 Scalar
rhoO = FluidSystem::referenceDensity(FluidSystem::oilPhaseIdx,
regionIdx);
282 Scalar
mMGas = use_mole_fraction_? FluidSystem::molarMass(FluidSystem::gasCompIdx,
regionIdx) : 1;
283 Scalar
rhoG = FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx,
regionIdx);
286 static Scalar toFractionGasWater (
unsigned regionIdx) {
287 Scalar
mMWater = use_mole_fraction_? FluidSystem::molarMass(FluidSystem::waterCompIdx,
regionIdx) : 1;
288 Scalar
rhoW = FluidSystem::referenceDensity(FluidSystem::waterPhaseIdx,
regionIdx);
289 Scalar
mMGas = use_mole_fraction_? FluidSystem::molarMass(FluidSystem::gasCompIdx,
regionIdx) : 1;
290 Scalar
rhoG = FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx,
regionIdx);
294 static bool use_mole_fraction_;
297template <
typename TypeTag>
299BlackOilDiffusionModule<TypeTag, true>::use_mole_fraction_;
308template <
class TypeTag,
bool enableDiffusion>
314template <
class TypeTag>
328 throw std::logic_error(
"Method tortuosity() does not make sense "
329 "if diffusion is disabled");
338 throw std::logic_error(
"Method diffusionCoefficient() does not "
339 "make sense if diffusion is disabled");
348 throw std::logic_error(
"Method effectiveDiffusionCoefficient() "
349 "does not make sense if diffusion is disabled");
357 template <
class Flu
idState>
360 const ElementContext&,
369template <
class TypeTag>
379 enum { numPhases = FluidSystem::numPhases };
380 enum { numComponents = FluidSystem::numComponents };
383 static constexpr unsigned waterPhaseIdx = FluidSystem::waterPhaseIdx;
395 if (
this == &rhs)
return *
this;
397 if (FluidSystem::enableDiffusion()) {
398 std::copy(rhs.tortuosity_, rhs.tortuosity_ + numPhases, tortuosity_);
399 for (
size_t i = 0; i < numPhases; ++i) {
400 std::copy(rhs.diffusionCoefficient_[i],
401 rhs.diffusionCoefficient_[i]+numComponents,
402 diffusionCoefficient_[i]);
441 template <
class Flu
idState>
444 const ElementContext& elemCtx,
449 if(!FluidSystem::enableDiffusion())
451 const auto& intQuants = elemCtx.intensiveQuantities(dofIdx,
timeIdx);
455 template<
class Flu
idState>
456 void update_(FluidState& fluidState,
458 const IntensiveQuantities& intQuants) {
461 if constexpr(enableMICP) {
462 unsigned pvtRegionIndex = intQuants.fluidState().pvtRegionIndex();
463 diffusionCoefficient_[waterPhaseIdx][0] = MICPModule::microbialDiffusion(pvtRegionIndex);
464 diffusionCoefficient_[waterPhaseIdx][1] = MICPModule::oxygenDiffusion(pvtRegionIndex);
465 diffusionCoefficient_[waterPhaseIdx][2] = MICPModule::ureaDiffusion(pvtRegionIndex);
470 if (!FluidSystem::phaseIsActive(
phaseIdx)) {
475 if (!FluidSystem::enableDissolvedGasInWater() && FluidSystem::waterPhaseIdx ==
phaseIdx) {
482 constexpr double myeps = 0.0001;
483 const Evaluation&
base =
486 * intQuants.fluidState().saturation(
phaseIdx));
488 1.0 / (intQuants.porosity() * intQuants.porosity())
489 * Toolbox::pow(
base, 10.0/3.0);
493 FluidSystem::diffusionCoefficient(fluidState,
502 Evaluation tortuosity_[numPhases];
503 Evaluation diffusionCoefficient_[numPhases][numComponents];
512template <
class TypeTag,
bool enableDiffusion>
513class BlackOilDiffusionExtensiveQuantities;
518template <
class TypeTag>
535 template <
class Context,
class Flu
idState>
536 void updateBoundary_(
const Context&,
549 throw std::logic_error(
"The method diffusivity() does not "
550 "make sense if diffusion is disabled.");
563 throw std::logic_error(
"The method effectiveDiffusionCoefficient() "
564 "does not make sense if diffusion is disabled.");
572template <
class TypeTag>
583 enum { dimWorld = GridView::dimensionworld };
588 static constexpr unsigned waterPhaseIdx = FluidSystem::waterPhaseIdx;
590 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
591 using DimEvalVector = Dune::FieldVector<Evaluation, dimWorld>;
593 using EvaluationArray = Evaluation[numPhases][numComponents];
602 if(!FluidSystem::enableDiffusion())
605 const auto& stencil = elemCtx.stencil(
timeIdx);
606 const auto&
face = stencil.interiorFace(faceIdx);
611 const Scalar diffusivity = elemCtx.problem().diffusivity(elemCtx,
face.interiorIndex(),
face.exteriorIndex());
612 const Scalar faceArea =
face.area();
613 diffusivity_ = diffusivity / faceArea;
615 Valgrind::CheckDefined(diffusivity_);
619 static void update(EvaluationArray& effectiveDiffusionCoefficient,
625 if constexpr(enableMICP) {
627 effectiveDiffusionCoefficient[waterPhaseIdx][
compIdx] = 0.5 *
635 if (!FluidSystem::phaseIsActive(
phaseIdx)) {
639 if (!FluidSystem::enableDissolvedGasInWater() && FluidSystem::waterPhaseIdx ==
phaseIdx) {
646 Valgrind::CheckDefined(effectiveDiffusionCoefficient[
phaseIdx][
compIdx]);
651 template <
class Context,
class Flu
idState>
652 void updateBoundary_(
const Context&,
657 throw std::runtime_error(
"Not implemented: Diffusion across boundary not implemented for blackoil");
668 {
return diffusivity_; }
680 const auto& effectiveDiffusionCoefficient()
const{
681 return effectiveDiffusionCoefficient_;
686 EvaluationArray effectiveDiffusionCoefficient_;
Contains the classes required to extend the black-oil model by MICP.
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:560
const Scalar & diffusivity() const
The diffusivity the face.
Definition blackoildiffusionmodule.hh:547
void update_(const ElementContext &, unsigned, unsigned)
Update the quantities required to calculate the diffusive mass fluxes.
Definition blackoildiffusionmodule.hh:530
Provides the quantities required to calculate diffusive mass fluxes.
Definition blackoildiffusionmodule.hh:574
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:677
void update_(const ElementContext &elemCtx, unsigned faceIdx, unsigned timeIdx)
Update the quantities required to calculate the diffusive mass fluxes.
Definition blackoildiffusionmodule.hh:599
const Scalar & diffusivity() const
The diffusivity of the face.
Definition blackoildiffusionmodule.hh:667
Provides the quantities required to calculate diffusive mass fluxes.
Definition blackoildiffusionmodule.hh:52
Scalar effectiveDiffusionCoefficient(unsigned, unsigned) const
Returns the effective molecular diffusion coefficient of the porous medium for a component in a phase...
Definition blackoildiffusionmodule.hh:346
Scalar diffusionCoefficient(unsigned, unsigned) const
Returns the molecular diffusion coefficient for a component in a phase.
Definition blackoildiffusionmodule.hh:336
void update_(FluidState &, typename FluidSystem::template ParameterCache< typename FluidState::Scalar > &, const ElementContext &, unsigned, unsigned)
Update the quantities required to calculate diffusive mass fluxes.
Definition blackoildiffusionmodule.hh:358
Scalar tortuosity(unsigned) const
Returns the tortuousity of the sub-domain of a fluid phase in the porous medium.
Definition blackoildiffusionmodule.hh:326
void update_(FluidState &fluidState, typename FluidSystem::template ParameterCache< typename FluidState::Scalar > ¶mCache, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Update the quantities required to calculate diffusive mass fluxes.
Definition blackoildiffusionmodule.hh:442
Evaluation diffusionCoefficient(unsigned phaseIdx, unsigned compIdx) const
Returns the molecular diffusion coefficient for a component in a phase.
Definition blackoildiffusionmodule.hh:411
Evaluation tortuosity(unsigned phaseIdx) const
Returns the tortuousity of the sub-domain of a fluid phase in the porous medium.
Definition blackoildiffusionmodule.hh:418
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:425
Provides the volumetric quantities required for the calculation of molecular diffusive fluxes.
Definition blackoildiffusionmodule.hh:309
static void registerParameters()
Register all run-time parameters for the diffusion module.
Definition blackoildiffusionmodule.hh:78
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:86
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:162
static void registerParameters()
Register all run-time parameters for the diffusion module.
Definition blackoildiffusionmodule.hh:136
Provides the auxiliary methods required for consideration of the diffusion equation.
Definition blackoildiffusionmodule.hh:49
Contains the high level supplements required to extend the black oil model by MICP.
Definition blackoilmicpmodules.hh:49
Declare the properties used by the infrastructure code of the finite volume discretizations.
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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