blackoilboundaryratevector.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_BLACK_OIL_BOUNDARY_RATE_VECTOR_HH
29 #define EWOMS_BLACK_OIL_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 BlackOilBoundaryRateVector : 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 { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
55  enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
56  enum { conti0EqIdx = Indices::conti0EqIdx };
57 
58 public:
62  BlackOilBoundaryRateVector() : ParentType()
63  {}
64 
68  BlackOilBoundaryRateVector(Scalar value) : ParentType(value)
69  {}
70 
75  BlackOilBoundaryRateVector& operator=(const BlackOilBoundaryRateVector& value) = default;
76 
80  template <class Context, class FluidState>
81  void setFreeFlow(const Context& context,
82  unsigned bfIdx,
83  unsigned timeIdx,
84  const FluidState& fluidState)
85  {
86  typename FluidSystem::ParameterCache paramCache;
87  paramCache.updateAll(fluidState);
88 
89  ExtensiveQuantities extQuants;
90  extQuants.updateBoundary(context, bfIdx, timeIdx, fluidState, paramCache);
91  const auto& insideIntQuants = context.intensiveQuantities(bfIdx, timeIdx);
92 
94  // advective fluxes of all components in all phases
96  (*this) = 0.0;
97  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
98  Scalar meanMBoundary = 0;
99  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
100  meanMBoundary += fluidState.moleFraction(phaseIdx, compIdx)
101  * FluidSystem::molarMass(compIdx);
102 
103  Scalar density;
104  if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
105  density = FluidSystem::density(fluidState, paramCache, phaseIdx);
106  else
107  density = insideIntQuants.fluidState().density(phaseIdx);
108 
109  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
110  Scalar molarity;
111  if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
112  molarity =
113  fluidState.moleFraction(phaseIdx, compIdx)
114  * density
115  / meanMBoundary;
116  else
117  molarity = insideIntQuants.fluidState().molarity(phaseIdx, compIdx);
118 
119  // add advective flux of current component in current
120  // phase
121  (*this)[conti0EqIdx + compIdx] += extQuants.volumeFlux(phaseIdx) * molarity;
122  }
123  }
124 
125 #ifndef NDEBUG
126  for (unsigned i = 0; i < numEq; ++i) {
127  Opm::Valgrind::CheckDefined((*this)[i]);
128  }
129  Opm::Valgrind::CheckDefined(*this);
130 #endif
131  }
132 
136  template <class Context, class FluidState>
137  void setInFlow(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 direction opposite to the outer
145  // unit normal
146  for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
147  Scalar& val = this->operator[](eqIdx);
148  val = std::min<Scalar>(0.0, val);
149  }
150  }
151 
155  template <class Context, class FluidState>
156  void setOutFlow(const Context& context,
157  unsigned bfIdx,
158  unsigned timeIdx,
159  const FluidState& fluidState)
160  {
161  this->setFreeFlow(context, bfIdx, timeIdx, fluidState);
162 
163  // we only allow fluxes in the same direction as the outer
164  // unit normal
165  for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
166  Scalar& val = this->operator[](eqIdx);
167  val = std::max( Scalar(0), val);
168  }
169  }
170 
174  void setNoFlow()
175  { (*this) = Scalar(0); }
176 };
177 
178 } // namespace Ewoms
179 
180 #endif
Definition: baseauxiliarymodule.hh:37
Contains the quantities which are are constant within a finite volume in the black-oil model...
#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
BlackOilBoundaryRateVector()
Default constructor.
Definition: blackoilboundaryratevector.hh:62
void setNoFlow()
Specify a no-flow boundary for all conserved quantities.
Definition: blackoilboundaryratevector.hh:174
Implements a boundary vector for the fully implicit black-oil model.
Definition: blackoilboundaryratevector.hh:44
void setFreeFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify a free-flow boundary.
Definition: blackoilboundaryratevector.hh:81
BlackOilBoundaryRateVector(Scalar value)
Definition: blackoilboundaryratevector.hh:68
void setOutFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an outflow boundary.
Definition: blackoilboundaryratevector.hh:156
void setInFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an inflow boundary.
Definition: blackoilboundaryratevector.hh:137