31#ifndef EWOMS_TRANS_FLUX_MODULE_HH
32#define EWOMS_TRANS_FLUX_MODULE_HH
36#include <opm/material/common/Valgrind.hpp>
37#include <dune/common/fvector.hh>
38#include <dune/common/fmatrix.hh>
42template <
class TypeTag>
43class TransIntensiveQuantities;
45template <
class TypeTag>
46class TransExtensiveQuantities;
48template <
class TypeTag>
49class TransBaseProblem;
54template <
class TypeTag>
71template <
class TypeTag>
78template <
class TypeTag>
83 void update_(
const ElementContext&,
unsigned,
unsigned)
90template <
class TypeTag>
103 enum { dimWorld = GridView::dimensionworld };
104 enum { numPhases = FluidSystem::numPhases };
107 typedef Dune::FieldVector<Scalar, dimWorld> DimVector;
108 typedef Dune::FieldVector<Evaluation, dimWorld> EvalDimVector;
109 typedef Dune::FieldMatrix<Scalar, dimWorld, dimWorld> DimMatrix;
117 throw std::logic_error(
"The ECL transmissibility module does not provide an explicit intrinsic permeability");
128 throw std::logic_error(
"The ECL transmissibility module does not provide explicit potential gradients");
138 {
return pressureDifference_[
phaseIdx]; }
148 throw std::logic_error(
"The ECL transmissibility module does not provide explicit filter velocities");
192 void updateSolvent(
const ElementContext& elemCtx,
unsigned scvfIdx,
unsigned timeIdx)
195 void updatePolymer(
const ElementContext& elemCtx,
unsigned scvfIdx,
unsigned timeIdx)
203 Valgrind::SetUndefined(*
this);
206 static const bool isEcfv = std::is_same<Discretization, EcfvDiscretization<TypeTag> >::value;
210 const auto& stencil = elemCtx.stencil(
timeIdx);
213 interiorDofIdx_ =
scvf.interiorIndex();
214 exteriorDofIdx_ =
scvf.exteriorIndex();
215 assert(interiorDofIdx_ != exteriorDofIdx_);
217 unsigned I = stencil.globalSpaceIndex(interiorDofIdx_);
218 unsigned J = stencil.globalSpaceIndex(exteriorDofIdx_);
225 Scalar
g = elemCtx.problem().gravity()[dimWorld - 1];
230 Scalar
zIn = dofCenterDepth_(elemCtx, interiorDofIdx_,
timeIdx);
231 Scalar
zEx = dofCenterDepth_(elemCtx, exteriorDofIdx_,
timeIdx);
237 if (!FluidSystem::phaseIsActive(
phaseIdx))
247 pressureDifference_[
phaseIdx] = 0.0;
269 if (pressureDifference_[
phaseIdx] > 0.0) {
273 else if (pressureDifference_[
phaseIdx] < 0.0) {
280 Scalar Vin = elemCtx.dofVolume(interiorDofIdx_, 0);
281 Scalar Vex = elemCtx.dofVolume(exteriorDofIdx_, 0);
286 else if (Vin < Vex) {
323 template <
class Flu
idState>
327 const FluidState& exFluidState)
329 const auto& stencil = elemCtx.stencil(
timeIdx);
332 interiorDofIdx_ =
scvf.interiorIndex();
334 Scalar trans = transmissibilityBoundary_(elemCtx,
scvfIdx,
timeIdx);
339 Scalar
g = elemCtx.problem().gravity()[dimWorld - 1];
348 Scalar
zIn = dofCenterDepth_(elemCtx, interiorDofIdx_,
timeIdx);
349 Scalar
zEx =
scvf.integrationPos()[dimWorld - 1];
356 if (!FluidSystem::phaseIsActive(
phaseIdx))
375 if (pressureDifference_[
phaseIdx] > 0.0) {
400 elemCtx.problem().materialLawParams(elemCtx,
403 std::array<typename FluidState::Scalar,numPhases>
kr;
404 MaterialLaw::relativePermeabilities(
kr,
matParams, exFluidState);
420 void calculateBoundaryFluxes_(
const ElementContext&,
unsigned,
unsigned)
425 Scalar transmissibility_(
const ElementContext& elemCtx,
unsigned scvfIdx,
unsigned timeIdx)
const
427 const auto& stencil = elemCtx.stencil(
timeIdx);
429 const auto&
interiorPos = stencil.subControlVolume(
face.interiorIndex()).globalPos();
430 const auto&
exteriorPos = stencil.subControlVolume(
face.exteriorIndex()).globalPos();
438 const auto&
K0mat = elemCtx.problem().intrinsicPermeability(elemCtx,
face.interiorIndex(),
timeIdx);
439 const auto&
K1mat = elemCtx.problem().intrinsicPermeability(elemCtx,
face.exteriorIndex(),
timeIdx);
445 for (
unsigned i = 0; i < dimWorld; ++ i){
446 if (std::abs(
face.normal()[i]) >
val) {
447 val = std::abs(
face.normal()[i]);
457 Scalar transmissibilityBoundary_(
const ElementContext& elemCtx,
unsigned scvfIdx,
unsigned timeIdx)
const
459 const auto& stencil = elemCtx.stencil(
timeIdx);
461 const auto&
interiorPos = stencil.subControlVolume(
face.interiorIndex()).globalPos();
465 const auto&
K0mat = elemCtx.problem().intrinsicPermeability(elemCtx,
face.interiorIndex(),
timeIdx);
471 for (
unsigned i = 0; i < dimWorld; ++ i){
472 if (std::abs(
face.normal()[i]) >
val) {
473 val = std::abs(
face.normal()[i]);
482 template <
class Context>
483 Scalar dofCenterDepth_(
const Context& context,
unsigned spaceIdx,
unsigned timeIdx)
const
486 return pos[dimWorld-1];
489 Implementation& asImp_()
490 {
return *
static_cast<Implementation*
>(
this); }
492 const Implementation& asImp_()
const
493 {
return *
static_cast<const Implementation*
>(
this); }
496 Evaluation volumeFlux_[numPhases];
500 Evaluation pressureDifference_[numPhases];
503 unsigned short interiorDofIdx_{};
504 unsigned short exteriorDofIdx_{};
505 short upIdx_[numPhases]{};
506 short dnIdx_[numPhases]{};
Provides the defaults for the parameters required by the transmissibility based volume flux calculati...
Definition transfluxmodule.hh:73
Provides the transmissibility based flux module.
Definition transfluxmodule.hh:92
const Evaluation & volumeFlux(unsigned phaseIdx) const
Return the volume flux of a fluid phase at the face's integration point .
Definition transfluxmodule.hh:160
unsigned downstreamIndex_(unsigned phaseIdx) const
Returns the local index of the degree of freedom in which is in downstream direction.
Definition transfluxmodule.hh:185
const EvalDimVector & filterVelocity(unsigned) const
Return the filter velocity of a fluid phase at the face's integration point [m/s].
Definition transfluxmodule.hh:146
void calculateGradients_(const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Update the required gradients for interior faces.
Definition transfluxmodule.hh:201
const Evaluation & pressureDifference(unsigned phaseIdx) const
Return the gravity corrected pressure difference between the interior and the exterior of a face.
Definition transfluxmodule.hh:137
void calculateFluxes_(const ElementContext &, unsigned, unsigned)
Update the volumetric fluxes for all fluid phases on the interior faces of the context.
Definition transfluxmodule.hh:417
const EvalDimVector & potentialGrad(unsigned) const
Return the pressure potential gradient of a fluid phase at the face's integration point [Pa/m].
Definition transfluxmodule.hh:126
void calculateBoundaryGradients_(const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx, const FluidState &exFluidState)
Update the required gradients for boundary faces.
Definition transfluxmodule.hh:324
const DimMatrix & intrinsicPermeability() const
Return the intrinsic permeability tensor at a face [m^2].
Definition transfluxmodule.hh:115
unsigned upstreamIndex_(unsigned phaseIdx) const
Returns the local index of the degree of freedom in which is in upstream direction.
Definition transfluxmodule.hh:171
Provides the intensive quantities for the transmissibility based flux module.
Definition transfluxmodule.hh:80
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
Specifies a flux module which uses transmissibilities.
Definition transfluxmodule.hh:56
static void registerParameters()
Register all run-time parameters for the flux module.
Definition transfluxmodule.hh:63