26 #ifndef EWOMS_DISCRETE_FRACTURE_INTENSIVE_QUANTITIES_HH
27 #define EWOMS_DISCRETE_FRACTURE_INTENSIVE_QUANTITIES_HH
43 template <
class TypeTag>
47 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
48 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
49 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
50 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
51 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
53 enum { numPhases = FluidSystem::numPhases };
54 enum { dimWorld = GridView::dimensionworld };
56 static_assert(dimWorld == 2,
"The fracture module currently is only "
57 "implemented for the 2D case!");
58 static_assert(numPhases == 2,
"The fracture module currently is only "
59 "implemented for two fluid phases!");
62 enum { wettingPhaseIdx = MaterialLaw::wettingPhaseIdx };
63 enum { nonWettingPhaseIdx = MaterialLaw::nonWettingPhaseIdx };
64 typedef Dune::FieldMatrix<Scalar, dimWorld, dimWorld> DimMatrix;
65 typedef Opm::ImmiscibleFluidState<Scalar, FluidSystem,
66 enableEnergy> FluidState;
72 void update(
const ElementContext &elemCtx,
int vertexIdx,
int timeIdx)
76 const auto &problem = elemCtx.problem();
77 const auto &fractureMapper = problem.fractureMapper();
78 int globalVertexIdx = elemCtx.globalSpaceIndex(vertexIdx, timeIdx);
80 Valgrind::SetUndefined(fractureFluidState_);
81 Valgrind::SetUndefined(fractureVolume_);
82 Valgrind::SetUndefined(fracturePorosity_);
83 Valgrind::SetUndefined(fractureIntrinsicPermeability_);
84 Valgrind::SetUndefined(fractureRelativePermeabilities_);
87 if (!fractureMapper.isFractureVertex(globalVertexIdx)) {
94 Scalar SwMatrix = std::min(1.0, this->
fluidState_.saturation(wettingPhaseIdx));
95 this->
fluidState_.setSaturation(wettingPhaseIdx, SwMatrix);
96 this->
fluidState_.setSaturation(nonWettingPhaseIdx, 1 - SwMatrix);
101 problem.fracturePorosity(elemCtx, vertexIdx, timeIdx);
102 fractureIntrinsicPermeability_ =
103 problem.fractureIntrinsicPermeability(elemCtx, vertexIdx, timeIdx);
109 const auto &vertexPos = elemCtx.pos(vertexIdx, timeIdx);
110 for (
int vertex2Idx = 0; vertex2Idx < elemCtx.numDof(0); ++ vertex2Idx) {
111 int globalVertex2Idx = elemCtx.globalSpaceIndex(vertex2Idx, timeIdx);
113 if (vertexIdx == vertex2Idx ||
114 !fractureMapper.isFractureEdge(globalVertexIdx, globalVertex2Idx))
117 Scalar fractureWidth =
118 problem.fractureWidth(elemCtx, vertexIdx, vertex2Idx, timeIdx);
120 auto distVec = elemCtx.pos(vertex2Idx, timeIdx);
121 distVec -= vertexPos;
123 Scalar edgeLength = distVec.two_norm();
131 fractureVolume_ += (fractureWidth / 2) * (edgeLength / 2);
134 if (fractureVolume_ <= 0.0)
148 const auto &fractureMatParams =
149 problem.fractureMaterialLawParams(elemCtx, vertexIdx, timeIdx);
153 Scalar saturations[numPhases];
154 MaterialLaw::saturations(saturations, fractureMatParams,
155 fractureFluidState_);
156 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
157 fractureFluidState_.setSaturation(phaseIdx, saturations[phaseIdx]);
161 Scalar SwFracture = std::max(0.0, fractureFluidState_.saturation(wettingPhaseIdx));
162 fractureFluidState_.setSaturation(wettingPhaseIdx, SwFracture);
163 fractureFluidState_.setSaturation(nonWettingPhaseIdx, 1 - SwFracture);
166 MaterialLaw::relativePermeabilities(fractureRelativePermeabilities_,
168 fractureFluidState_);
172 fractureFluidState_.checkDefined();
182 Scalar fractureRelativePermeability(
int phaseIdx)
const
183 {
return fractureRelativePermeabilities_[phaseIdx]; }
191 Scalar fractureMobility(
int phaseIdx)
const
193 return fractureRelativePermeabilities_[phaseIdx]
194 / fractureFluidState_.viscosity(phaseIdx);
200 Scalar fracturePorosity()
const
201 {
return fracturePorosity_; }
207 const DimMatrix &fractureIntrinsicPermeability()
const
208 {
return fractureIntrinsicPermeability_; }
214 Scalar fractureVolume()
const
215 {
return fractureVolume_; }
221 const FluidState &fractureFluidState()
const
222 {
return fractureFluidState_; }
225 FluidState fractureFluidState_;
226 Scalar fractureVolume_;
227 Scalar fracturePorosity_;
228 DimMatrix fractureIntrinsicPermeability_;
229 Scalar fractureRelativePermeabilities_[numPhases];
void update(const ElementContext &elemCtx, int dofIdx, int timeIdx)
Definition: immiscibleintensivequantities.hh:82
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:468
Contains the quantities which are are constant within a finite volume in the discret fracture immisci...
Definition: discretefractureintensivequantities.hh:44
Definition: baseauxiliarymodule.hh:35
FluidState fluidState_
Definition: immiscibleintensivequantities.hh:181
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 ...
Definition: immiscibleintensivequantities.hh:46
Contains the quantities which are are constant within a finite volume for the immiscible multi-phase ...