26 #ifndef EWOMS_DISCRETE_FRACTURE_LOCAL_RESIDUAL_BASE_HH
27 #define EWOMS_DISCRETE_FRACTURE_LOCAL_RESIDUAL_BASE_HH
39 template <
class TypeTag>
44 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
45 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
46 typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
47 typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
48 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
50 enum { conti0EqIdx = Indices::conti0EqIdx };
66 const ElementContext &elemCtx,
71 EqVector phaseStorage(0.0);
74 const auto &problem = elemCtx.problem();
75 const auto &fractureMapper = problem.fractureMapper();
76 int globalIdx = elemCtx.globalSpaceIndex(dofIdx, timeIdx);
78 if (!fractureMapper.isFractureVertex(globalIdx)) {
81 storage += phaseStorage;
85 const auto &intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
86 const auto &scv = elemCtx.stencil(timeIdx).subControlVolume(dofIdx);
89 phaseStorage *= 1 - intQuants.fractureVolume()/scv.volume();
92 const auto &fsFracture = intQuants.fractureFluidState();
94 phaseStorage[conti0EqIdx + phaseIdx] +=
95 intQuants.fracturePorosity()*
96 fsFracture.saturation(phaseIdx) *
97 fsFracture.density(phaseIdx) *
98 intQuants.fractureVolume()/scv.volume();
100 EnergyModule::addFracturePhaseStorage(phaseStorage, intQuants, scv,
104 storage += phaseStorage;
111 int scvfIdx,
int timeIdx)
const
115 const auto &extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
117 int i = extQuants.interiorIndex();
118 int j = extQuants.exteriorIndex();
119 int I = elemCtx.globalSpaceIndex(i, timeIdx);
120 int J = elemCtx.globalSpaceIndex(j, timeIdx);
121 const auto &fractureMapper = elemCtx.problem().fractureMapper();
122 if (!fractureMapper.isFractureEdge(I, J))
127 const auto &scvf = elemCtx.stencil(timeIdx).interiorFace(scvfIdx);
128 Scalar scvfArea = scvf.area();
131 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
132 if (!elemCtx.model().phaseIsConsidered(phaseIdx))
139 flux[conti0EqIdx + phaseIdx] *=
140 1 - extQuants.fractureWidth() / (2 * scvfArea);
143 int upIdx = extQuants.upstreamIndex(phaseIdx);
144 const auto &up = elemCtx.intensiveQuantities(upIdx, timeIdx);
145 flux[conti0EqIdx + phaseIdx] +=
146 extQuants.fractureVolumeFlux(phaseIdx) * up.fractureFluidState().density(phaseIdx);
149 EnergyModule::handleFractureFlux(flux, elemCtx, scvfIdx, timeIdx);
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:468
Calculates the local residual of the immiscible multi-phase model.
Definition: immisciblelocalresidual.hh:41
void computeFlux(RateVector &flux, const ElementContext &elemCtx, int scvfIdx, int timeIdx) const
Evaluates the total mass flux of all conservation quantities over a face of a sub-control volume...
Definition: immisciblelocalresidual.hh:111
Provides the auxiliary methods required for consideration of the energy equation. ...
Definition: energymodule.hh:54
Definition: baseauxiliarymodule.hh:35
void addPhaseStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, int dofIdx, int timeIdx, int phaseIdx) const
Adds the amount all conservation quantities (e.g. phase mass) within a single fluid phase...
Definition: immisciblelocalresidual.hh:71
Calculates the local residual of the discrete fracture immiscible multi-phase model.
Definition: discretefracturelocalresidual.hh:40
void computeFlux(RateVector &flux, const ElementContext &elemCtx, int scvfIdx, int timeIdx) const
Evaluates the total mass flux of all conservation quantities over a face of a sub-control volume...
Definition: discretefracturelocalresidual.hh:110
Calculates the local residual of the immiscible multi-phase model.
void addPhaseStorage(EqVector &storage, const ElementContext &elemCtx, int dofIdx, int timeIdx, int phaseIdx) const
Adds the amount all conservation quantities (e.g. phase mass) within a single fluid phase...
Definition: discretefracturelocalresidual.hh:65