30#ifndef EWOMS_DARCY_FLUX_MODULE_HH
31#define EWOMS_DARCY_FLUX_MODULE_HH
33#include <dune/common/fvector.hh>
34#include <dune/common/fmatrix.hh>
36#include <opm/common/Exceptions.hpp>
38#include <opm/material/common/Valgrind.hpp>
46template <
class TypeTag>
47class DarcyIntensiveQuantities;
49template <
class TypeTag>
50class DarcyExtensiveQuantities;
52template <
class TypeTag>
53class DarcyBaseProblem;
59template <
class TypeTag>
78template <
class TypeTag>
86template <
class TypeTag>
91 void update_(
const ElementContext&,
113template <
class TypeTag>
124 enum { dimWorld = GridView::dimensionworld };
128 using ParameterCache =
typename FluidSystem::template ParameterCache<Evaluation>;
129 using EvalDimVector = Dune::FieldVector<Evaluation, dimWorld>;
130 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
131 using DimMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
148 {
return potentialGrad_[
phaseIdx]; }
157 {
return filterVelocity_[
phaseIdx]; }
172 short upstreamIndex_(
unsigned phaseIdx)
const
173 {
return upstreamDofIdx_[
phaseIdx]; }
175 short downstreamIndex_(
unsigned phaseIdx)
const
176 {
return downstreamDofIdx_[
phaseIdx]; }
193 unsigned i =
scvf.interiorIndex();
194 unsigned j =
scvf.exteriorIndex();
195 interiorDofIdx_ =
static_cast<short>(i);
196 exteriorDofIdx_ =
static_cast<short>(j);
202 Valgrind::SetUndefined(potentialGrad_[
phaseIdx]);
211 Valgrind::CheckDefined(potentialGrad_[
phaseIdx]);
215 if (Parameters::Get<Parameters::EnableGravity>()) {
244 if (std::is_same<Scalar, Evaluation>::value ||
259 if (std::is_same<Scalar, Evaluation>::value ||
284 + std::string(FluidSystem::phaseName(
phaseIdx))+
"'");
290 Valgrind::SetUndefined(K_);
292 Valgrind::CheckDefined(K_);
296 Valgrind::SetUndefined(potentialGrad_[
phaseIdx]);
301 Evaluation tmp = 0.0;
306 upstreamDofIdx_[
phaseIdx] = exteriorDofIdx_;
307 downstreamDofIdx_[
phaseIdx] = interiorDofIdx_;
310 upstreamDofIdx_[
phaseIdx] = interiorDofIdx_;
311 downstreamDofIdx_[
phaseIdx] = exteriorDofIdx_;
330 template <
class Flu
idState>
334 const FluidState& fluidState)
342 Valgrind::SetUndefined(potentialGrad_[
phaseIdx]);
351 Valgrind::CheckDefined(potentialGrad_[
phaseIdx]);
355 auto i =
scvf.interiorIndex();
356 interiorDofIdx_ =
static_cast<short>(i);
357 exteriorDofIdx_ = -1;
365 if (Parameters::Get<Parameters::EnableGravity>()) {
386 Valgrind::CheckDefined(
pStatIn);
399 Valgrind::CheckDefined(potentialGrad_[
phaseIdx]);
403 + std::string(FluidSystem::phaseName(
phaseIdx))+
"'");
414 Scalar
kr[numPhases];
415 MaterialLaw::relativePermeabilities(
kr,
matParams, fluidState);
416 Valgrind::CheckDefined(
kr);
422 Evaluation tmp = 0.0;
427 upstreamDofIdx_[
phaseIdx] = exteriorDofIdx_;
428 downstreamDofIdx_[
phaseIdx] = interiorDofIdx_;
431 upstreamDofIdx_[
phaseIdx] = interiorDofIdx_;
432 downstreamDofIdx_[
phaseIdx] = exteriorDofIdx_;
436 if (upstreamDofIdx_[
phaseIdx] < 0) {
443 / Toolbox::value(fluidState.viscosity(
phaseIdx));
449 Valgrind::CheckDefined(mobility_[
phaseIdx]);
462 const DimVector& normal =
scvf.normal();
463 Valgrind::CheckDefined(normal);
471 asImp_().calculateFilterVelocity_(
phaseIdx);
472 Valgrind::CheckDefined(filterVelocity_[
phaseIdx]);
475 for (
unsigned i = 0; i < normal.size(); ++i)
491 const DimVector& normal =
scvf.normal();
492 Valgrind::CheckDefined(normal);
501 asImp_().calculateFilterVelocity_(
phaseIdx);
502 Valgrind::CheckDefined(filterVelocity_[
phaseIdx]);
504 for (
unsigned i = 0; i < normal.size(); ++i)
509 void calculateFilterVelocity_(
unsigned phaseIdx)
513 for (
unsigned i = 0; i < K_.M(); ++ i)
514 for (
unsigned j = 0; j < K_.N(); ++ j)
515 assert(std::isfinite(K_[i][j]));
522 for (
unsigned i = 0; i < filterVelocity_[
phaseIdx].size(); ++ i)
528 Implementation& asImp_()
529 {
return *
static_cast<Implementation*
>(
this); }
531 const Implementation& asImp_()
const
532 {
return *
static_cast<const Implementation*
>(
this); }
539 Evaluation mobility_[numPhases];
542 EvalDimVector filterVelocity_[numPhases];
546 Evaluation volumeFlux_[numPhases];
549 EvalDimVector potentialGrad_[numPhases];
552 short upstreamDofIdx_[numPhases];
553 short downstreamDofIdx_[numPhases];
554 short interiorDofIdx_;
555 short exteriorDofIdx_;
Callback class for a phase pressure.
Definition quantitycallbacks.hh:133
Provides the defaults for the parameters required by the Darcy velocity approach.
Definition darcyfluxmodule.hh:80
Provides the Darcy flux module.
Definition darcyfluxmodule.hh:115
void calculateBoundaryGradients_(const ElementContext &elemCtx, unsigned boundaryFaceIdx, unsigned timeIdx, const FluidState &fluidState)
Calculate the gradients at the grid boundary which are required to determine the volumetric fluxes.
Definition darcyfluxmodule.hh:331
const EvalDimVector & filterVelocity(unsigned phaseIdx) const
Return the filter velocity of a fluid phase at the face's integration point [m/s].
Definition darcyfluxmodule.hh:156
void calculateBoundaryFluxes_(const ElementContext &elemCtx, unsigned boundaryFaceIdx, unsigned timeIdx)
Calculate the volumetric fluxes at a boundary face of all fluid phases.
Definition darcyfluxmodule.hh:486
const EvalDimVector & potentialGrad(unsigned phaseIdx) const
Return the pressure potential gradient of a fluid phase at the face's integration point [Pa/m].
Definition darcyfluxmodule.hh:147
void calculateFluxes_(const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Calculate the volumetric fluxes of all phases.
Definition darcyfluxmodule.hh:459
const Evaluation & volumeFlux(unsigned phaseIdx) const
Return the volume flux of a fluid phase at the face's integration point .
Definition darcyfluxmodule.hh:168
void calculateGradients_(const ElementContext &elemCtx, unsigned faceIdx, unsigned timeIdx)
Calculate the gradients which are required to determine the volumetric fluxes.
Definition darcyfluxmodule.hh:183
const DimMatrix & intrinsicPermability() const
Returns the intrinsic permeability tensor for a given sub-control volume face.
Definition darcyfluxmodule.hh:138
Provides the intensive quantities for the Darcy flux module.
Definition darcyfluxmodule.hh:88
Callback class for a phase pressure.
Definition quantitycallbacks.hh:84
Defines the common parameters for the porous medium multi-phase models.
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
This method contains all callback classes for quantities that are required by some extensive quantiti...
Specifies a flux module which uses the Darcy relation.
Definition darcyfluxmodule.hh:61
static void registerParameters()
Register all run-time parameters for the flux module.
Definition darcyfluxmodule.hh:69