flash/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/constraintsolvers/NcpFlash.hpp>
32
33#include <opm/material/densead/Math.hpp>
34
35#include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
36#include <opm/material/fluidmatrixinteractions/NullMaterial.hpp>
37
40
50
54
55#include <sstream>
56#include <string>
57
58namespace Opm {
59template <class TypeTag>
60class FlashModel;
61}
62
63namespace Opm::Properties {
64
65namespace TTag {
67struct FlashModel { using InheritsFrom = std::tuple<MultiPhaseBaseModel>; };
68} // namespace TTag
69
71template<class TypeTag>
73
75template<class TypeTag>
76struct FlashSolver<TypeTag, TTag::FlashModel>
77{ using type = Opm::NcpFlash<GetPropType<TypeTag, Properties::Scalar>,
79
81template<class TypeTag>
82struct Model<TypeTag, TTag::FlashModel> { using type = Opm::FlashModel<TypeTag>; };
83
85template<class TypeTag>
87
89template<class TypeTag>
90struct RateVector<TypeTag, TTag::FlashModel> { using type = Opm::FlashRateVector<TypeTag>; };
91
93template<class TypeTag>
95
97template<class TypeTag>
99
101template<class TypeTag>
103
105template<class TypeTag>
106struct Indices<TypeTag, TTag::FlashModel> { using type = Opm::FlashIndices<TypeTag, /*PVIdx=*/0>; };
107
108// disable molecular diffusion by default
109template<class TypeTag>
110struct EnableDiffusion<TypeTag, TTag::FlashModel> { static constexpr bool value = false; };
111
113template<class TypeTag>
114struct EnableEnergy<TypeTag, TTag::FlashModel> { static constexpr bool value = false; };
115
116} // namespace Opm::Properties
117
118namespace Opm {
119
179template <class TypeTag>
180class FlashModel
181 : public MultiPhaseBaseModel<TypeTag>
182{
183 using ParentType = MultiPhaseBaseModel<TypeTag>;
184
188
190
191 enum { numComponents = getPropValue<TypeTag, Properties::NumComponents>() };
192 enum { enableDiffusion = getPropValue<TypeTag, Properties::EnableDiffusion>() };
193 enum { enableEnergy = getPropValue<TypeTag, Properties::EnableEnergy>() };
194
195
196 using EnergyModule = Opm::EnergyModule<TypeTag, enableEnergy>;
197
198public:
199 FlashModel(Simulator& simulator)
200 : ParentType(simulator)
201 {}
202
206 static void registerParameters()
207 {
209
210 // register runtime parameters of the VTK output modules
212
213 if (enableDiffusion)
215
216 if (enableEnergy)
218
219 Parameters::Register<Parameters::FlashTolerance<Scalar>>
220 ("The maximum tolerance for the flash solver to "
221 "consider the solution converged");
222
223 // The updates of intensive quantities tend to be _very_ expensive for this
224 // model, so let's try to minimize the number of required ones
225 Parameters::SetDefault<Parameters::EnableIntensiveQuantityCache>(true);
226
227 // since thermodynamic hints are basically free if the cache for intensive quantities is
228 // enabled, and this model usually shows quite a performance improvment if they are
229 // enabled, let's enable them by default.
230 Parameters::SetDefault<Parameters::EnableThermodynamicHints>(true);
231 }
232
236 static std::string name()
237 { return "flash"; }
238
242 std::string primaryVarName(unsigned pvIdx) const
243 {
244 const std::string& tmp = EnergyModule::primaryVarName(pvIdx);
245 if (tmp != "")
246 return tmp;
247
248 std::ostringstream oss;
249 if (Indices::cTot0Idx <= pvIdx && pvIdx < Indices::cTot0Idx
250 + numComponents)
251 oss << "c_tot," << FluidSystem::componentName(/*compIdx=*/pvIdx
252 - Indices::cTot0Idx);
253 else
254 assert(false);
255
256 return oss.str();
257 }
258
262 std::string eqName(unsigned eqIdx) const
263 {
264 const std::string& tmp = EnergyModule::eqName(eqIdx);
265 if (tmp != "")
266 return tmp;
267
268 std::ostringstream oss;
269 if (Indices::conti0EqIdx <= eqIdx && eqIdx < Indices::conti0EqIdx
270 + numComponents) {
271 unsigned compIdx = eqIdx - Indices::conti0EqIdx;
272 oss << "continuity^" << FluidSystem::componentName(compIdx);
273 }
274 else
275 assert(false);
276
277 return oss.str();
278 }
279
283 Scalar primaryVarWeight(unsigned globalDofIdx, unsigned pvIdx) const
284 {
285 Scalar tmp = EnergyModule::primaryVarWeight(*this, globalDofIdx, pvIdx);
286 if (tmp > 0)
287 return tmp;
288
289 unsigned compIdx = pvIdx - Indices::cTot0Idx;
290
291 // make all kg equal. also, divide the weight of all total
292 // compositions by 100 to make the relative errors more
293 // comparable to the ones of the other models (at 10% porosity
294 // the medium is fully saturated with water at atmospheric
295 // conditions if 100 kg/m^3 are present!)
296 return FluidSystem::molarMass(compIdx) / 100.0;
297 }
298
302 Scalar eqWeight(unsigned globalDofIdx, unsigned eqIdx) const
303 {
304 Scalar tmp = EnergyModule::eqWeight(*this, globalDofIdx, eqIdx);
305 if (tmp > 0)
306 return tmp;
307
308 unsigned compIdx = eqIdx - Indices::conti0EqIdx;
309
310 // make all kg equal
311 return FluidSystem::molarMass(compIdx);
312 }
313
315 {
317
318 // add the VTK output modules which are meaningful for the model
319 this->addOutputModule(new Opm::VtkCompositionModule<TypeTag>(this->simulator_));
320 if (enableDiffusion)
321 this->addOutputModule(new Opm::VtkDiffusionModule<TypeTag>(this->simulator_));
322 if (enableEnergy)
323 this->addOutputModule(new Opm::VtkEnergyModule<TypeTag>(this->simulator_));
324 }
325};
326
327} // namespace Opm
328
329#endif
Provides the auxiliary methods required for consideration of the energy equation.
Definition: energymodule.hh:50
Implements a boundary vector for the fully implicit compositional multi-phase model which is based on...
Definition: flashboundaryratevector.hh:45
This template class contains the data which is required to calculate all fluxes of components over a ...
Definition: flashextensivequantities.hh:54
Defines the primary variable and equation indices for the compositional multi-phase model based on fl...
Definition: flash/flashindices.hh:47
Contains the intensive quantities of the flash-based compositional multi-phase model.
Definition: flash/flashintensivequantities.hh:58
Calculates the local residual of the compositional multi-phase model based on flash calculations.
Definition: flash/flashlocalresidual.hh:46
A compositional multi-phase model based on flash-calculations.
Definition: ptflash/flashmodel.hh:187
FlashModel(Simulator &simulator)
Definition: flash/flashmodel.hh:199
static void registerParameters()
Register all run-time parameters for the immiscible model.
Definition: flash/flashmodel.hh:206
Scalar primaryVarWeight(unsigned globalDofIdx, unsigned pvIdx) const
Returns the relative weight of a primary variable for calculating relative errors.
Definition: flash/flashmodel.hh:283
static std::string name()
Definition: flash/flashmodel.hh:236
Scalar eqWeight(unsigned globalDofIdx, unsigned eqIdx) const
Returns the relative weight of an equation.
Definition: flash/flashmodel.hh:302
std::string primaryVarName(unsigned pvIdx) const
Given an primary variable index, return a human readable name.
Definition: flash/flashmodel.hh:242
void registerOutputModules_()
Definition: flash/flashmodel.hh:314
std::string eqName(unsigned eqIdx) const
Given an equation index, return a human readable name.
Definition: flash/flashmodel.hh:262
Represents the primary variables used by the compositional flow model based on flash calculations.
Definition: flash/flashprimaryvariables.hh:57
Definition: flashratevector.hh:48
A base class for fully-implicit multi-phase porous-media flow models which assume multiple fluid phas...
Definition: multiphasebasemodel.hh:153
static void registerParameters()
Register all run-time parameters for the immiscible model.
Definition: multiphasebasemodel.hh:179
void registerOutputModules_()
Definition: multiphasebasemodel.hh:254
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:102
VTK output module for the fluid composition.
Definition: vtkcompositionmodule.hh:69
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: vtkcompositionmodule.hh:96
VTK output module for quantities which make sense for models which incorperate molecular diffusion.
Definition: vtkdiffusionmodule.hh:65
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: vtkdiffusionmodule.hh:92
VTK output module for quantities which make sense for models which assume thermal equilibrium.
Definition: vtkenergymodule.hh:66
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: vtkenergymodule.hh:93
Contains the classes required to consider energy as a conservation quantity in a multi-phase module.
Declares the parameters for the compositional multi-phase model based on flash calculations.
Declares the properties required by the compositional multi-phase model based on flash calculations.
Definition: blackoilmodel.hh:72
Definition: blackoilboundaryratevector.hh:37
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:235
Type of object for specifying boundary conditions.
Definition: fvbaseproperties.hh:119
Enable diffusive fluxes?
Definition: multiphasebaseproperties.hh:79
Specify whether energy should be considered as a conservation quantity or not.
Definition: multiphasebaseproperties.hh:76
Data required to calculate a flux over a face.
Definition: fvbaseproperties.hh:149
Opm::NcpFlash< GetPropType< TypeTag, Properties::Scalar >, GetPropType< TypeTag, Properties::FluidSystem > > type
Definition: flash/flashmodel.hh:78
The type of the flash constraint solver.
Definition: flashproperties.hh:39
Enumerations used by the model.
Definition: multiphasebaseproperties.hh:48
The secondary variables within a sub-control volume.
Definition: fvbaseproperties.hh:133
The type of the local residual function.
Definition: fvbaseproperties.hh:94
The type of the model.
Definition: basicproperties.hh:88
A vector of primary variables within a sub-control volume.
Definition: fvbaseproperties.hh:130
Vector containing volumetric or areal rates of quantities.
Definition: fvbaseproperties.hh:116
The type tag for the isothermal single phase problems.
Definition: flash/flashmodel.hh:67
std::tuple< MultiPhaseBaseModel > InheritsFrom
Definition: flash/flashmodel.hh:67