31#ifndef OPM_NEWTRAN_FLUX_MODULE_HPP
32#define OPM_NEWTRAN_FLUX_MODULE_HPP
34#include <dune/common/fvector.hh>
35#include <dune/common/fmatrix.hh>
37#include <opm/common/OpmLog/OpmLog.hpp>
39#include <opm/input/eclipse/EclipseState/Grid/FaceDir.hpp>
41#include <opm/material/common/MathToolbox.hpp>
42#include <opm/material/common/Valgrind.hpp>
53template <
class TypeTag>
54class NewTranIntensiveQuantities;
56template <
class TypeTag>
57class NewTranExtensiveQuantities;
59template <
class TypeTag>
60class NewTranBaseProblem;
66template <
class TypeTag>
85template <
class TypeTag>
93template <
class TypeTag>
98 void update_(
const ElementContext&,
unsigned,
unsigned)
106template <
class TypeTag>
119 enum { dimWorld = GridView::dimensionworld };
120 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
121 enum { numPhases = FluidSystem::numPhases };
130 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
131 using EvalDimVector = Dune::FieldVector<Evaluation, dimWorld>;
132 using DimMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
142 throw std::invalid_argument(
"The ECL transmissibility module does not provide an explicit intrinsic permeability");
153 throw std::invalid_argument(
"The ECL transmissibility module does not provide explicit potential gradients");
163 {
return pressureDifference_[
phaseIdx]; }
173 throw std::invalid_argument(
"The ECL transmissibility module does not provide explicit filter velocities");
225 static void volumeAndPhasePressureDifferences(std::array<short, numPhases>&
upIdx,
226 std::array<short, numPhases>&
dnIdx,
233 const auto& problem =
elemCtx.problem();
244 Scalar faceArea =
scvf.area();
245 Scalar thpres = problem.thresholdPressure(
I,
J);
250 Scalar
g =
elemCtx.problem().gravity()[dimWorld - 1];
271 if (!FluidSystem::phaseIsActive(
phaseIdx))
287 problem.moduleParams());
294 const auto& materialLawManager = problem.materialLawManager();
295 FaceDir::DirEnum
facedir = FaceDir::DirEnum::Unknown;
296 if (materialLawManager->hasDirectionalRelperms()) {
309 template<
class EvalType>
310 static void calculatePhasePressureDiff_(
short&
upIdx,
324 const ModuleParams& moduleParams)
344 if constexpr(enableConvectiveMixing) {
377 else if (Vin < Vex) {
419 Valgrind::SetUndefined(*
this);
421 volumeAndPhasePressureDifferences(upIdx_ , dnIdx_, volumeFlux_, pressureDifference_,
elemCtx,
scvfIdx,
timeIdx);
427 template <
class Flu
idState>
431 const FluidState& exFluidState)
434 const Scalar faceArea =
scvf.area();
435 const Scalar
zEx =
scvf.integrationPos()[dimWorld - 1];
436 const auto& problem =
elemCtx.problem();
450 pressureDifference_);
455 if constexpr (enableSolvent) {
456 if (upIdx_[gasPhaseIdx] == 0) {
462 asImp_().setSolventVolumeFlux(0.0);
471 template <
class Problem,
class Flu
idState,
class EvaluationContainer>
475 const unsigned bfIdx,
476 const double faceArea,
478 const FluidState& exFluidState,
479 std::array<short, numPhases>&
upIdx,
480 std::array<short, numPhases>&
dnIdx,
494 Scalar
g = problem.gravity()[dimWorld - 1];
508 if (!FluidSystem::phaseIsActive(
phaseIdx))
546 const Scalar
transMult = Toolbox::value(
up.rockCompTransMultiplier());
556 std::array<typename FluidState::Scalar,numPhases>
kr;
557 MaterialLaw::relativePermeabilities(
kr,
matParams, exFluidState);
574 void calculateBoundaryFluxes_(
const ElementContext&,
unsigned,
unsigned)
578 Implementation& asImp_()
579 {
return *
static_cast<Implementation*
>(
this); }
581 const Implementation& asImp_()
const
582 {
return *
static_cast<const Implementation*
>(
this); }
585 Evaluation volumeFlux_[numPhases];
589 Evaluation pressureDifference_[numPhases];
592 std::array<short, numPhases> upIdx_;
593 std::array<short, numPhases> dnIdx_;
Calculates the local residual of the black oil model.
Declares the properties required by the black oil model.
Definition blackoilconvectivemixingmodule.hh:58
Provides the defaults for the parameters required by the transmissibility based volume flux calculati...
Definition NewTranFluxModule.hpp:87
Provides the ECL flux module.
Definition NewTranFluxModule.hpp:108
const Evaluation & volumeFlux(unsigned phaseIdx) const
Return the volume flux of a fluid phase at the face's integration point .
Definition NewTranFluxModule.hpp:185
const DimMatrix & intrinsicPermeability() const
Return the intrinsic permeability tensor at a face [m^2].
Definition NewTranFluxModule.hpp:140
void calculateBoundaryGradients_(const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx, const FluidState &exFluidState)
Update the required gradients for boundary faces.
Definition NewTranFluxModule.hpp:428
static void calculateBoundaryGradients_(const Problem &problem, const unsigned globalSpaceIdx, const IntensiveQuantities &intQuantsIn, const unsigned bfIdx, const double faceArea, const double zEx, const FluidState &exFluidState, std::array< short, numPhases > &upIdx, std::array< short, numPhases > &dnIdx, EvaluationContainer &volumeFlux, EvaluationContainer &pressureDifference)
Update the required gradients for boundary faces.
Definition NewTranFluxModule.hpp:472
unsigned upstreamIndex_(unsigned phaseIdx) const
Returns the local index of the degree of freedom in which is in upstream direction.
Definition NewTranFluxModule.hpp:196
void calculateFluxes_(const ElementContext &, unsigned, unsigned)
Update the volumetric fluxes for all fluid phases on the interior faces of the context.
Definition NewTranFluxModule.hpp:571
const EvalDimVector & potentialGrad(unsigned) const
Return the pressure potential gradient of a fluid phase at the face's integration point [Pa/m].
Definition NewTranFluxModule.hpp:151
unsigned downstreamIndex_(unsigned phaseIdx) const
Returns the local index of the degree of freedom in which is in downstream direction.
Definition NewTranFluxModule.hpp:210
const EvalDimVector & filterVelocity(unsigned) const
Return the filter velocity of a fluid phase at the face's integration point [m/s].
Definition NewTranFluxModule.hpp:171
const Evaluation & pressureDifference(unsigned phaseIdx) const
Return the gravity corrected pressure difference between the interior and the exterior of a face.
Definition NewTranFluxModule.hpp:162
void calculateGradients_(const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Update the required gradients for interior faces.
Definition NewTranFluxModule.hpp:417
Provides the intensive quantities for the ECL flux module.
Definition NewTranFluxModule.hpp:95
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
Definition blackoillocalresidualtpfa.hh:133
Specifies a flux module which uses ECL transmissibilities.
Definition NewTranFluxModule.hpp:68
static void registerParameters()
Register all run-time parameters for the flux module.
Definition NewTranFluxModule.hpp:76