immiscibleboundaryratevector.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_IMMISCIBLE_BOUNDARY_RATE_VECTOR_HH
29 #define EWOMS_IMMISCIBLE_BOUNDARY_RATE_VECTOR_HH
30 
31 #include <opm/common/Valgrind.hpp>
32 #include <opm/material/constraintsolvers/NcpFlash.hpp>
33 
35 
36 namespace Ewoms {
37 
44 template <class TypeTag>
45 class ImmiscibleBoundaryRateVector : public GET_PROP_TYPE(TypeTag, RateVector)
46 {
47  typedef typename GET_PROP_TYPE(TypeTag, RateVector) ParentType;
48  typedef typename GET_PROP_TYPE(TypeTag, ExtensiveQuantities) ExtensiveQuantities;
49  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
50  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
51  typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
52  typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
53 
54  enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
55  enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
56  enum { conti0EqIdx = Indices::conti0EqIdx };
57  enum { enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy) };
58 
59  typedef Opm::MathToolbox<Evaluation> Toolbox;
61 
62 public:
64  : ParentType()
65  {}
66 
73  ImmiscibleBoundaryRateVector(const Evaluation& value)
74  : ParentType(value)
75  {}
76 
83 
84  ImmiscibleBoundaryRateVector& operator=(const ImmiscibleBoundaryRateVector& value) = default;
85 
97  template <class Context, class FluidState>
98  void setFreeFlow(const Context& context, unsigned bfIdx, unsigned timeIdx, const FluidState& fluidState)
99  {
100  typename FluidSystem::template ParameterCache<typename FluidState::Scalar> paramCache;
101  paramCache.updateAll(fluidState);
102 
103  ExtensiveQuantities extQuants;
104  extQuants.updateBoundary(context, bfIdx, timeIdx, fluidState, paramCache);
105  const auto& insideIntQuants = context.intensiveQuantities(bfIdx, timeIdx);
106 
108  // advective fluxes of all components in all phases
110  (*this) = Evaluation(0.0);
111  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
112  Evaluation density;
113  if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
114  density = FluidSystem::density(fluidState, paramCache, phaseIdx);
115  else
116  density = insideIntQuants.fluidState().density(phaseIdx);
117 
118  Opm::Valgrind::CheckDefined(density);
119  Opm::Valgrind::CheckDefined(extQuants.volumeFlux(phaseIdx));
120 
121  // add advective flux of current component in current
122  // phase
123  (*this)[conti0EqIdx + phaseIdx] += extQuants.volumeFlux(phaseIdx)*density;
124 
125  if (enableEnergy) {
126  Evaluation specificEnthalpy;
127  Scalar pBoundary = fluidState.pressure(phaseIdx);
128  const Evaluation& pElement = insideIntQuants.fluidState().pressure(phaseIdx);
129  if (pBoundary > pElement)
130  specificEnthalpy = FluidSystem::enthalpy(fluidState, paramCache, phaseIdx);
131  else
132  specificEnthalpy = insideIntQuants.fluidState().enthalpy(phaseIdx);
133 
134  // currently we neglect heat conduction!
135  Evaluation enthalpyRate = density*extQuants.volumeFlux(phaseIdx)*specificEnthalpy;
136  EnergyModule::addToEnthalpyRate(*this, enthalpyRate);
137  }
138  }
139 
140  EnergyModule::addToEnthalpyRate(*this, EnergyModule::heatConductionRate(extQuants));
141 
142 #ifndef NDEBUG
143  for (unsigned i = 0; i < numEq; ++i)
144  Opm::Valgrind::CheckDefined((*this)[i]);
145  Opm::Valgrind::CheckDefined(*this);
146 #endif
147  }
148 
158  template <class Context, class FluidState>
159  void setInFlow(const Context& context,
160  unsigned bfIdx,
161  unsigned timeIdx,
162  const FluidState& fluidState)
163  {
164  this->setFreeFlow(context, bfIdx, timeIdx, fluidState);
165 
166  // we only allow fluxes in the direction opposite to the outer unit normal
167  for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
168  Evaluation& val = this->operator[](eqIdx);
169  val = Toolbox::min(0.0, val);
170  }
171  }
172 
182  template <class Context, class FluidState>
183  void setOutFlow(const Context& context,
184  unsigned bfIdx,
185  unsigned timeIdx,
186  const FluidState& fluidState)
187  {
188  this->setFreeFlow(context, bfIdx, timeIdx, fluidState);
189 
190  // we only allow fluxes in the same direction as the outer unit normal
191  for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
192  Evaluation& val = this->operator[](eqIdx);
193  val = Toolbox::max(0.0, val);
194  }
195  }
196 
200  void setNoFlow()
201  { (*this) = Evaluation(0.0); }
202 };
203 
204 } // namespace Ewoms
205 
206 #endif
Contains the quantities which are are constant within a finite volume for the immiscible multi-phase ...
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
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
void setNoFlow()
Specify a no-flow boundary for all conserved quantities.
Definition: immiscibleboundaryratevector.hh:200
Implements a boundary vector for the fully implicit multi-phase model which assumes immiscibility...
Definition: immiscibleboundaryratevector.hh:45
void setOutFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an outflow boundary.
Definition: immiscibleboundaryratevector.hh:183
void setInFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an inflow boundary.
Definition: immiscibleboundaryratevector.hh:159
ImmiscibleBoundaryRateVector(const Evaluation &value)
Constructor that assigns all entries to a scalar value.
Definition: immiscibleboundaryratevector.hh:73
void setFreeFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify a free-flow boundary.
Definition: immiscibleboundaryratevector.hh:98