28 #ifndef EWOMS_DISCRETE_FRACTURE_INTENSIVE_QUANTITIES_HH 29 #define EWOMS_DISCRETE_FRACTURE_INTENSIVE_QUANTITIES_HH 35 #include <opm/common/Valgrind.hpp> 47 template <
class TypeTag>
51 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
52 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
53 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
54 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
55 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
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 typedef Dune::FieldMatrix<Scalar, dimWorld, dimWorld> DimMatrix;
69 typedef Opm::ImmiscibleFluidState<Scalar, FluidSystem,
70 enableEnergy> FluidState;
83 void update(
const ElementContext& elemCtx,
unsigned vertexIdx,
unsigned timeIdx)
87 const auto& problem = elemCtx.problem();
88 const auto& fractureMapper = problem.fractureMapper();
89 unsigned globalVertexIdx = elemCtx.globalSpaceIndex(vertexIdx, timeIdx);
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_);
98 if (!fractureMapper.isFractureVertex(globalVertexIdx)) {
106 std::min<Scalar>(1.0, this->fluidState_.saturation(wettingPhaseIdx));
107 this->fluidState_.setSaturation(wettingPhaseIdx, SwMatrix);
108 this->fluidState_.setSaturation(nonWettingPhaseIdx, 1 - SwMatrix);
113 problem.fracturePorosity(elemCtx, vertexIdx, timeIdx);
114 fractureIntrinsicPermeability_ =
115 problem.fractureIntrinsicPermeability(elemCtx, vertexIdx, timeIdx);
121 const auto& vertexPos = elemCtx.pos(vertexIdx, timeIdx);
122 for (
unsigned vertex2Idx = 0; vertex2Idx < elemCtx.numDof(0); ++ vertex2Idx) {
123 unsigned globalVertex2Idx = elemCtx.globalSpaceIndex(vertex2Idx, timeIdx);
125 if (vertexIdx == vertex2Idx ||
126 !fractureMapper.isFractureEdge(globalVertexIdx, globalVertex2Idx))
129 Scalar fractureWidth =
130 problem.fractureWidth(elemCtx, vertexIdx, vertex2Idx, timeIdx);
132 auto distVec = elemCtx.pos(vertex2Idx, timeIdx);
133 distVec -= vertexPos;
135 Scalar edgeLength = distVec.two_norm();
143 fractureVolume_ += (fractureWidth / 2) * (edgeLength / 2);
153 fractureFluidState_.assign(this->fluidState_);
157 const auto& fractureMatParams =
158 problem.fractureMaterialLawParams(elemCtx, vertexIdx, timeIdx);
162 Scalar saturations[numPhases];
163 MaterialLaw::saturations(saturations, fractureMatParams,
164 fractureFluidState_);
165 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
166 fractureFluidState_.setSaturation(phaseIdx, saturations[phaseIdx]);
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];
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 ...
void update(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Definition: immiscibleintensivequantities.hh:93
Scalar fractureMobility(unsigned phaseIdx) const
Returns the effective mobility of a given phase within the control volume.
Definition: discretefractureintensivequantities.hh:201
Contains the quantities which are are constant within a finite volume for the immiscible multi-phase ...
Definition: immiscibleintensivequantities.hh:50
Definition: baseauxiliarymodule.hh:37
Contains the quantities which are are constant within a finite volume in the discret fracture immisci...
Definition: discretefractureintensivequantities.hh:48
Scalar fracturePorosity() const
Returns the average porosity within the fracture.
Definition: discretefractureintensivequantities.hh:210
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
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
void update(const ElementContext &elemCtx, unsigned vertexIdx, unsigned timeIdx)
Definition: discretefractureintensivequantities.hh:83
const FluidState & fractureFluidState() const
Returns a fluid state object which represents the thermodynamic state of the fluids within the fractu...
Definition: discretefractureintensivequantities.hh:231
const DimMatrix & fractureIntrinsicPermeability() const
Returns the average intrinsic permeability within the fracture.
Definition: discretefractureintensivequantities.hh:217