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
31#include <dune/common/fvector.hh>
32
33#include <opm/material/common/MathToolbox.hpp>
34
38
39namespace Opm {
40
45template <class TypeTag>
46class RichardsLocalResidual : public GetPropType<TypeTag, Properties::DiscLocalResidual>
47{
54
55 enum { contiEqIdx = Indices::contiEqIdx };
56 enum { liquidPhaseIdx = getPropValue<TypeTag, Properties::LiquidPhaseIndex>() };
57 enum { numEq = getPropValue<TypeTag, Properties::NumEq>() };
58 using Toolbox = MathToolbox<Evaluation>;
59
60public:
64 template <class LhsEval>
65 void computeStorage(Dune::FieldVector<LhsEval, numEq>& storage,
66 const ElementContext& elemCtx,
67 unsigned dofIdx,
68 unsigned timeIdx) const
69 {
70 const IntensiveQuantities& intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
71
72 // partial time derivative of the wetting phase mass
73 storage[contiEqIdx] =
74 Toolbox::template decay<LhsEval>(intQuants.fluidState().density(liquidPhaseIdx)) *
75 Toolbox::template decay<LhsEval>(intQuants.fluidState().saturation(liquidPhaseIdx)) *
76 Toolbox::template decay<LhsEval>(intQuants.porosity());
77 }
78
82 void computeFlux(RateVector& flux,
83 const ElementContext& elemCtx,
84 unsigned scvfIdx,
85 unsigned timeIdx) const
86 {
87 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
88
89 const unsigned focusDofIdx = elemCtx.focusDofIndex();
90 const unsigned upIdx = static_cast<unsigned>(extQuants.upstreamIndex(liquidPhaseIdx));
91
92 const IntensiveQuantities& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
93
94 // compute advective mass flux of the liquid phase. This is slightly hacky
95 // because it is specific to the element-centered finite volume method.
96 const Evaluation& rho = up.fluidState().density(liquidPhaseIdx);
97 if (focusDofIdx == upIdx) {
98 flux[contiEqIdx] = extQuants.volumeFlux(liquidPhaseIdx) * rho;
99 }
100 else {
101 flux[contiEqIdx] = extQuants.volumeFlux(liquidPhaseIdx) * Toolbox::value(rho);
102 }
103 }
104
108 void computeSource(RateVector& source,
109 const ElementContext& elemCtx,
110 unsigned dofIdx,
111 unsigned timeIdx) const
112 { elemCtx.problem().source(source, elemCtx, dofIdx, timeIdx); }
113};
114
115} // namespace Opm
116
117#endif
Element-wise calculation of the residual for the Richards model.
Definition: richardslocalresidual.hh:47
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:82
void computeStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Evaluate the amount all conservation quantities (e.g. phase mass) within a finite sub-control volume.
Definition: richardslocalresidual.hh:65
void computeSource(RateVector &source, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Calculate the source term of the equation.
Definition: richardslocalresidual.hh:108
Declare the properties used by the infrastructure code of the finite volume discretizations.
Defines the common properties required by the porous medium multi-phase models.
Definition: blackoilboundaryratevector.hh:39
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:233
Contains the property declarations for the Richards model.