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