28#ifndef EWOMS_DISPERSION_MODULE_HH
29#define EWOMS_DISPERSION_MODULE_HH
31#include <dune/common/fvector.hh>
36#include <opm/material/common/MathToolbox.hpp>
37#include <opm/material/common/Valgrind.hpp>
49template <
class TypeTag,
bool enableDispersion>
52template <
class TypeTag,
bool enableDispersion>
58template <
class TypeTag>
67 enum { numPhases = FluidSystem::numPhases };
82 template <
class Context>
89 template<
class IntensiveQuantities,
class Scalar>
90 static void addDispersiveFlux(RateVector&,
91 const IntensiveQuantities&,
92 const IntensiveQuantities&,
101template <
class TypeTag>
116 enum { numPhases = FluidSystem::numPhases };
117 enum { numComponents = FluidSystem::numComponents };
118 enum { conti0EqIdx = Indices::conti0EqIdx };
122 static constexpr unsigned contiMicrobialEqIdx = Indices::contiMicrobialEqIdx;
123 static constexpr unsigned contiOxygenEqIdx = Indices::contiOxygenEqIdx;
124 static constexpr unsigned waterPhaseIdx = FluidSystem::waterPhaseIdx;
125 static constexpr unsigned contiUreaEqIdx = Indices::contiUreaEqIdx;
134 if (!eclState.getSimulationConfig().rock_config().dispersion()) {
138 if (eclState.getSimulationConfig().hasVAPWAT() || eclState.getSimulationConfig().hasVAPOIL()) {
139 OpmLog::warning(
"Dispersion is activated in combination with VAPWAT/VAPOIL. \n"
140 "Water/oil is still allowed to vaporize, but dispersion in the "
141 "gas phase is ignored.");
149 template <
class Context>
154 if (!context.simulator().vanguard().eclState().getSimulationConfig().rock_config().dispersion()) {
160 const auto& dispersivity =
extQuants.dispersivity();
161 const auto& normVelocityAvg =
extQuants.normVelocityAvg();
162 addDispersiveFlux(
flux,
inIq,
exIq, dispersivity, normVelocityAvg);
185 template<
class IntensiveQuantities,
class Scalar>
187 const IntensiveQuantities&
inIq,
188 const IntensiveQuantities&
exIq,
189 const Evaluation& dispersivity,
190 const Scalar& normVelocityAvg)
192 const auto&
inFs =
inIq.fluidState();
193 const auto&
exFs =
exIq.fluidState();
194 Evaluation
diffR = 0.0;
195 if constexpr(enableMICP) {
197 Evaluation
bAvg = (
inFs.invB(waterPhaseIdx) + Toolbox::value(
exFs.invB(waterPhaseIdx))) / 2;
198 diffR =
inIq.microbialConcentration()- Toolbox::value(
exIq.microbialConcentration());
199 flux[contiMicrobialEqIdx] +=
201 * normVelocityAvg[waterPhaseIdx]
204 diffR =
inIq.oxygenConcentration()- Toolbox::value(
exIq.oxygenConcentration());
205 flux[contiOxygenEqIdx] +=
207 * normVelocityAvg[waterPhaseIdx]
210 diffR =
inIq.ureaConcentration()- Toolbox::value(
exIq.ureaConcentration());
211 flux[contiUreaEqIdx] +=
213 * normVelocityAvg[waterPhaseIdx]
219 unsigned pvtRegionIndex =
inFs.pvtRegionIndex();
221 if (!FluidSystem::phaseIsActive(
phaseIdx)) {
226 if (!FluidSystem::enableDissolvedGasInWater() && FluidSystem::waterPhaseIdx ==
phaseIdx) {
233 if (FluidSystem::gasPhaseIdx ==
phaseIdx) {
243 if (FluidSystem::enableDissolvedGas() && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) &&
phaseIdx == FluidSystem::oilPhaseIdx) {
244 Evaluation
rsAvg = (
inFs.Rs() + Toolbox::value(
exFs.Rs())) / 2;
248 if (FluidSystem::enableVaporizedOil() && FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) &&
phaseIdx == FluidSystem::gasPhaseIdx) {
249 Evaluation
rvAvg = (
inFs.Rv() + Toolbox::value(
exFs.Rv())) / 2;
250 convFactor = toMassFractionGasOil(pvtRegionIndex) / (1.0 +
rvAvg*toMassFractionGasOil(pvtRegionIndex));
253 if (FluidSystem::enableDissolvedGasInWater() &&
phaseIdx == FluidSystem::waterPhaseIdx) {
254 Evaluation
rsAvg = (
inFs.Rsw() + Toolbox::value(
exFs.Rsw())) / 2;
258 if (FluidSystem::enableVaporizedWater() &&
phaseIdx == FluidSystem::gasPhaseIdx) {
259 Evaluation
rvAvg = (
inFs.Rvw() + Toolbox::value(
exFs.Rvw())) / 2;
260 convFactor = toMassFractionGasWater(pvtRegionIndex)/ (1.0 +
rvAvg*toMassFractionGasWater(pvtRegionIndex));
288 static Scalar toMassFractionGasOil (
unsigned regionIdx) {
289 Scalar
rhoO = FluidSystem::referenceDensity(FluidSystem::oilPhaseIdx,
regionIdx);
290 Scalar
rhoG = FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx,
regionIdx);
293 static Scalar toMassFractionGasWater (
unsigned regionIdx) {
294 Scalar
rhoW = FluidSystem::referenceDensity(FluidSystem::waterPhaseIdx,
regionIdx);
295 Scalar
rhoG = FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx,
regionIdx);
307template <
class TypeTag,
bool enableDispersion>
313template <
class TypeTag>
326 throw std::logic_error(
"Method normVelocityCell() "
327 "does not make sense if dispersion is disabled");
335 template<
class ElementContext>
345template <
class TypeTag>
354 enum { numPhases = FluidSystem::numPhases };
355 enum { numComponents = FluidSystem::numComponents };
356 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
357 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
358 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
359 enum { gasCompIdx = FluidSystem::gasCompIdx };
360 enum { oilCompIdx = FluidSystem::oilCompIdx };
361 enum { waterCompIdx = FluidSystem::waterCompIdx };
362 enum { conti0EqIdx = Indices::conti0EqIdx };
385 template<
class ElementContext>
389 if (!elemCtx.simulator().vanguard().eclState().getSimulationConfig().rock_config().dispersion()) {
392 const auto& problem = elemCtx.simulator().problem();
393 if (problem.model().linearizer().getVelocityInfo().empty()) {
396 const std::array<int, 3>
phaseIdxs = { gasPhaseIdx, oilPhaseIdx, waterPhaseIdx };
397 const std::array<int, 3>
compIdxs = { gasCompIdx, oilCompIdx, waterCompIdx };
398 const auto&
velocityInf = problem.model().linearizer().getVelocityInfo();
399 unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx,
timeIdx);
401 for (
unsigned i = 0; i <
phaseIdxs.size(); ++i) {
402 normVelocityCell_[i] = 0;
405 for (
unsigned i = 0; i <
phaseIdxs.size(); ++i) {
406 if (FluidSystem::phaseIsActive(
phaseIdxs[i])) {
409 + Indices::canonicalToActiveComponentIndex(
compIdxs[i])] ));
416 Scalar normVelocityCell_[numPhases];
425template <
class TypeTag,
bool enableDispersion>
426class BlackOilDispersionExtensiveQuantities;
431template <
class TypeTag>
440 enum { numPhases = FluidSystem::numPhases };
452 template <
class Context,
class Flu
idState>
453 void updateBoundary_(
const Context&,
460 using ScalarArray = Scalar[numPhases];
462 static void update(ScalarArray&,
463 const IntensiveQuantities&,
464 const IntensiveQuantities&)
473 throw std::logic_error(
"The method dispersivity() does not "
474 "make sense if dispersion is disabled.");
486 throw std::logic_error(
"The method normVelocityAvg() "
487 "does not make sense if dispersion is disabled.");
495template <
class TypeTag>
506 enum { dimWorld = GridView::dimensionworld };
510 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
511 using DimEvalVector = Dune::FieldVector<Evaluation, dimWorld>;
513 using ScalarArray = Scalar[numPhases];
514 static void update(ScalarArray& normVelocityAvg,
519 if (!FluidSystem::phaseIsActive(
phaseIdx)) {
523 if (!FluidSystem::enableDissolvedGasInWater() && FluidSystem::waterPhaseIdx ==
phaseIdx) {
529 if (FluidSystem::gasPhaseIdx ==
phaseIdx) {
537 Valgrind::CheckDefined(normVelocityAvg[
phaseIdx]);
541 template <
class Context,
class Flu
idState>
542 void updateBoundary_(
const Context&,
547 throw std::runtime_error(
"Not implemented: Dispersion across boundary not implemented for blackoil");
558 {
return dispersivity_; }
568 {
return normVelocityAvg_[
phaseIdx]; }
570 const auto& normVelocityAvg()
const{
571 return normVelocityAvg_;
575 Scalar dispersivity_;
576 ScalarArray normVelocityAvg_;
Definition blackoildispersionmodule.hh:433
void update_(const ElementContext &, unsigned, unsigned)
Update the quantities required to calculate the dispersive fluxes.
Definition blackoildispersionmodule.hh:447
const Scalar & dispersivity() const
The dispersivity the face.
Definition blackoildispersionmodule.hh:471
const Scalar & normVelocityAvg(unsigned) const
The effective filter velocity coefficient in a fluid phase at the face's integration point.
Definition blackoildispersionmodule.hh:484
Provides the quantities required to calculate dispersive mass fluxes.
Definition blackoildispersionmodule.hh:497
const Scalar & dispersivity() const
The dispersivity of the face.
Definition blackoildispersionmodule.hh:557
const Scalar & normVelocityAvg(unsigned phaseIdx) const
The effective velocity coefficient in a fluid phase at the face's integration point.
Definition blackoildispersionmodule.hh:567
Provides the quantities required to calculate dispersive mass fluxes.
Definition blackoildispersionmodule.hh:53
Scalar normVelocityCell(unsigned, unsigned) const
Returns the max.
Definition blackoildispersionmodule.hh:324
void update_(ElementContext &, unsigned, unsigned)
Update the quantities required to calculate dispersive fluxes.
Definition blackoildispersionmodule.hh:336
void update_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Update the quantities required to calculate dispersive mass fluxes.
Definition blackoildispersionmodule.hh:386
Scalar normVelocityCell(unsigned phaseIdx) const
Returns the max.
Definition blackoildispersionmodule.hh:369
Provides the volumetric quantities required for the calculation of dispersive fluxes.
Definition blackoildispersionmodule.hh:308
static void addDispersiveFlux(RateVector &, const Context &, unsigned, unsigned)
Adds the dispersive flux to the flux vector over a flux integration point.
Definition blackoildispersionmodule.hh:83
static void addDispersiveFlux(RateVector &flux, const IntensiveQuantities &inIq, const IntensiveQuantities &exIq, const Evaluation &dispersivity, const Scalar &normVelocityAvg)
Adds the mass flux due to dispersion to the flux vector over the integration point.
Definition blackoildispersionmodule.hh:186
static void addDispersiveFlux(RateVector &flux, const Context &context, unsigned spaceIdx, unsigned timeIdx)
Adds the mass flux due to dispersion to the flux vector over the flux integration point.
Definition blackoildispersionmodule.hh:150
Provides the auxiliary methods required for consideration of the dispersion equation.
Definition blackoildispersionmodule.hh:50
Declare the properties used by the infrastructure code of the finite volume discretizations.
Defines the common properties required by the porous medium multi-phase models.
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