All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
pvslocalresidual.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_PVS_LOCAL_RESIDUAL_HH
29 #define EWOMS_PVS_LOCAL_RESIDUAL_HH
30 
31 #include "pvsproperties.hh"
32 
35 
36 #include <opm/common/Valgrind.hpp>
37 
38 namespace Ewoms {
39 
46 template <class TypeTag>
47 class PvsLocalResidual : public GET_PROP_TYPE(TypeTag, DiscLocalResidual)
48 {
49  typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
50  typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
51  typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
52  typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
53  typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
54  typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
55 
56  enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
57  enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
58  enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
59  enum { conti0EqIdx = Indices::conti0EqIdx };
60 
61  enum { enableDiffusion = GET_PROP_VALUE(TypeTag, EnableDiffusion) };
63 
64  enum { enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy) };
66 
67  typedef Opm::MathToolbox<Evaluation> Toolbox;
68 
69 public:
73  template <class LhsEval>
74  void addPhaseStorage(Dune::FieldVector<LhsEval, numEq>& storage,
75  const ElementContext& elemCtx,
76  unsigned dofIdx,
77  unsigned timeIdx,
78  unsigned phaseIdx) const
79  {
80  const IntensiveQuantities& intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
81  const auto& fs = intQuants.fluidState();
82 
83  // compute storage term of all components within all phases
84  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
85  unsigned eqIdx = conti0EqIdx + compIdx;
86  storage[eqIdx] +=
87  Toolbox::template decay<LhsEval>(fs.molarity(phaseIdx, compIdx))
88  * Toolbox::template decay<LhsEval>(fs.saturation(phaseIdx))
89  * Toolbox::template decay<LhsEval>(intQuants.porosity());
90  }
91 
92  EnergyModule::addPhaseStorage(storage, elemCtx.intensiveQuantities(dofIdx, timeIdx), phaseIdx);
93  }
94 
98  template <class LhsEval>
99  void computeStorage(Dune::FieldVector<LhsEval, numEq>& storage,
100  const ElementContext& elemCtx,
101  unsigned dofIdx,
102  unsigned timeIdx) const
103  {
104  storage = 0.0;
105  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
106  addPhaseStorage(storage, elemCtx, dofIdx, timeIdx, phaseIdx);
107 
108  EnergyModule::addSolidHeatStorage(storage, elemCtx.intensiveQuantities(dofIdx, timeIdx));
109  }
110 
114  void computeFlux(RateVector& flux,
115  const ElementContext& elemCtx,
116  unsigned scvfIdx,
117  unsigned timeIdx) const
118  {
119  flux = 0.0;
120  addAdvectiveFlux(flux, elemCtx, scvfIdx, timeIdx);
121  Opm::Valgrind::CheckDefined(flux);
122 
123  addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
124  Opm::Valgrind::CheckDefined(flux);
125  }
126 
130  void addAdvectiveFlux(RateVector& flux,
131  const ElementContext& elemCtx,
132  unsigned scvfIdx,
133  unsigned timeIdx) const
134  {
135  const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
136 
137  unsigned focusDofIdx = elemCtx.focusDofIndex();
138  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
139  // data attached to upstream and the downstream DOFs
140  // of the current phase
141  unsigned upIdx = static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
142  const IntensiveQuantities& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
143 
144  // this is a bit hacky because it is specific to the element-centered
145  // finite volume scheme. (N.B. that if finite differences are used to
146  // linearize the system of equations, it does not matter.)
147  if (upIdx == focusDofIdx) {
148  Evaluation tmp =
149  up.fluidState().molarDensity(phaseIdx)
150  * extQuants.volumeFlux(phaseIdx);
151 
152  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
153  flux[conti0EqIdx + compIdx] +=
154  tmp*up.fluidState().moleFraction(phaseIdx, compIdx);
155  }
156  }
157  else {
158  Evaluation tmp =
159  Toolbox::value(up.fluidState().molarDensity(phaseIdx))
160  * extQuants.volumeFlux(phaseIdx);
161 
162  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
163  flux[conti0EqIdx + compIdx] +=
164  tmp*Toolbox::value(up.fluidState().moleFraction(phaseIdx, compIdx));
165  }
166  }
167  }
168 
169  EnergyModule::addAdvectiveFlux(flux, elemCtx, scvfIdx, timeIdx);
170  }
171 
175  void addDiffusiveFlux(RateVector& flux,
176  const ElementContext& elemCtx,
177  unsigned scvfIdx,
178  unsigned timeIdx) const
179  {
180  DiffusionModule::addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
181  EnergyModule::addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
182  }
183 
187  void computeSource(RateVector& source,
188  const ElementContext& elemCtx,
189  unsigned dofIdx,
190  unsigned timeIdx) const
191  {
192  Opm::Valgrind::SetUndefined(source);
193  elemCtx.problem().source(source, elemCtx, dofIdx, timeIdx);
194  Opm::Valgrind::CheckDefined(source);
195  }
196 };
197 
198 } // namespace Ewoms
199 
200 #endif
void computeSource(RateVector &source, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Calculate the source term of the equation.
Definition: pvslocalresidual.hh:187
void addPhaseStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx, unsigned phaseIdx) const
Adds the amount all conservation quantities (e.g.
Definition: pvslocalresidual.hh:74
void computeStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Evaluate the amount all conservation quantities (e.g.
Definition: pvslocalresidual.hh:99
Provides the auxiliary methods required for consideration of the energy equation. ...
Definition: energymodule.hh:59
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
void addAdvectiveFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx) const
Add the advective mass flux at a given flux integration point.
Definition: pvslocalresidual.hh:130
Contains the classes required to consider energy as a conservation quantity in a multi-phase module...
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
Element-wise calculation of the local residual for the compositional multi-phase primary variable swi...
Definition: pvslocalresidual.hh:47
Provides the auxiliary methods required for consideration of the diffusion equation.
Definition: diffusionmodule.hh:53
void addDiffusiveFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx) const
Adds the diffusive flux at a given flux integration point.
Definition: pvslocalresidual.hh:175
Declares the properties required for the compositional multi-phase primary variable switching model...
Classes required for molecular diffusion.
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: pvslocalresidual.hh:114