All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
multiphasebaseextensivequantities.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_MULTI_PHASE_BASE_EXTENSIVE_QUANTITIES_HH
29 #define EWOMS_MULTI_PHASE_BASE_EXTENSIVE_QUANTITIES_HH
30 
32 
36 
37 #include <opm/common/Valgrind.hpp>
38 #include <opm/common/Unused.hpp>
39 
40 #include <dune/common/fvector.hh>
41 
42 namespace Ewoms {
49 template <class TypeTag>
51  : public GET_PROP_TYPE(TypeTag, DiscExtensiveQuantities)
52  , public GET_PROP_TYPE(TypeTag, FluxModule)::FluxExtensiveQuantities
53 {
54  typedef typename GET_PROP_TYPE(TypeTag, DiscExtensiveQuantities) ParentType;
55  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
56  typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
57  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
58 
59  enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
60 
61  typedef typename GET_PROP_TYPE(TypeTag, FluxModule) FluxModule;
62  typedef typename FluxModule::FluxExtensiveQuantities FluxExtensiveQuantities;
63 
64 public:
68  static void registerParameters()
69  {
70  FluxModule::registerParameters();
71  }
72 
81  void update(const ElementContext& elemCtx, unsigned scvfIdx, unsigned timeIdx)
82  {
83  ParentType::update(elemCtx, scvfIdx, timeIdx);
84 
85  // compute the pressure potential gradients
86  FluxExtensiveQuantities::calculateGradients_(elemCtx, scvfIdx, timeIdx);
87 
88  // Check whether the pressure potential is in the same direction as the face
89  // normal or in the opposite one
90  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
91  if (!elemCtx.model().phaseIsConsidered(phaseIdx)) {
92  Opm::Valgrind::SetUndefined(upstreamScvIdx_[phaseIdx]);
93  Opm::Valgrind::SetUndefined(downstreamScvIdx_[phaseIdx]);
94  continue;
95  }
96 
97  upstreamScvIdx_[phaseIdx] = FluxExtensiveQuantities::upstreamIndex_(phaseIdx);
98  downstreamScvIdx_[phaseIdx] = FluxExtensiveQuantities::downstreamIndex_(phaseIdx);
99  }
100 
101  FluxExtensiveQuantities::calculateFluxes_(elemCtx, scvfIdx, timeIdx);
102  }
103 
104 
115  template <class Context, class FluidState>
116  void updateBoundary(const Context& context,
117  unsigned bfIdx,
118  unsigned timeIdx,
119  const FluidState& fluidState,
120  typename FluidSystem::template ParameterCache<typename FluidState::Scalar>& paramCache)
121  {
122  ParentType::updateBoundary(context, bfIdx, timeIdx, fluidState, paramCache);
123 
124  FluxExtensiveQuantities::calculateBoundaryGradients_(context.elementContext(),
125  bfIdx,
126  timeIdx,
127  fluidState,
128  paramCache);
129  FluxExtensiveQuantities::calculateBoundaryFluxes_(context.elementContext(),
130  bfIdx,
131  timeIdx);
132  }
133 
141  short upstreamIndex(unsigned phaseIdx) const
142  { return upstreamScvIdx_[phaseIdx]; }
143 
151  short downstreamIndex(unsigned phaseIdx) const
152  { return downstreamScvIdx_[phaseIdx]; }
153 
160  Scalar upstreamWeight(unsigned phaseIdx OPM_UNUSED) const
161  { return 1.0; }
162 
169  Scalar downstreamWeight(unsigned phaseIdx) const
170  { return 1.0 - upstreamWeight(phaseIdx); }
171 
172 private:
173  short upstreamScvIdx_[numPhases];
174  short downstreamScvIdx_[numPhases];
175 };
176 
177 } // namespace Ewoms
178 
179 #endif
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
Defines the common properties required by the porous medium multi-phase models.
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
short downstreamIndex(unsigned phaseIdx) const
Return the local index of the downstream control volume for a given phase as a function of the normal...
Definition: multiphasebaseextensivequantities.hh:151
This file provides the infrastructure to retrieve run-time parameters.
void updateBoundary(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState, typename FluidSystem::template ParameterCache< typename FluidState::Scalar > &paramCache)
Update the extensive quantities for a given boundary face.
Definition: multiphasebaseextensivequantities.hh:116
This method contains all callback classes for quantities that are required by some extensive quantiti...
void update(const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Update the extensive quantities for a given sub-control-volume-face.
Definition: multiphasebaseextensivequantities.hh:81
static void registerParameters()
Register all run-time parameters for the extensive quantities.
Definition: multiphasebaseextensivequantities.hh:68
Scalar upstreamWeight(unsigned phaseIdx OPM_UNUSED) const
Return the weight of the upstream control volume for a given phase as a function of the normal flux...
Definition: multiphasebaseextensivequantities.hh:160
short upstreamIndex(unsigned phaseIdx) const
Return the local index of the upstream control volume for a given phase as a function of the normal f...
Definition: multiphasebaseextensivequantities.hh:141
This class calculates the pressure potential gradients and the filter velocities for multi-phase flow...
Definition: multiphasebaseextensivequantities.hh:50
Provide the properties at a face which make sense indepentently of the conserved quantities.
Scalar downstreamWeight(unsigned phaseIdx) const
Return the weight of the downstream control volume for a given phase as a function of the normal flux...
Definition: multiphasebaseextensivequantities.hh:169