30#ifndef EWOMS_FORCHHEIMER_FLUX_MODULE_HH
31#define EWOMS_FORCHHEIMER_FLUX_MODULE_HH
35#include <opm/common/Exceptions.hpp>
39#include <opm/material/common/Valgrind.hpp>
41#include <dune/common/fvector.hh>
42#include <dune/common/fmatrix.hh>
47template <
class TypeTag>
48class ForchheimerIntensiveQuantities;
50template <
class TypeTag>
51class ForchheimerExtensiveQuantities;
53template <
class TypeTag>
54class ForchheimerBaseProblem;
60template <
class TypeTag>
79template <
class TypeTag>
95 template <
class Context>
100 throw std::logic_error(
"Not implemented: Problem::ergunCoefficient()");
114 template <
class Context>
128template <
class TypeTag>
146 {
return ergunCoefficient_; }
152 {
return mobilityPassabilityRatio_[
phaseIdx]; }
157 const auto& problem =
elemCtx.problem();
161 mobilityPassabilityRatio_[
phaseIdx] =
162 problem.mobilityPassabilityRatio(
elemCtx,
169 Evaluation ergunCoefficient_;
170 Evaluation mobilityPassabilityRatio_[numPhases];
212template <
class TypeTag>
225 enum { dimWorld = GridView::dimensionworld };
230 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
231 using DimEvalVector = Dune::FieldVector<Evaluation, dimWorld>;
232 using DimMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
233 using DimEvalMatrix = Dune::FieldMatrix<Evaluation, dimWorld, dimWorld>;
240 {
return ergunCoefficient_; }
249 {
return mobilityPassabilityRatio_[
phaseIdx]; }
252 void calculateGradients_(
const ElementContext&
elemCtx,
259 unsigned i =
static_cast<unsigned>(this->interiorDofIdx_);
260 unsigned j =
static_cast<unsigned>(this->exteriorDofIdx_);
290 unsigned upIdx =
static_cast<unsigned>(this->upstreamIndex_(
phaseIdx));
296 mobilityPassabilityRatio_[
phaseIdx] =
302 mobilityPassabilityRatio_[
phaseIdx] =
308 template <
class Flu
idState>
309 void calculateBoundaryGradients_(
const ElementContext&
elemCtx,
312 const FluidState& fluidState)
320 unsigned i =
static_cast<unsigned>(this->interiorDofIdx_);
326 ergunCoefficient_ =
intQuantsIn.ergunCoefficient();
347 mobilityPassabilityRatio_[
phaseIdx] =
362 auto i = asImp_().interiorIndex();
363 auto j = asImp_().exteriorIndex();
368 const auto& normal =
scvf.normal();
369 Valgrind::CheckDefined(normal);
391 this->filterVelocity_[
phaseIdx] = 0.0;
396 calculateForchheimerFlux_(
phaseIdx);
413 const auto& normal = boundaryFace.normal();
420 this->filterVelocity_[
phaseIdx] = 0.0;
425 calculateForchheimerFlux_(
phaseIdx);
434 void calculateForchheimerFlux_(
unsigned phaseIdx)
437 DimEvalVector& velocity = this->filterVelocity_[
phaseIdx];
444 DimEvalVector residual;
450 while (
deltaV.one_norm() > 1
e-11) {
465 void forchheimerResid_(DimEvalVector& residual,
unsigned phaseIdx)
const
467 const DimEvalVector& velocity = this->filterVelocity_[
phaseIdx];
470 const auto& mobility = this->mobility_[
phaseIdx];
471 const auto& density = density_[
phaseIdx];
509 Valgrind::CheckDefined(residual);
512 void gradForchheimerResid_(DimEvalVector& residual,
517 DimEvalVector& velocity = this->filterVelocity_[
phaseIdx];
518 forchheimerResid_(residual,
phaseIdx);
522 for (
unsigned i = 0; i < dimWorld; ++i) {
523 Scalar
coordEps = std::max(
eps, Toolbox::scalarValue(velocity[i]) * (1 +
eps));
542 for (
unsigned i = 0; i < dimWorld; i++) {
543 for (
unsigned j = 0; j < dimWorld; j++) {
547 if (std::abs(
K[i][j]) > 1
e-25)
555 Implementation& asImp_()
556 {
return *
static_cast<Implementation *
>(
this); }
558 const Implementation& asImp_()
const
559 {
return *
static_cast<const Implementation *
>(
this); }
566 Evaluation ergunCoefficient_;
569 Evaluation mobilityPassabilityRatio_[numPhases];
572 Evaluation density_[numPhases];
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
void calculateGradients_(const ElementContext &elemCtx, unsigned faceIdx, unsigned timeIdx)
Calculate the gradients which are required to determine the volumetric fluxes.
Definition darcyfluxmodule.hh:183
Provides the defaults for the parameters required by the Forchheimer velocity approach.
Definition forchheimerfluxmodule.hh:81
Evaluation mobilityPassabilityRatio(Context &context, unsigned spaceIdx, unsigned timeIdx, unsigned phaseIdx) const
Returns the ratio between the phase mobility and its passability for a given fluid phase .
Definition forchheimerfluxmodule.hh:115
Scalar ergunCoefficient(const Context &, unsigned, unsigned) const
Returns the Ergun coefficient.
Definition forchheimerfluxmodule.hh:96
Provides the Forchheimer flux module.
Definition forchheimerfluxmodule.hh:215
const Evaluation & ergunCoefficient() const
Return the Ergun coefficent at the face's integration point.
Definition forchheimerfluxmodule.hh:239
void calculateBoundaryFluxes_(const ElementContext &elemCtx, unsigned bfIdx, unsigned timeIdx)
Calculate the volumetric flux rates of all phases at the domain boundary.
Definition forchheimerfluxmodule.hh:408
void calculateFluxes_(const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Calculate the volumetric fluxes of all phases.
Definition forchheimerfluxmodule.hh:359
bool isDiagonal_(const DimMatrix &K) const
Check whether all off-diagonal entries of a tensor are zero.
Definition forchheimerfluxmodule.hh:540
Evaluation & mobilityPassabilityRatio(unsigned phaseIdx) const
Return the ratio of the mobility divided by the passability at the face's integration point for a giv...
Definition forchheimerfluxmodule.hh:248
Provides the intensive quantities for the Forchheimer module.
Definition forchheimerfluxmodule.hh:130
const Evaluation & ergunCoefficient() const
Returns the Ergun coefficient.
Definition forchheimerfluxmodule.hh:145
const Evaluation & mobilityPassabilityRatio(unsigned phaseIdx) const
Returns the passability of a phase.
Definition forchheimerfluxmodule.hh:151
This file contains the necessary classes to calculate the volumetric fluxes out of a pressure potenti...
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
Specifies a flux module which uses the Forchheimer relation.
Definition forchheimerfluxmodule.hh:62
static void registerParameters()
Register all run-time parameters for the flux module.
Definition forchheimerfluxmodule.hh:70