28#ifndef EWOMS_DISCRETE_FRACTURE_INTENSIVE_QUANTITIES_HH
29#define EWOMS_DISCRETE_FRACTURE_INTENSIVE_QUANTITIES_HH
35#include <opm/material/common/Valgrind.hpp>
47template <
class TypeTag>
57 enum { numPhases = FluidSystem::numPhases };
58 enum { dimWorld = GridView::dimensionworld };
60 static_assert(dimWorld == 2,
"The fracture module currently is only "
61 "implemented for the 2D case!");
62 static_assert(numPhases == 2,
"The fracture module currently is only "
63 "implemented for two fluid phases!");
66 enum { wettingPhaseIdx = MaterialLaw::wettingPhaseIdx };
67 enum { nonWettingPhaseIdx = MaterialLaw::nonWettingPhaseIdx };
68 using DimMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
87 const auto& problem =
elemCtx.problem();
88 const auto& fractureMapper = problem.fractureMapper();
91 Opm::Valgrind::SetUndefined(fractureFluidState_);
92 Opm::Valgrind::SetUndefined(fractureVolume_);
93 Opm::Valgrind::SetUndefined(fracturePorosity_);
94 Opm::Valgrind::SetUndefined(fractureIntrinsicPermeability_);
95 Opm::Valgrind::SetUndefined(fractureRelativePermeabilities_);
106 std::min<Scalar>(1.0, this->fluidState_.saturation(wettingPhaseIdx));
107 this->fluidState_.setSaturation(wettingPhaseIdx,
SwMatrix);
108 this->fluidState_.setSaturation(nonWettingPhaseIdx, 1 -
SwMatrix);
114 fractureIntrinsicPermeability_ =
129 Scalar fractureWidth =
143 fractureVolume_ += (fractureWidth / 2) * (
edgeLength / 2);
153 fractureFluidState_.assign(this->fluidState_);
164 fractureFluidState_);
171 std::max<Scalar>(0.0, fractureFluidState_.saturation(wettingPhaseIdx));
172 fractureFluidState_.setSaturation(wettingPhaseIdx,
SwFracture);
173 fractureFluidState_.setSaturation(nonWettingPhaseIdx, 1 -
SwFracture);
176 MaterialLaw::relativePermeabilities(fractureRelativePermeabilities_,
178 fractureFluidState_);
182 fractureFluidState_.checkDefined();
193 {
return fractureRelativePermeabilities_[
phaseIdx]; }
203 return fractureRelativePermeabilities_[
phaseIdx]
204 / fractureFluidState_.viscosity(
phaseIdx);
211 {
return fracturePorosity_; }
218 {
return fractureIntrinsicPermeability_; }
225 {
return fractureVolume_; }
232 {
return fractureFluidState_; }
235 FluidState fractureFluidState_;
236 Scalar fractureVolume_;
237 Scalar fracturePorosity_;
238 DimMatrix fractureIntrinsicPermeability_;
239 Scalar fractureRelativePermeabilities_[numPhases];
Contains the quantities which are are constant within a finite volume in the discret fracture immisci...
Definition discretefractureintensivequantities.hh:49
Scalar fracturePorosity() const
Returns the average porosity within the fracture.
Definition discretefractureintensivequantities.hh:210
const FluidState & fractureFluidState() const
Returns a fluid state object which represents the thermodynamic state of the fluids within the fractu...
Definition discretefractureintensivequantities.hh:231
Scalar fractureMobility(unsigned phaseIdx) const
Returns the effective mobility of a given phase within the control volume.
Definition discretefractureintensivequantities.hh:201
Scalar fractureRelativePermeability(unsigned phaseIdx) const
Returns the effective mobility of a given phase within the control volume.
Definition discretefractureintensivequantities.hh:192
Scalar fractureVolume() const
Return the volume [m^2] occupied by fractures within the given sub-control volume.
Definition discretefractureintensivequantities.hh:224
const DimMatrix & fractureIntrinsicPermeability() const
Returns the average intrinsic permeability within the fracture.
Definition discretefractureintensivequantities.hh:217
void update(const ElementContext &elemCtx, unsigned vertexIdx, unsigned timeIdx)
Definition discretefractureintensivequantities.hh:83
Contains the quantities which are are constant within a finite volume for the immiscible multi-phase ...
Definition immiscibleintensivequantities.hh:54
void update(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Definition immiscibleintensivequantities.hh:93
Defines the properties required for the immiscible multi-phase model which considers discrete fractur...
Contains the quantities which are are constant within a finite volume for the immiscible multi-phase ...
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