26 #ifndef EWOMS_DISCRETE_FRACTURE_EXTENSIVE_QUANTITIES_HH
27 #define EWOMS_DISCRETE_FRACTURE_EXTENSIVE_QUANTITIES_HH
31 #include <dune/common/fvector.hh>
32 #include <dune/common/fmatrix.hh>
42 template <
class TypeTag>
47 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
48 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
49 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
50 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
52 enum { dimWorld = GridView::dimensionworld };
53 enum { numPhases = FluidSystem::numPhases };
55 typedef Dune::FieldMatrix<Scalar, dimWorld, dimWorld> DimMatrix;
56 typedef Dune::FieldVector<Scalar, dimWorld> DimVector;
62 void update(
const ElementContext &elemCtx,
int scvfIdx,
int timeIdx)
66 const auto &extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
67 const auto &stencil = elemCtx.stencil(timeIdx);
68 const auto &scvf = stencil.interiorFace(scvfIdx);
69 int insideScvIdx = scvf.interiorIndex();
70 int outsideScvIdx = scvf.exteriorIndex();
72 int globalI = elemCtx.globalSpaceIndex(insideScvIdx, timeIdx);
73 int globalJ = elemCtx.globalSpaceIndex(outsideScvIdx, timeIdx);
74 const auto &fractureMapper = elemCtx.problem().fractureMapper();
75 if (!fractureMapper.isFractureEdge(globalI, globalJ))
80 elemCtx.problem().fractureFaceIntrinsicPermeability(fractureIntrinsicPermeability_,
81 elemCtx, scvfIdx, timeIdx);
83 auto distDirection = elemCtx.pos(outsideScvIdx, timeIdx);
84 distDirection -= elemCtx.pos(insideScvIdx, timeIdx);
85 distDirection /= distDirection.two_norm();
87 const auto &problem = elemCtx.problem();
88 fractureWidth_ = problem.fractureWidth(elemCtx, insideScvIdx,
89 outsideScvIdx, timeIdx);
90 assert(fractureWidth_ < scvf.area());
92 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
93 const auto &pGrad = extQuants.potentialGrad(phaseIdx);
95 int upstreamIdx = extQuants.upstreamIndex(phaseIdx);
96 const auto &up = elemCtx.intensiveQuantities(upstreamIdx, timeIdx);
99 fractureIntrinsicPermeability_.mv(pGrad,
100 fractureFilterVelocity_[phaseIdx]);
101 fractureFilterVelocity_[phaseIdx] *= -up.fractureMobility(phaseIdx);
106 fractureVolumeFlux_[phaseIdx] = 0;
107 for (
int dimIdx = 0; dimIdx < dimWorld; ++dimIdx)
108 fractureVolumeFlux_[phaseIdx] +=
109 (fractureFilterVelocity_[phaseIdx][dimIdx] * distDirection[dimIdx])
110 * (fractureWidth_ / 2.0) / scvf.area();
116 {
return fractureIntrinsicPermeability_; }
119 {
return fractureVolumeFlux_[phaseIdx]; }
122 {
return fractureWidth_; }
125 {
return fractureFilterVelocity_[phaseIdx]; }
128 DimMatrix fractureIntrinsicPermeability_;
129 DimVector fractureFilterVelocity_[numPhases];
130 Scalar fractureVolumeFlux_[numPhases];
131 Scalar fractureWidth_;
This class provides the data all quantities that are required to calculate the fluxes of the fluid ph...
void update(const ElementContext &elemCtx, int scvfIdx, int timeIdx)
Update the extensive quantities for a given sub-control-volume-face.
Definition: immiscibleextensivequantities.hh:71
const DimMatrix & fractureIntrinsicPermeability() const
Definition: discretefractureextensivequantities.hh:115
const DimVector & fractureFilterVelocity(int phaseIdx) const
Definition: discretefractureextensivequantities.hh:124
This class expresses all intensive quantities of the discrete fracture model.
Definition: discretefractureextensivequantities.hh:43
This class provides the data all quantities that are required to calculate the fluxes of the fluid ph...
Definition: immiscibleextensivequantities.hh:47
Definition: baseauxiliarymodule.hh:35
Scalar fractureVolumeFlux(int phaseIdx) const
Definition: discretefractureextensivequantities.hh:118
void update(const ElementContext &elemCtx, int scvfIdx, int timeIdx)
Update the extensive quantities for a given sub-control-volume-face.
Definition: discretefractureextensivequantities.hh:62
Scalar fractureWidth() const
Definition: discretefractureextensivequantities.hh:121