28#ifndef EWOMS_DISPERSION_MODULE_HH
29#define EWOMS_DISPERSION_MODULE_HH
33#include <opm/material/common/Valgrind.hpp>
35#include <dune/common/fvector.hh>
47template <
class TypeTag,
bool enableDispersion>
50template <
class TypeTag,
bool enableDispersion>
56template <
class TypeTag>
64 enum { numPhases = FluidSystem::numPhases };
79 template <
class Context>
86 template<
class Flu
idState,
class Scalar>
87 static void addDispersiveFlux(RateVector&,
98template <
class TypeTag>
113 enum { numPhases = FluidSystem::numPhases };
114 enum { numComponents = FluidSystem::numComponents };
115 enum { conti0EqIdx = Indices::conti0EqIdx };
125 if (!eclState.getSimulationConfig().rock_config().dispersion()) {
129 if (eclState.getSimulationConfig().hasVAPWAT() || eclState.getSimulationConfig().hasVAPOIL()) {
130 OpmLog::warning(
"Dispersion is activated in combination with VAPWAT/VAPOIL. \n"
131 "Water/oil is still allowed to vaporize, but dispersion in the "
132 "gas phase is ignored.");
140 template <
class Context>
145 if (!context.simulator().vanguard().eclState().getSimulationConfig().rock_config().dispersion()) {
151 const auto& dispersivity =
extQuants.dispersivity();
152 const auto& normVelocityAvg =
extQuants.normVelocityAvg();
176 template<
class Flu
idState,
class Scalar>
180 const Evaluation& dispersivity,
181 const Scalar& normVelocityAvg)
183 unsigned pvtRegionIndex =
fluidStateI.pvtRegionIndex();
185 if (!FluidSystem::phaseIsActive(
phaseIdx)) {
190 if (!FluidSystem::enableDissolvedGasInWater() && FluidSystem::waterPhaseIdx ==
phaseIdx) {
197 if (FluidSystem::gasPhaseIdx ==
phaseIdx) {
202 if ((!FluidSystem::enableVaporizedWater() && !FluidSystem::enableVaporizedOil()) && FluidSystem::gasPhaseIdx ==
phaseIdx) {
212 Evaluation
diffR = 0.0;
213 if (FluidSystem::enableDissolvedGas() && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) &&
phaseIdx == FluidSystem::oilPhaseIdx) {
218 if (FluidSystem::enableVaporizedOil() && FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) &&
phaseIdx == FluidSystem::gasPhaseIdx) {
220 convFactor = toMassFractionGasOil(pvtRegionIndex) / (1.0 +
rvAvg*toMassFractionGasOil(pvtRegionIndex));
223 if (FluidSystem::enableDissolvedGasInWater() &&
phaseIdx == FluidSystem::waterPhaseIdx) {
228 if (FluidSystem::enableVaporizedWater() &&
phaseIdx == FluidSystem::gasPhaseIdx) {
230 convFactor = toMassFractionGasWater(pvtRegionIndex)/ (1.0 +
rvAvg*toMassFractionGasWater(pvtRegionIndex));
258 static Scalar toMassFractionGasOil (
unsigned regionIdx) {
259 Scalar
rhoO = FluidSystem::referenceDensity(FluidSystem::oilPhaseIdx,
regionIdx);
260 Scalar
rhoG = FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx,
regionIdx);
263 static Scalar toMassFractionGasWater (
unsigned regionIdx) {
264 Scalar
rhoW = FluidSystem::referenceDensity(FluidSystem::waterPhaseIdx,
regionIdx);
265 Scalar
rhoG = FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx,
regionIdx);
277template <
class TypeTag,
bool enableDispersion>
283template <
class TypeTag>
296 throw std::logic_error(
"Method normVelocityCell() "
297 "does not make sense if dispersion is disabled");
305 template<
class ElementContext>
315template <
class TypeTag>
324 enum { numPhases = FluidSystem::numPhases };
325 enum { numComponents = FluidSystem::numComponents };
326 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
327 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
328 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
329 enum { gasCompIdx = FluidSystem::gasCompIdx };
330 enum { oilCompIdx = FluidSystem::oilCompIdx };
331 enum { waterCompIdx = FluidSystem::waterCompIdx };
332 enum { conti0EqIdx = Indices::conti0EqIdx };
355 template<
class ElementContext>
359 if (!
elemCtx.simulator().vanguard().eclState().getSimulationConfig().rock_config().dispersion()) {
362 const auto& problem =
elemCtx.simulator().problem();
363 if (problem.model().linearizer().getVelocityInfo().empty()) {
366 const std::array<int, 3>
phaseIdxs = { gasPhaseIdx, oilPhaseIdx, waterPhaseIdx };
367 const std::array<int, 3>
compIdxs = { gasCompIdx, oilCompIdx, waterCompIdx };
368 const auto&
velocityInf = problem.model().linearizer().getVelocityInfo();
371 for (
unsigned i = 0; i <
phaseIdxs.size(); ++i) {
372 normVelocityCell_[i] = 0;
375 for (
unsigned i = 0; i <
phaseIdxs.size(); ++i) {
376 if (FluidSystem::phaseIsActive(
phaseIdxs[i])) {
379 + Indices::canonicalToActiveComponentIndex(
compIdxs[i])] ));
386 Scalar normVelocityCell_[numPhases];
395template <
class TypeTag,
bool enableDispersion>
396class BlackOilDispersionExtensiveQuantities;
401template <
class TypeTag>
410 enum { numPhases = FluidSystem::numPhases };
422 template <
class Context,
class Flu
idState>
423 void updateBoundary_(
const Context&,
430 using ScalarArray = Scalar[numPhases];
432 static void update(ScalarArray&,
433 const IntensiveQuantities&,
434 const IntensiveQuantities&)
443 throw std::logic_error(
"The method dispersivity() does not "
444 "make sense if dispersion is disabled.");
456 throw std::logic_error(
"The method normVelocityAvg() "
457 "does not make sense if dispersion is disabled.");
465template <
class TypeTag>
476 enum { dimWorld = GridView::dimensionworld };
480 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
481 using DimEvalVector = Dune::FieldVector<Evaluation, dimWorld>;
483 using ScalarArray = Scalar[numPhases];
484 static void update(ScalarArray& normVelocityAvg,
489 if (!FluidSystem::phaseIsActive(
phaseIdx)) {
493 if (!FluidSystem::enableDissolvedGasInWater() && FluidSystem::waterPhaseIdx ==
phaseIdx) {
497 if ((!FluidSystem::enableVaporizedWater() && !FluidSystem::enableVaporizedOil()) && FluidSystem::gasPhaseIdx ==
phaseIdx) {
505 Valgrind::CheckDefined(normVelocityAvg[
phaseIdx]);
509 template <
class Context,
class Flu
idState>
510 void updateBoundary_(
const Context&,
515 throw std::runtime_error(
"Not implemented: Dispersion across boundary not implemented for blackoil");
526 {
return dispersivity_; }
536 {
return normVelocityAvg_[
phaseIdx]; }
538 const auto& normVelocityAvg()
const{
539 return normVelocityAvg_;
543 Scalar dispersivity_;
544 ScalarArray normVelocityAvg_;
Definition blackoildispersionmodule.hh:403
void update_(const ElementContext &, unsigned, unsigned)
Update the quantities required to calculate the dispersive fluxes.
Definition blackoildispersionmodule.hh:417
const Scalar & dispersivity() const
The dispersivity the face.
Definition blackoildispersionmodule.hh:441
const Scalar & normVelocityAvg(unsigned) const
The effective filter velocity coefficient in a fluid phase at the face's integration point.
Definition blackoildispersionmodule.hh:454
Provides the quantities required to calculate dispersive mass fluxes.
Definition blackoildispersionmodule.hh:467
const Scalar & dispersivity() const
The dispersivity of the face.
Definition blackoildispersionmodule.hh:525
const Scalar & normVelocityAvg(unsigned phaseIdx) const
The effective velocity coefficient in a fluid phase at the face's integration point.
Definition blackoildispersionmodule.hh:535
Provides the quantities required to calculate dispersive mass fluxes.
Definition blackoildispersionmodule.hh:51
Scalar normVelocityCell(unsigned, unsigned) const
Returns the max.
Definition blackoildispersionmodule.hh:294
void update_(ElementContext &, unsigned, unsigned)
Update the quantities required to calculate dispersive fluxes.
Definition blackoildispersionmodule.hh:306
void update_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Update the quantities required to calculate dispersive mass fluxes.
Definition blackoildispersionmodule.hh:356
Scalar normVelocityCell(unsigned phaseIdx) const
Returns the max.
Definition blackoildispersionmodule.hh:339
Provides the volumetric quantities required for the calculation of dispersive fluxes.
Definition blackoildispersionmodule.hh:278
static void addDispersiveFlux(RateVector &, const Context &, unsigned, unsigned)
Adds the dispersive flux to the flux vector over a flux integration point.
Definition blackoildispersionmodule.hh:80
static void addDispersiveFlux(RateVector &flux, const FluidState &fluidStateI, const FluidState &fluidStateJ, const Evaluation &dispersivity, const Scalar &normVelocityAvg)
Adds the mass flux due to dispersion to the flux vector over the integration point.
Definition blackoildispersionmodule.hh:177
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:141
Provides the auxiliary methods required for consideration of the dispersion equation.
Definition blackoildispersionmodule.hh:48
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