pvsboundaryratevector.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_PVS_BOUNDARY_RATE_VECTOR_HH
29 #define EWOMS_PVS_BOUNDARY_RATE_VECTOR_HH
30 
31 #include "pvsproperties.hh"
32 
34 #include <opm/common/Valgrind.hpp>
35 
36 namespace Ewoms {
37 
45 template <class TypeTag>
46 class PvsBoundaryRateVector : public GET_PROP_TYPE(TypeTag, RateVector)
47 {
48  typedef typename GET_PROP_TYPE(TypeTag, RateVector) ParentType;
49  typedef typename GET_PROP_TYPE(TypeTag, ExtensiveQuantities) ExtensiveQuantities;
50  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
51  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
52  typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
53  typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
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  enum { enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy) };
60 
62  typedef Opm::MathToolbox<Evaluation> Toolbox;
63 
64 public:
66  : ParentType()
67  {}
68 
73  PvsBoundaryRateVector(const Evaluation& value)
74  : ParentType(value)
75  {}
76 
81  PvsBoundaryRateVector(const PvsBoundaryRateVector& value) = default;
82  PvsBoundaryRateVector& operator=(const PvsBoundaryRateVector& value) = default;
83 
87  template <class Context, class FluidState>
88  void setFreeFlow(const Context& context, unsigned bfIdx, unsigned timeIdx, const FluidState& fluidState)
89  {
90  typename FluidSystem::template ParameterCache<typename FluidState::Scalar> paramCache;
91  paramCache.updateAll(fluidState);
92 
93  ExtensiveQuantities extQuants;
94  extQuants.updateBoundary(context, bfIdx, timeIdx, fluidState, paramCache);
95  const auto& insideIntQuants = context.intensiveQuantities(bfIdx, timeIdx);
96 
98  // advective fluxes of all components in all phases
100  (*this) = Evaluation(0.0);
101  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
102  Evaluation meanMBoundary = 0;
103  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
104  meanMBoundary +=
105  fluidState.moleFraction(phaseIdx, compIdx)*FluidSystem::molarMass(compIdx);
106 
107  Evaluation density;
108  if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
109  density = FluidSystem::density(fluidState, paramCache, phaseIdx);
110  else
111  density = insideIntQuants.fluidState().density(phaseIdx);
112 
113  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
114  Evaluation molarity;
115  if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
116  molarity = fluidState.moleFraction(phaseIdx, compIdx)*density/meanMBoundary;
117  else
118  molarity = insideIntQuants.fluidState().molarity(phaseIdx, compIdx);
119 
120  // add advective flux of current component in current
121  // phase
122  (*this)[conti0EqIdx + compIdx] += extQuants.volumeFlux(phaseIdx)*molarity;
123  }
124 
125  if (enableEnergy) {
126  Evaluation specificEnthalpy;
127  if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
128  specificEnthalpy = FluidSystem::enthalpy(fluidState, paramCache, phaseIdx);
129  else
130  specificEnthalpy = insideIntQuants.fluidState().enthalpy(phaseIdx);
131 
132  // currently we neglect heat conduction!
133  Evaluation enthalpyRate = density*extQuants.volumeFlux(phaseIdx)*specificEnthalpy;
134  EnergyModule::addToEnthalpyRate(*this, enthalpyRate);
135  }
136  }
137 
138  EnergyModule::addToEnthalpyRate(*this, EnergyModule::heatConductionRate(extQuants));
139 
140 #ifndef NDEBUG
141  for (unsigned i = 0; i < numEq; ++i)
142  Opm::Valgrind::CheckDefined((*this)[i]);
143 #endif
144  }
145 
149  template <class Context, class FluidState>
150  void setInFlow(const Context& context,
151  unsigned bfIdx,
152  unsigned timeIdx,
153  const FluidState& fluidState)
154  {
155  this->setFreeFlow(context, bfIdx, timeIdx, fluidState);
156 
157  // we only allow fluxes in the direction opposite to the outer unit normal
158  for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
159  Evaluation& val = this->operator[](eqIdx);
160  val = Toolbox::min(0.0, val);
161  }
162  }
163 
167  template <class Context, class FluidState>
168  void setOutFlow(const Context& context,
169  unsigned bfIdx,
170  unsigned timeIdx,
171  const FluidState& fluidState)
172  {
173  this->setFreeFlow(context, bfIdx, timeIdx, fluidState);
174 
175  // we only allow fluxes in the same direction as the outer unit normal
176  for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
177  Evaluation& val = this->operator[](eqIdx);
178  val = Toolbox::max(0.0, val);
179  }
180  }
181 
185  void setNoFlow()
186  { (*this) = Evaluation(0.0); }
187 };
188 
189 } // namespace Ewoms
190 
191 #endif
PvsBoundaryRateVector(const Evaluation &value)
Definition: pvsboundaryratevector.hh:73
Definition: baseauxiliarymodule.hh:37
void setNoFlow()
Specify a no-flow boundary for all conserved quantities.
Definition: pvsboundaryratevector.hh:185
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 setOutFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an outflow boundary.
Definition: pvsboundaryratevector.hh:168
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
void setFreeFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify a free-flow boundary.
Definition: pvsboundaryratevector.hh:88
void setInFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an inflow boundary.
Definition: pvsboundaryratevector.hh:150
Declares the properties required for the compositional multi-phase primary variable switching model...
Implements a rate vector on the boundary for the fully implicit compositional multi-phase primary var...
Definition: pvsboundaryratevector.hh:46