richardsboundaryratevector.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_BOUNDARY_RATE_VECTOR_HH
29 #define EWOMS_RICHARDS_BOUNDARY_RATE_VECTOR_HH
30 
31 #include <opm/common/Valgrind.hpp>
32 #include <opm/material/constraintsolvers/NcpFlash.hpp>
33 
35 
36 namespace Ewoms {
37 
43 template <class TypeTag>
44 class RichardsBoundaryRateVector : public GET_PROP_TYPE(TypeTag, RateVector)
45 {
46  typedef typename GET_PROP_TYPE(TypeTag, RateVector) ParentType;
47  typedef typename GET_PROP_TYPE(TypeTag, ExtensiveQuantities) ExtensiveQuantities;
48  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
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, Indices) Indices;
52 
53  enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
54  enum { contiEqIdx = Indices::contiEqIdx };
55  enum { liquidPhaseIdx = GET_PROP_VALUE(TypeTag, LiquidPhaseIndex) };
56 
57  typedef Opm::MathToolbox<Evaluation> Toolbox;
58 
59 public:
60  RichardsBoundaryRateVector() : ParentType()
61  {}
62 
67  RichardsBoundaryRateVector(const Evaluation& value)
68  : ParentType(value)
69  {}
70 
76  RichardsBoundaryRateVector& operator=(const RichardsBoundaryRateVector& value) = default;
77 
81  template <class Context, class FluidState>
82  void setFreeFlow(const Context& context, unsigned bfIdx, unsigned timeIdx, const FluidState& fluidState)
83  {
84  typename FluidSystem::template ParameterCache<typename FluidState::Scalar> paramCache;
85  paramCache.updateAll(fluidState);
86 
87  ExtensiveQuantities extQuants;
88  extQuants.updateBoundary(context, bfIdx, timeIdx, fluidState, paramCache);
89  const auto& insideIntQuants = context.intensiveQuantities(bfIdx, timeIdx);
90 
92  // advective fluxes of all components in all phases
94  (*this) = Evaluation(0.0);
95 
96  unsigned phaseIdx = liquidPhaseIdx;
97  Evaluation density;
98  if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
99  density = FluidSystem::density(fluidState, paramCache, phaseIdx);
100  else
101  density = insideIntQuants.fluidState().density(phaseIdx);
102 
103  // add advective flux of current component in current
104  // phase
105  (*this)[contiEqIdx] += extQuants.volumeFlux(phaseIdx) * density;
106 
107 #ifndef NDEBUG
108  for (unsigned i = 0; i < numEq; ++i) {
109  Opm::Valgrind::CheckDefined((*this)[i]);
110  }
111  Opm::Valgrind::CheckDefined(*this);
112 #endif
113  }
114 
118  template <class Context, class FluidState>
119  void setInFlow(const Context& context,
120  unsigned bfIdx,
121  unsigned timeIdx,
122  const FluidState& fluidState)
123  {
124  this->setFreeFlow(context, bfIdx, timeIdx, fluidState);
125 
126  // we only allow fluxes in the direction opposite to the outer unit normal
127  for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
128  Evaluation& val = this->operator[](eqIdx);
129  val = Toolbox::min(0.0, val);
130  }
131  }
132 
136  template <class Context, class FluidState>
137  void setOutFlow(const Context& context,
138  unsigned bfIdx,
139  unsigned timeIdx,
140  const FluidState& fluidState)
141  {
142  this->setFreeFlow(context, bfIdx, timeIdx, fluidState);
143 
144  // we only allow fluxes in the same direction as the outer unit normal
145  for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
146  Evaluation& val = this->operator[](eqIdx);
147  val = Toolbox::max(0.0, val);
148  }
149  }
150 
154  void setNoFlow()
155  { (*this) = Evaluation(0.0); }
156 };
157 
158 } // namespace Ewoms
159 
160 #endif
Intensive quantities required by the Richards model.
Definition: baseauxiliarymodule.hh:37
RichardsBoundaryRateVector(const Evaluation &value)
Definition: richardsboundaryratevector.hh:67
#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 setOutFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an outflow boundary.
Definition: richardsboundaryratevector.hh:137
void setInFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an inflow boundary.
Definition: richardsboundaryratevector.hh:119
void setFreeFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify a free-flow boundary.
Definition: richardsboundaryratevector.hh:82
Implements a boundary vector for the fully implicit Richards model.
Definition: richardsboundaryratevector.hh:44
void setNoFlow()
Specify a no-flow boundary for all conserved quantities.
Definition: richardsboundaryratevector.hh:154