discretefractureextensivequantities.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
28#ifndef EWOMS_DISCRETE_FRACTURE_EXTENSIVE_QUANTITIES_HH
29#define EWOMS_DISCRETE_FRACTURE_EXTENSIVE_QUANTITIES_HH
30
32
33#include <dune/common/fvector.hh>
34#include <dune/common/fmatrix.hh>
35
36namespace Opm {
37
44template <class TypeTag>
46{
48
53
54 enum { dimWorld = GridView::dimensionworld };
55 enum { numPhases = FluidSystem::numPhases };
56
57 using DimMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
58 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
59
60public:
64 void update(const ElementContext& elemCtx, unsigned scvfIdx, unsigned timeIdx)
65 {
66 ParentType::update(elemCtx, scvfIdx, timeIdx);
67
68 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
69 const auto& stencil = elemCtx.stencil(timeIdx);
70 const auto& scvf = stencil.interiorFace(scvfIdx);
71 unsigned insideScvIdx = scvf.interiorIndex();
72 unsigned outsideScvIdx = scvf.exteriorIndex();
73
74 unsigned globalI = elemCtx.globalSpaceIndex(insideScvIdx, timeIdx);
75 unsigned globalJ = elemCtx.globalSpaceIndex(outsideScvIdx, timeIdx);
76 const auto& fractureMapper = elemCtx.problem().fractureMapper();
77 if (!fractureMapper.isFractureEdge(globalI, globalJ))
78 // do nothing if no fracture goes though the current edge
79 return;
80
81 // average the intrinsic permeability of the fracture
82 elemCtx.problem().fractureFaceIntrinsicPermeability(fractureIntrinsicPermeability_,
83 elemCtx, scvfIdx, timeIdx);
84
85 auto distDirection = elemCtx.pos(outsideScvIdx, timeIdx);
86 distDirection -= elemCtx.pos(insideScvIdx, timeIdx);
87 distDirection /= distDirection.two_norm();
88
89 const auto& problem = elemCtx.problem();
90 fractureWidth_ = problem.fractureWidth(elemCtx, insideScvIdx,
91 outsideScvIdx, timeIdx);
92 assert(fractureWidth_ < scvf.area());
93
94 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
95 const auto& pGrad = extQuants.potentialGrad(phaseIdx);
96
97 unsigned upstreamIdx = static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
98 const auto& up = elemCtx.intensiveQuantities(upstreamIdx, timeIdx);
99
100 // multiply with the fracture mobility of the upstream vertex
101 fractureIntrinsicPermeability_.mv(pGrad,
102 fractureFilterVelocity_[phaseIdx]);
103 fractureFilterVelocity_[phaseIdx] *= -up.fractureMobility(phaseIdx);
104
105 // divide the volume flux by two. This is required because
106 // a fracture is always shared by two sub-control-volume
107 // faces.
108 fractureVolumeFlux_[phaseIdx] = 0;
109 for (unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx)
110 fractureVolumeFlux_[phaseIdx] +=
111 (fractureFilterVelocity_[phaseIdx][dimIdx] * distDirection[dimIdx])
112 * (fractureWidth_ / 2.0) / scvf.area();
113 }
114 }
115
116public:
117 const DimMatrix& fractureIntrinsicPermeability() const
118 { return fractureIntrinsicPermeability_; }
119
120 Scalar fractureVolumeFlux(unsigned phaseIdx) const
121 { return fractureVolumeFlux_[phaseIdx]; }
122
123 Scalar fractureWidth() const
124 { return fractureWidth_; }
125
126 const DimVector& fractureFilterVelocity(unsigned phaseIdx) const
127 { return fractureFilterVelocity_[phaseIdx]; }
128
129private:
130 DimMatrix fractureIntrinsicPermeability_;
131 DimVector fractureFilterVelocity_[numPhases];
132 Scalar fractureVolumeFlux_[numPhases];
133 Scalar fractureWidth_;
134};
135
136} // namespace Opm
137
138#endif
This class expresses all intensive quantities of the discrete fracture model.
Definition: discretefractureextensivequantities.hh:46
Scalar fractureVolumeFlux(unsigned phaseIdx) const
Definition: discretefractureextensivequantities.hh:120
Scalar fractureWidth() const
Definition: discretefractureextensivequantities.hh:123
const DimVector & fractureFilterVelocity(unsigned phaseIdx) const
Definition: discretefractureextensivequantities.hh:126
const DimMatrix & fractureIntrinsicPermeability() const
Definition: discretefractureextensivequantities.hh:117
void update(const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Update the extensive quantities for a given sub-control-volume-face.
Definition: discretefractureextensivequantities.hh:64
This class provides the data all quantities that are required to calculate the fluxes of the fluid ph...
Definition: immiscibleextensivequantities.hh:52
void update(const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Update the extensive quantities for a given sub-control-volume-face.
Definition: immiscibleextensivequantities.hh:75
Definition: blackoilboundaryratevector.hh:37
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