All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
richardslocalresidual.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_RICHARDS_LOCAL_RESIDUAL_HH
29 #define EWOMS_RICHARDS_LOCAL_RESIDUAL_HH
30 
32 
34 
35 namespace Ewoms {
36 
41 template <class TypeTag>
42 class RichardsLocalResidual : public GET_PROP_TYPE(TypeTag, DiscLocalResidual)
43 {
44  typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
45  typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
46  typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
47  typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
48  typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
49  typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
50 
51  enum { contiEqIdx = Indices::contiEqIdx };
52  enum { liquidPhaseIdx = GET_PROP_VALUE(TypeTag, LiquidPhaseIndex) };
53  enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
54  typedef Opm::MathToolbox<Evaluation> Toolbox;
55 
56 public:
60  template <class LhsEval>
61  void computeStorage(Dune::FieldVector<LhsEval, numEq>& storage,
62  const ElementContext& elemCtx,
63  unsigned dofIdx,
64  unsigned timeIdx) const
65  {
66  const IntensiveQuantities& intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
67 
68  // partial time derivative of the wetting phase mass
69  storage[contiEqIdx] =
70  Toolbox::template decay<LhsEval>(intQuants.fluidState().density(liquidPhaseIdx))
71  *Toolbox::template decay<LhsEval>(intQuants.fluidState().saturation(liquidPhaseIdx))
72  *Toolbox::template decay<LhsEval>(intQuants.porosity());
73  }
74 
78  void computeFlux(RateVector& flux,
79  const ElementContext& elemCtx,
80  unsigned scvfIdx,
81  unsigned timeIdx) const
82  {
83  const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
84 
85  unsigned focusDofIdx = elemCtx.focusDofIndex();
86  unsigned upIdx = static_cast<unsigned>(extQuants.upstreamIndex(liquidPhaseIdx));
87 
88  const IntensiveQuantities& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
89 
90  // compute advective mass flux of the liquid phase. This is slightly hacky
91  // because it is specific to the element-centered finite volume method.
92  const Evaluation& rho = up.fluidState().density(liquidPhaseIdx);
93  if (focusDofIdx == upIdx)
94  flux[contiEqIdx] = extQuants.volumeFlux(liquidPhaseIdx)*rho;
95  else
96  flux[contiEqIdx] = extQuants.volumeFlux(liquidPhaseIdx)*Toolbox::value(rho);
97  }
98 
102  void computeSource(RateVector& source,
103  const ElementContext& elemCtx,
104  unsigned dofIdx,
105  unsigned timeIdx) const
106  { elemCtx.problem().source(source, elemCtx, dofIdx, timeIdx); }
107 };
108 
109 } // namespace Ewoms
110 
111 #endif
Intensive quantities required by the Richards model.
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
void computeStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Evaluate the amount all conservation quantities (e.g.
Definition: richardslocalresidual.hh:61
Element-wise calculation of the residual for the Richards model.
Definition: richardslocalresidual.hh:42
Calculates and stores the data which is required to calculate the flux of fluid over a face of a fini...
void computeSource(RateVector &source, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Calculate the source term of the equation.
Definition: richardslocalresidual.hh:102
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: richardslocalresidual.hh:78