discretefractureintensivequantities.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_INTENSIVE_QUANTITIES_HH
29 #define EWOMS_DISCRETE_FRACTURE_INTENSIVE_QUANTITIES_HH
30 
32 
34 
35 #include <opm/common/Valgrind.hpp>
36 
37 namespace Ewoms {
38 
47 template <class TypeTag>
49 {
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;
56 
57  enum { numPhases = FluidSystem::numPhases };
58  enum { dimWorld = GridView::dimensionworld };
59 
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!");
64 
65  enum { enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy) };
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  /*storeEnthalpy=*/enableEnergy> FluidState;
71 
72 public:
74  { }
75 
77 
79 
83  void update(const ElementContext& elemCtx, unsigned vertexIdx, unsigned timeIdx)
84  {
85  ParentType::update(elemCtx, vertexIdx, timeIdx);
86 
87  const auto& problem = elemCtx.problem();
88  const auto& fractureMapper = problem.fractureMapper();
89  unsigned globalVertexIdx = elemCtx.globalSpaceIndex(vertexIdx, timeIdx);
90 
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_);
96 
97  // do nothing if there is no fracture within the current degree of freedom
98  if (!fractureMapper.isFractureVertex(globalVertexIdx)) {
99  fractureVolume_ = 0;
100  return;
101  }
102 
103  // Make sure that the wetting saturation in the matrix fluid
104  // state does not get larger than 1
105  Scalar SwMatrix =
106  std::min<Scalar>(1.0, this->fluidState_.saturation(wettingPhaseIdx));
107  this->fluidState_.setSaturation(wettingPhaseIdx, SwMatrix);
108  this->fluidState_.setSaturation(nonWettingPhaseIdx, 1 - SwMatrix);
109 
110  // retrieve the facture width and intrinsic permeability from
111  // the problem
112  fracturePorosity_ =
113  problem.fracturePorosity(elemCtx, vertexIdx, timeIdx);
114  fractureIntrinsicPermeability_ =
115  problem.fractureIntrinsicPermeability(elemCtx, vertexIdx, timeIdx);
116 
117  // compute the fracture volume for the current sub-control
118  // volume. note, that we don't take overlaps of fractures into
119  // account for this.
120  fractureVolume_ = 0;
121  const auto& vertexPos = elemCtx.pos(vertexIdx, timeIdx);
122  for (unsigned vertex2Idx = 0; vertex2Idx < elemCtx.numDof(/*timeIdx=*/0); ++ vertex2Idx) {
123  unsigned globalVertex2Idx = elemCtx.globalSpaceIndex(vertex2Idx, timeIdx);
124 
125  if (vertexIdx == vertex2Idx ||
126  !fractureMapper.isFractureEdge(globalVertexIdx, globalVertex2Idx))
127  continue;
128 
129  Scalar fractureWidth =
130  problem.fractureWidth(elemCtx, vertexIdx, vertex2Idx, timeIdx);
131 
132  auto distVec = elemCtx.pos(vertex2Idx, timeIdx);
133  distVec -= vertexPos;
134 
135  Scalar edgeLength = distVec.two_norm();
136 
137  // the fracture is always adjacent to two sub-control
138  // volumes of the control volume, so when calculating the
139  // volume of the fracture which gets attributed to one
140  // SCV, the fracture width needs to divided by 2. Also,
141  // only half of the edge is located in the current control
142  // volume, so its length also needs to divided by 2.
143  fractureVolume_ += (fractureWidth / 2) * (edgeLength / 2);
144  }
145 
147  // set the fluid state for the fracture.
149 
150  // start with the same fluid state as in the matrix. This
151  // implies equal saturations, pressures, temperatures,
152  // enthalpies, etc.
153  fractureFluidState_.assign(this->fluidState_);
154 
155  // ask the problem for the material law parameters of the
156  // fracture.
157  const auto& fractureMatParams =
158  problem.fractureMaterialLawParams(elemCtx, vertexIdx, timeIdx);
159 
160  // calculate the fracture saturations which would be required
161  // to be consistent with the pressures
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]);
167 
168  // Make sure that the wetting saturation in the fracture does
169  // not get negative
170  Scalar SwFracture =
171  std::max<Scalar>(0.0, fractureFluidState_.saturation(wettingPhaseIdx));
172  fractureFluidState_.setSaturation(wettingPhaseIdx, SwFracture);
173  fractureFluidState_.setSaturation(nonWettingPhaseIdx, 1 - SwFracture);
174 
175  // calculate the relative permeabilities of the fracture
176  MaterialLaw::relativePermeabilities(fractureRelativePermeabilities_,
177  fractureMatParams,
178  fractureFluidState_);
179 
180  // make sure that valgrind complains if the fluid state is not
181  // fully defined.
182  fractureFluidState_.checkDefined();
183  }
184 
185 public:
192  Scalar fractureRelativePermeability(unsigned phaseIdx) const
193  { return fractureRelativePermeabilities_[phaseIdx]; }
194 
201  Scalar fractureMobility(unsigned phaseIdx) const
202  {
203  return fractureRelativePermeabilities_[phaseIdx]
204  / fractureFluidState_.viscosity(phaseIdx);
205  }
206 
210  Scalar fracturePorosity() const
211  { return fracturePorosity_; }
212 
217  const DimMatrix& fractureIntrinsicPermeability() const
218  { return fractureIntrinsicPermeability_; }
219 
224  Scalar fractureVolume() const
225  { return fractureVolume_; }
226 
231  const FluidState& fractureFluidState() const
232  { return fractureFluidState_; }
233 
234 protected:
235  FluidState fractureFluidState_;
236  Scalar fractureVolume_;
237  Scalar fracturePorosity_;
238  DimMatrix fractureIntrinsicPermeability_;
239  Scalar fractureRelativePermeabilities_[numPhases];
240 };
241 
242 } // namespace Ewoms
243 
244 #endif
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