flashintensivequantities.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_INTENSIVE_QUANTITIES_HH
29 #define EWOMS_FLASH_INTENSIVE_QUANTITIES_HH
30 
31 #include "flashproperties.hh"
32 #include "flashindices.hh"
33 
36 
37 #include <opm/material/fluidstates/CompositionalFluidState.hpp>
38 #include <opm/common/Valgrind.hpp>
39 
40 #include <dune/common/fvector.hh>
41 #include <dune/common/fmatrix.hh>
42 
43 namespace Ewoms {
44 
51 template <class TypeTag>
53  : public GET_PROP_TYPE(TypeTag, DiscIntensiveQuantities)
54  , public DiffusionIntensiveQuantities<TypeTag, GET_PROP_VALUE(TypeTag, EnableDiffusion) >
55  , public EnergyIntensiveQuantities<TypeTag, GET_PROP_VALUE(TypeTag, EnableEnergy) >
56  , public GET_PROP_TYPE(TypeTag, FluxModule)::FluxIntensiveQuantities
57 {
58  typedef typename GET_PROP_TYPE(TypeTag, DiscIntensiveQuantities) ParentType;
59 
60  typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
61  typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
62  typedef typename GET_PROP_TYPE(TypeTag, MaterialLawParams) MaterialLawParams;
63  typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
64  typedef typename GET_PROP_TYPE(TypeTag, FluxModule) FluxModule;
65  typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
66  typedef typename GET_PROP_TYPE(TypeTag, ThreadManager) ThreadManager;
67 
68  // primary variable indices
69  enum { cTot0Idx = Indices::cTot0Idx };
70  enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
71  enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
72  enum { enableDiffusion = GET_PROP_VALUE(TypeTag, EnableDiffusion) };
73  enum { enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy) };
74  enum { dimWorld = GridView::dimensionworld };
75 
76  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
77  typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
78  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
79  typedef typename GET_PROP_TYPE(TypeTag, FlashSolver) FlashSolver;
80 
81  typedef Dune::FieldVector<Evaluation, numComponents> ComponentVector;
82  typedef Dune::FieldMatrix<Scalar, dimWorld, dimWorld> DimMatrix;
83 
84  typedef typename FluxModule::FluxIntensiveQuantities FluxIntensiveQuantities;
87 
88 public:
90  typedef Opm::CompositionalFluidState<Evaluation, FluidSystem, enableEnergy> FluidState;
91 
93  { }
94 
95  FlashIntensiveQuantities(const FlashIntensiveQuantities& other) = default;
96 
97  FlashIntensiveQuantities& operator=(const FlashIntensiveQuantities& other) = default;
98 
102  void update(const ElementContext& elemCtx, unsigned dofIdx, unsigned timeIdx)
103  {
104  ParentType::update(elemCtx, dofIdx, timeIdx);
105  EnergyIntensiveQuantities::updateTemperatures_(fluidState_, elemCtx, dofIdx, timeIdx);
106 
107  const auto& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
108  const auto& problem = elemCtx.problem();
109  Scalar flashTolerance = EWOMS_GET_PARAM(TypeTag, Scalar, FlashTolerance);
110 
111  // extract the total molar densities of the components
112  ComponentVector cTotal;
113  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
114  cTotal[compIdx] = priVars.makeEvaluation(cTot0Idx + compIdx, timeIdx);
115 
116  const auto *hint = elemCtx.thermodynamicHint(dofIdx, timeIdx);
117  if (hint) {
118  // use the same fluid state as the one of the hint, but
119  // make sure that we don't overwrite the temperature
120  // specified by the primary variables
121  Evaluation T = fluidState_.temperature(/*phaseIdx=*/0);
122  fluidState_.assign(hint->fluidState());
123  fluidState_.setTemperature(T);
124  }
125  else
126  FlashSolver::guessInitial(fluidState_, cTotal);
127 
128  // compute the phase compositions, densities and pressures
129  typename FluidSystem::template ParameterCache<Evaluation> paramCache;
130  const MaterialLawParams& materialParams =
131  problem.materialLawParams(elemCtx, dofIdx, timeIdx);
132  FlashSolver::template solve<MaterialLaw>(fluidState_,
133  materialParams,
134  paramCache,
135  cTotal,
136  flashTolerance);
137 
138  // calculate relative permeabilities
139  MaterialLaw::relativePermeabilities(relativePermeability_,
140  materialParams, fluidState_);
141  Opm::Valgrind::CheckDefined(relativePermeability_);
142 
143  // set the phase viscosities
144  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
145  paramCache.updatePhase(fluidState_, phaseIdx);
146 
147  const Evaluation& mu = FluidSystem::viscosity(fluidState_, paramCache, phaseIdx);
148  fluidState_.setViscosity(phaseIdx, mu);
149 
150  mobility_[phaseIdx] = relativePermeability_[phaseIdx] / mu;
151  Opm::Valgrind::CheckDefined(mobility_[phaseIdx]);
152  }
153 
155  // calculate the remaining quantities
157 
158  // porosity
159  porosity_ = problem.porosity(elemCtx, dofIdx, timeIdx);
160  Opm::Valgrind::CheckDefined(porosity_);
161 
162  // intrinsic permeability
163  intrinsicPerm_ = problem.intrinsicPermeability(elemCtx, dofIdx, timeIdx);
164 
165  // update the quantities specific for the velocity model
166  FluxIntensiveQuantities::update_(elemCtx, dofIdx, timeIdx);
167 
168  // energy related quantities
169  EnergyIntensiveQuantities::update_(fluidState_, paramCache, elemCtx, dofIdx, timeIdx);
170 
171  // update the diffusion specific quantities of the intensive quantities
172  DiffusionIntensiveQuantities::update_(fluidState_, paramCache, elemCtx, dofIdx, timeIdx);
173  }
174 
178  const FluidState& fluidState() const
179  { return fluidState_; }
180 
184  const DimMatrix& intrinsicPermeability() const
185  { return intrinsicPerm_; }
186 
190  const Evaluation& relativePermeability(unsigned phaseIdx) const
191  { return relativePermeability_[phaseIdx]; }
192 
196  const Evaluation& mobility(unsigned phaseIdx) const
197  {
198  return mobility_[phaseIdx];
199  }
200 
204  const Evaluation& porosity() const
205  { return porosity_; }
206 
207 private:
208  DimMatrix intrinsicPerm_;
209  FluidState fluidState_;
210  Evaluation porosity_;
211  Evaluation relativePermeability_[numPhases];
212  Evaluation mobility_[numPhases];
213 };
214 
215 } // namespace Ewoms
216 
217 #endif
Simplifies multi-threaded capabilities.
Definition: threadmanager.hh:52
Provides the volumetric quantities required for the calculation of molecular diffusive fluxes...
Definition: diffusionmodule.hh:146
Definition: baseauxiliarymodule.hh:37
#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
Contains the intensive quantities of the flash-based compositional multi-phase model.
Definition: flashintensivequantities.hh:52
const FluidState & fluidState() const
Returns the phase state for the control-volume.
Definition: flashintensivequantities.hh:178
Declares the properties required by the compositional multi-phase model based on flash calculations...
Opm::CompositionalFluidState< Evaluation, FluidSystem, enableEnergy > FluidState
The type of the object returned by the fluidState() method.
Definition: flashintensivequantities.hh:90
Provides the volumetric quantities required for the energy equation.
Definition: energymodule.hh:539
#define EWOMS_GET_PARAM(TypeTag, ParamType, ParamName)
Retrieve a runtime parameter.
Definition: parametersystem.hh:99
const DimMatrix & intrinsicPermeability() const
Returns the intrinsic permeability tensor a degree of freedom.
Definition: flashintensivequantities.hh:184
const Evaluation & relativePermeability(unsigned phaseIdx) const
Returns the relative permeability of a given phase within the control volume.
Definition: flashintensivequantities.hh:190
void update(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Definition: flashintensivequantities.hh:102
Defines the primary variable and equation indices for the compositional multi-phase model based on fl...
const Evaluation & porosity() const
Returns the average porosity within the control volume.
Definition: flashintensivequantities.hh:204
const Evaluation & mobility(unsigned phaseIdx) const
Returns the effective mobility of a given phase within the control volume.
Definition: flashintensivequantities.hh:196
Classes required for molecular diffusion.