All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
flashmodel.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_MODEL_HH
29 #define EWOMS_FLASH_MODEL_HH
30 
31 #include <opm/material/densead/Math.hpp>
32 
33 #include "flashproperties.hh"
34 #include "flashprimaryvariables.hh"
35 #include "flashlocalresidual.hh"
36 #include "flashratevector.hh"
40 #include "flashindices.hh"
41 
47 #include <opm/material/fluidmatrixinteractions/NullMaterial.hpp>
48 #include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
49 #include <opm/material/constraintsolvers/NcpFlash.hpp>
50 
51 #include <sstream>
52 #include <string>
53 
54 namespace Ewoms {
55 template <class TypeTag>
56 class FlashModel;
57 }
58 
59 namespace Ewoms {
60 namespace Properties {
63  VtkComposition,
64  VtkEnergy,
65  VtkDiffusion));
66 
68 SET_TYPE_PROP(FlashModel, LocalResidual,
70 
72 SET_TYPE_PROP(FlashModel, FlashSolver,
73  Opm::NcpFlash<typename GET_PROP_TYPE(TypeTag, Scalar),
74  typename GET_PROP_TYPE(TypeTag, FluidSystem)>);
75 
77 SET_SCALAR_PROP(FlashModel, FlashTolerance, -1.0);
78 
81 
84 
87 
90 
93 
96 
99 
100 // The updates of intensive quantities tend to be _very_ expensive for this
101 // model, so let's try to minimize the number of required ones
102 SET_BOOL_PROP(FlashModel, EnableIntensiveQuantityCache, true);
103 
104 // since thermodynamic hints are basically free if the cache for intensive quantities is
105 // enabled, and this model usually shows quite a performance improvment if they are
106 // enabled, let's enable them by default.
107 SET_BOOL_PROP(FlashModel, EnableThermodynamicHints, true);
108 
109 // disable molecular diffusion by default
110 SET_BOOL_PROP(FlashModel, EnableDiffusion, false);
111 
113 SET_BOOL_PROP(FlashModel, EnableEnergy, false);
114 } // namespace Properties
115 
174 template <class TypeTag>
175 class FlashModel
176  : public MultiPhaseBaseModel<TypeTag>
177 {
178  typedef MultiPhaseBaseModel<TypeTag> ParentType;
179 
180  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
181  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
182  typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
183 
184  typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
185 
186  enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
187  enum { enableDiffusion = GET_PROP_VALUE(TypeTag, EnableDiffusion) };
188  enum { enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy) };
189 
190 
191  typedef Ewoms::EnergyModule<TypeTag, enableEnergy> EnergyModule;
192 
193 public:
194  FlashModel(Simulator& simulator)
195  : ParentType(simulator)
196  {}
197 
201  static void registerParameters()
202  {
204 
205  // register runtime parameters of the VTK output modules
207 
208  if (enableDiffusion)
210 
211  if (enableEnergy)
213 
214  EWOMS_REGISTER_PARAM(TypeTag, Scalar, FlashTolerance,
215  "The maximum tolerance for the flash solver to "
216  "consider the solution converged");
217  }
218 
222  static std::string name()
223  { return "flash"; }
224 
228  std::string primaryVarName(unsigned pvIdx) const
229  {
230  const std::string& tmp = EnergyModule::primaryVarName(pvIdx);
231  if (tmp != "")
232  return tmp;
233 
234  std::ostringstream oss;
235  if (Indices::cTot0Idx <= pvIdx && pvIdx < Indices::cTot0Idx
236  + numComponents)
237  oss << "c_tot," << FluidSystem::componentName(/*compIdx=*/pvIdx
238  - Indices::cTot0Idx);
239  else
240  assert(false);
241 
242  return oss.str();
243  }
244 
248  std::string eqName(unsigned eqIdx) const
249  {
250  const std::string& tmp = EnergyModule::eqName(eqIdx);
251  if (tmp != "")
252  return tmp;
253 
254  std::ostringstream oss;
255  if (Indices::conti0EqIdx <= eqIdx && eqIdx < Indices::conti0EqIdx
256  + numComponents) {
257  unsigned compIdx = eqIdx - Indices::conti0EqIdx;
258  oss << "continuity^" << FluidSystem::componentName(compIdx);
259  }
260  else
261  assert(false);
262 
263  return oss.str();
264  }
265 
269  Scalar primaryVarWeight(unsigned globalDofIdx, unsigned pvIdx) const
270  {
271  Scalar tmp = EnergyModule::primaryVarWeight(*this, globalDofIdx, pvIdx);
272  if (tmp > 0)
273  return tmp;
274 
275  unsigned compIdx = pvIdx - Indices::cTot0Idx;
276 
277  // make all kg equal. also, divide the weight of all total
278  // compositions by 100 to make the relative errors more
279  // comparable to the ones of the other models (at 10% porosity
280  // the medium is fully saturated with water at atmospheric
281  // conditions if 100 kg/m^3 are present!)
282  return FluidSystem::molarMass(compIdx) / 100.0;
283  }
284 
288  Scalar eqWeight(unsigned globalDofIdx, unsigned eqIdx) const
289  {
290  Scalar tmp = EnergyModule::eqWeight(*this, globalDofIdx, eqIdx);
291  if (tmp > 0)
292  return tmp;
293 
294  unsigned compIdx = eqIdx - Indices::conti0EqIdx;
295 
296  // make all kg equal
297  return FluidSystem::molarMass(compIdx);
298  }
299 
300  void registerOutputModules_()
301  {
302  ParentType::registerOutputModules_();
303 
304  // add the VTK output modules which are meaningful for the model
305  this->addOutputModule(new Ewoms::VtkCompositionModule<TypeTag>(this->simulator_));
306  if (enableDiffusion)
307  this->addOutputModule(new Ewoms::VtkDiffusionModule<TypeTag>(this->simulator_));
308  if (enableEnergy)
309  this->addOutputModule(new Ewoms::VtkEnergyModule<TypeTag>(this->simulator_));
310  }
311 };
312 
313 } // namespace Ewoms
314 
315 #endif
A compositional multi-phase model based on flash-calculations.
Definition: flashmodel.hh:56
Contains the intensive quantities of the flash-based compositional multi-phase model.
#define SET_BOOL_PROP(EffTypeTagName, PropTagName,...)
Set a property to a simple constant boolean value.
Definition: propertysystem.hh:361
std::string primaryVarName(unsigned pvIdx) const
Given an primary variable index, return a human readable name.
Definition: flashmodel.hh:228
This template class contains the data which is required to calculate all fluxes of components over a ...
Definition: flashextensivequantities.hh:52
static void registerParameters()
Register all run-time parameters for the immiscible model.
Definition: flashmodel.hh:201
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: vtkenergymodule.hh:101
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: vtkcompositionmodule.hh:105
This template class contains the data which is required to calculate all fluxes of components over a ...
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...
VTK output module for quantities which make sense for models which incorperate molecular diffusion...
#define NEW_TYPE_TAG(...)
Define a new type tag.
Definition: propertysystem.hh:169
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
Calculates the local residual of the compositional multi-phase model based on flash calculations...
Definition: flashlocalresidual.hh:46
Represents the primary variables used by the compositional flow model based on flash calculations...
#define INHERITS_FROM(...)
Syntactic sugar for NEW_TYPE_TAG.
Definition: propertysystem.hh:230
Contains the intensive quantities of the flash-based compositional multi-phase model.
Definition: flashintensivequantities.hh:52
#define EWOMS_REGISTER_PARAM(TypeTag, ParamType, ParamName, Description)
Register a run-time parameter.
Definition: parametersystem.hh:68
Implements a vector representing rates of conserved quantities.
Definition: flashratevector.hh:47
Declares the properties required by the compositional multi-phase model based on flash calculations...
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: vtkdiffusionmodule.hh:101
VTK output module for quantities which make sense for models which assume thermal equilibrium...
VTK output module for the fluid composition.
Scalar eqWeight(unsigned globalDofIdx, unsigned eqIdx) const
Returns the relative weight of an equation.
Definition: flashmodel.hh:288
static std::string name()
Definition: flashmodel.hh:222
Scalar primaryVarWeight(unsigned globalDofIdx, unsigned pvIdx) const
Returns the relative weight of a primary variable for calculating relative errors.
Definition: flashmodel.hh:269
Implements a boundary vector for the fully implicit compositional multi-phase model which is based on...
VTK output module for quantities which make sense for models which incorperate molecular diffusion...
Definition: vtkdiffusionmodule.hh:73
std::string eqName(unsigned eqIdx) const
Given an equation index, return a human readable name.
Definition: flashmodel.hh:248
Defines the primary variable and equation indices for the compositional multi-phase model based on fl...
static void registerParameters()
Register all run-time parameters for the immiscible model.
Definition: multiphasebasemodel.hh:147
Represents the primary variables used by the compositional flow model based on flash calculations...
Definition: flashprimaryvariables.hh:58
Defines the primary variable and equation indices for the compositional multi-phase model based on fl...
Definition: flashindices.hh:45
A base class for fully-implicit multi-phase porous-media flow models which assume multiple fluid phas...
A base class for fully-implicit multi-phase porous-media flow models which assume multiple fluid phas...
Definition: multiphasebasemodel.hh:47
VTK output module for quantities which make sense for models which assume thermal equilibrium...
Definition: vtkenergymodule.hh:73
Calculates the local residual of the compositional multi-phase model based on flash calculations...
VTK output module for the fluid composition.
Definition: vtkcompositionmodule.hh:77
#define SET_TYPE_PROP(EffTypeTagName, PropTagName,...)
Set a property which defines a type.
Definition: propertysystem.hh:377
Implements a vector representing rates of conserved quantities.
#define SET_SCALAR_PROP(EffTypeTagName, PropTagName,...)
Set a property to a simple constant scalar value.
Definition: propertysystem.hh:394
Implements a boundary vector for the fully implicit compositional multi-phase model which is based on...
Definition: flashboundaryratevector.hh:46