flashboundaryratevector.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_FLASH_BOUNDARY_RATE_VECTOR_HH
29 #define EWOMS_FLASH_BOUNDARY_RATE_VECTOR_HH
30 
31 #include "flashproperties.hh"
32 
34 #include <opm/common/Valgrind.hpp>
35 
36 namespace Ewoms {
37 
45 template <class TypeTag>
46 class FlashBoundaryRateVector : 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:
65  FlashBoundaryRateVector() : ParentType()
66  {}
67 
72  FlashBoundaryRateVector(const Evaluation& value) : ParentType(value)
73  {}
74 
79  FlashBoundaryRateVector(const FlashBoundaryRateVector& value) = default;
80  FlashBoundaryRateVector& operator=(const FlashBoundaryRateVector& value) = default;
81 
85  template <class Context, class FluidState>
86  void setFreeFlow(const Context& context,
87  unsigned bfIdx,
88  unsigned timeIdx,
89  const FluidState& fluidState)
90  {
91  typename FluidSystem::template ParameterCache<typename FluidState::Scalar> paramCache;
92  paramCache.updateAll(fluidState);
93 
94  ExtensiveQuantities extQuants;
95  extQuants.updateBoundary(context, bfIdx, timeIdx, fluidState, paramCache);
96  const auto& insideIntQuants = context.intensiveQuantities(bfIdx, timeIdx);
97 
99  // advective fluxes of all components in all phases
101  (*this) = Evaluation(0.0);
102  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
103  Evaluation meanMBoundary = 0;
104  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
105  meanMBoundary +=
106  fluidState.moleFraction(phaseIdx, compIdx)*FluidSystem::molarMass(compIdx);
107 
108  Evaluation density;
109  if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
110  density = FluidSystem::density(fluidState, paramCache, phaseIdx);
111  else
112  density = insideIntQuants.fluidState().density(phaseIdx);
113 
114  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
115  Evaluation molarity;
116  if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
117  molarity = fluidState.moleFraction(phaseIdx, compIdx)*density/meanMBoundary;
118  else
119  molarity = insideIntQuants.fluidState().molarity(phaseIdx, compIdx);
120 
121  // add advective flux of current component in current
122  // phase
123  (*this)[conti0EqIdx + compIdx] += extQuants.volumeFlux(phaseIdx)*molarity;
124  }
125 
126  if (enableEnergy) {
127  Evaluation specificEnthalpy;
128  if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
129  specificEnthalpy = FluidSystem::enthalpy(fluidState, paramCache, phaseIdx);
130  else
131  specificEnthalpy = insideIntQuants.fluidState().enthalpy(phaseIdx);
132 
133  // currently we neglect heat conduction!
134  Evaluation enthalpyRate = density*extQuants.volumeFlux(phaseIdx)*specificEnthalpy;
135  EnergyModule::addToEnthalpyRate(*this, enthalpyRate);
136  }
137  }
138 
139  EnergyModule::addToEnthalpyRate(*this, EnergyModule::heatConductionRate(extQuants));
140 
141 #ifndef NDEBUG
142  for (unsigned i = 0; i < numEq; ++i) {
143  Opm::Valgrind::CheckDefined((*this)[i]);
144  }
145 #endif
146  }
147 
151  template <class Context, class FluidState>
152  void setInFlow(const Context& context,
153  unsigned bfIdx,
154  unsigned timeIdx,
155  const FluidState& fluidState)
156  {
157  this->setFreeFlow(context, bfIdx, timeIdx, fluidState);
158 
159  // we only allow fluxes in the direction opposite to the outer unit normal
160  for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
161  Evaluation& val = this->operator[](eqIdx);
162  val = Toolbox::min(0.0, val);
163  }
164  }
165 
169  template <class Context, class FluidState>
170  void setOutFlow(const Context& context,
171  unsigned bfIdx,
172  unsigned timeIdx,
173  const FluidState& fluidState)
174  {
175  this->setFreeFlow(context, bfIdx, timeIdx, fluidState);
176 
177  // we only allow fluxes in the same direction as the outer unit normal
178  for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
179  Evaluation& val = this->operator[](eqIdx);
180  val = Toolbox::max(0.0, val);
181  }
182  }
183 
187  void setNoFlow()
188  { (*this) = Evaluation(0.0); }
189 };
190 
191 } // namespace Ewoms
192 
193 #endif
void setFreeFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify a free-flow boundary.
Definition: flashboundaryratevector.hh:86
FlashBoundaryRateVector(const Evaluation &value)
Definition: flashboundaryratevector.hh:72
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
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 setNoFlow()
Specify a no-flow boundary for all conserved quantities.
Definition: flashboundaryratevector.hh:187
void setInFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an inflow boundary.
Definition: flashboundaryratevector.hh:152
Declares the properties required by the compositional multi-phase model based on flash calculations...
void setOutFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an outflow boundary.
Definition: flashboundaryratevector.hh:170
Implements a boundary vector for the fully implicit compositional multi-phase model which is based on...
Definition: flashboundaryratevector.hh:46