All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
blackoillocalresidual.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_BLACK_OIL_LOCAL_RESIDUAL_HH
29 #define EWOMS_BLACK_OIL_LOCAL_RESIDUAL_HH
30 
31 #include "blackoilproperties.hh"
34 
35 
36 namespace Ewoms {
42 template <class TypeTag>
43 class BlackOilLocalResidual : public GET_PROP_TYPE(TypeTag, DiscLocalResidual)
44 {
45  typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
46  typedef typename GET_PROP_TYPE(TypeTag, ExtensiveQuantities) ExtensiveQuantities;
47  typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
48  typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
49  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
50  typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
51  typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
52  typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
53  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
54 
55  enum { conti0EqIdx = Indices::conti0EqIdx };
56  enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
57  enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
58  enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
59 
60  enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
61  enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
62  enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
63 
64  enum { gasCompIdx = FluidSystem::gasCompIdx };
65  enum { oilCompIdx = FluidSystem::oilCompIdx };
66  enum { waterCompIdx = FluidSystem::waterCompIdx };
67  enum { compositionSwitchIdx = Indices::compositionSwitchIdx };
68 
69  static const bool compositionSwitchEnabled = Indices::gasEnabled;
70  static const bool waterEnabled = Indices::waterEnabled;
71 
72  static constexpr bool blackoilConserveSurfaceVolume = GET_PROP_VALUE(TypeTag, BlackoilConserveSurfaceVolume);
73 
74  typedef Opm::MathToolbox<Evaluation> Toolbox;
77 
78 
79 public:
83  template <class LhsEval>
84  void computeStorage(Dune::FieldVector<LhsEval, numEq>& storage,
85  const ElementContext& elemCtx,
86  unsigned dofIdx,
87  unsigned timeIdx) const
88  {
89  // retrieve the intensive quantities for the SCV at the specified point in time
90  const IntensiveQuantities& intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
91  const auto& fs = intQuants.fluidState();
92 
93  storage = 0.0;
94 
95  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
96  if (!FluidSystem::phaseIsActive(phaseIdx))
97  continue;
98 
99  unsigned activeCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
100  LhsEval surfaceVolume =
101  Toolbox::template decay<LhsEval>(fs.saturation(phaseIdx))
102  * Toolbox::template decay<LhsEval>(fs.invB(phaseIdx))
103  * Toolbox::template decay<LhsEval>(intQuants.porosity());
104 
105  storage[conti0EqIdx + activeCompIdx] += surfaceVolume;
106 
107  // account for dissolved gas
108  if (phaseIdx == oilPhaseIdx && FluidSystem::enableDissolvedGas()) {
109  unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(gasCompIdx);
110  storage[conti0EqIdx + activeGasCompIdx] +=
111  Toolbox::template decay<LhsEval>(intQuants.fluidState().Rs())
112  * surfaceVolume;
113  }
114 
115  // account for vaporized oil
116  if (phaseIdx == gasPhaseIdx && FluidSystem::enableVaporizedOil()) {
117  unsigned activeOilCompIdx = Indices::canonicalToActiveComponentIndex(oilCompIdx);
118  storage[conti0EqIdx + activeOilCompIdx] +=
119  Toolbox::template decay<LhsEval>(intQuants.fluidState().Rv())
120  * surfaceVolume;
121  }
122  }
123 
124  // convert surface volumes to component masses
125  if (!blackoilConserveSurfaceVolume) {
126  unsigned pvtRegionIdx = intQuants.pvtRegionIndex();
127  if (waterEnabled) {
128  unsigned activeWaterCompIdx = Indices::canonicalToActiveComponentIndex(waterCompIdx);
129  storage[conti0EqIdx + activeWaterCompIdx] *=
130  FluidSystem::referenceDensity(waterPhaseIdx, pvtRegionIdx);
131  }
132  if (compositionSwitchEnabled) {
133  unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(gasCompIdx);
134  storage[conti0EqIdx + activeGasCompIdx] *=
135  FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx);
136  }
137  unsigned activeOilCompIdx = Indices::canonicalToActiveComponentIndex(oilCompIdx);
138  storage[conti0EqIdx + activeOilCompIdx] *=
139  FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx);
140  }
141  // deal with solvents (if present)
142  SolventModule::addStorage(storage, intQuants);
143 
144  // deal with polymer (if present)
145  PolymerModule::addStorage(storage, intQuants);
146 
147  }
148 
152  void computeFlux(RateVector& flux,
153  const ElementContext& elemCtx,
154  unsigned scvfIdx,
155  unsigned timeIdx) const
156  {
157  assert(timeIdx == 0);
158 
159  flux = 0.0;
160 
161  const ExtensiveQuantities& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
162  unsigned focusDofIdx = elemCtx.focusDofIndex();
163  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
164  if (!FluidSystem::phaseIsActive(phaseIdx))
165  continue;
166 
167  unsigned upIdx = static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
168  const IntensiveQuantities& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
169  if (upIdx == focusDofIdx)
170  evalPhaseFluxes_<Evaluation>(flux, phaseIdx, extQuants, up);
171  else
172  evalPhaseFluxes_<Scalar>(flux, phaseIdx, extQuants, up);
173  }
174 
175  if (!blackoilConserveSurfaceVolume) {
176  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
177  unsigned activeCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
178  unsigned upIdx = static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
179  const IntensiveQuantities& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
180  unsigned pvtRegionIdx = up.pvtRegionIndex();
181  Scalar refDensity = FluidSystem::referenceDensity(phaseIdx, pvtRegionIdx);
182  flux[conti0EqIdx + activeCompIdx] *= refDensity;
183  }
184  }
185 
186  // deal with solvents (if present)
187  SolventModule::computeFlux(flux, elemCtx, scvfIdx, timeIdx);
188 
189  // deal with polymer (if present)
190  PolymerModule::computeFlux(flux, elemCtx, scvfIdx, timeIdx);
191  }
192 
196  void computeSource(RateVector& source,
197  const ElementContext& elemCtx,
198  unsigned dofIdx,
199  unsigned timeIdx) const
200  {
201  // retrieve the source term intrinsic to the problem
202  elemCtx.problem().source(source, elemCtx, dofIdx, timeIdx);
203  }
204 
205 protected:
206  template <class UpEval>
207  void evalPhaseFluxes_(RateVector& flux,
208  unsigned phaseIdx,
209  const ExtensiveQuantities& extQuants,
210  const IntensiveQuantities& up) const
211  {
212  unsigned compIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
213  const auto& fs = up.fluidState();
214 
215  Evaluation surfaceVolumeFlux =
216  Toolbox::template decay<UpEval>(fs.invB(phaseIdx))
217  * extQuants.volumeFlux(phaseIdx);
218 
219  flux[conti0EqIdx + compIdx] +=
220  surfaceVolumeFlux;
221 
222  // dissolved gas (in the oil phase).
223  if (phaseIdx == oilPhaseIdx && FluidSystem::enableDissolvedGas()) {
224  flux[conti0EqIdx + Indices::canonicalToActiveComponentIndex(gasCompIdx)] +=
225  Toolbox::template decay<UpEval>(fs.Rs())
226  * surfaceVolumeFlux;
227 
228  }
229 
230  // vaporized oil (in the gas phase).
231  if (phaseIdx == gasPhaseIdx && FluidSystem::enableVaporizedOil()) {
232  flux[conti0EqIdx + Indices::canonicalToActiveComponentIndex(oilCompIdx)] +=
233  Toolbox::template decay<UpEval>(fs.Rv())
234  * surfaceVolumeFlux;
235 
236  }
237  }
238 
239 };
240 
241 } // namespace Ewoms
242 
243 #endif
void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx) const
Evaluates the total mass flux of all conservation quantities over a face of a sub-control volume...
Definition: blackoillocalresidual.hh:152
void computeSource(RateVector &source, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Calculate the source term of the equation.
Definition: blackoillocalresidual.hh:196
void computeStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Evaluate the amount all conservation quantities (e.g.
Definition: blackoillocalresidual.hh:84
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
Contains the classes required to extend the black-oil model by solvents.
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
Declares the properties required by the black oil model.
Contains the high level supplements required to extend the black oil model by solvents.
Definition: blackoilprimaryvariables.hh:47
Calculates the local residual of the black oil model.
Definition: blackoillocalresidual.hh:43
Contains the high level supplements required to extend the black oil model by polymer.
Definition: blackoilpolymermodules.hh:64
Contains the classes required to extend the black-oil model by polymer.