My Project
Loading...
Searching...
No Matches
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 OPM_PTFLASH_MODEL_HH
29#define OPM_PTFLASH_MODEL_HH
30
31#include <opm/material/constraintsolvers/PTFlash.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
45
50
57
58#include <sstream>
59#include <string>
60
61namespace Opm {
62
63template <class TypeTag>
64class FlashModel;
65
66}
67
68namespace Opm::Properties {
69
70namespace TTag {
72struct FlashModel { using InheritsFrom = std::tuple<MultiPhaseBaseModel>; };
73} // namespace TTag
74
76template<class TypeTag>
77struct LocalResidual<TypeTag, TTag::FlashModel>
78{ using type = Opm::FlashLocalResidual<TypeTag>; };
79
81template<class TypeTag>
82struct NewtonMethod<TypeTag, TTag::FlashModel>
84
86template<class TypeTag>
87struct FlashSolver<TypeTag, TTag::FlashModel>
88{
91};
92
94template<class TypeTag>
95struct Model<TypeTag, TTag::FlashModel>
96{ using type = Opm::FlashModel<TypeTag>; };
97
99template<class TypeTag>
100struct PrimaryVariables<TypeTag, TTag::FlashModel>
102
104template<class TypeTag>
105struct RateVector<TypeTag, TTag::FlashModel>
106{ using type = Opm::FlashRateVector<TypeTag>; };
107
109template<class TypeTag>
110struct BoundaryRateVector<TypeTag, TTag::FlashModel>
112
114template<class TypeTag>
115struct IntensiveQuantities<TypeTag, TTag::FlashModel>
117
119template<class TypeTag>
120struct ExtensiveQuantities<TypeTag, TTag::FlashModel>
122
124template<class TypeTag>
125struct Indices<TypeTag, TTag::FlashModel>
126{ using type = Opm::FlashIndices<TypeTag, /*PVIdx=*/0>; };
127
128// disable molecular diffusion by default
129template<class TypeTag>
130struct EnableDiffusion<TypeTag, TTag::FlashModel>
131{ static constexpr bool value = false; };
132
134template<class TypeTag>
135struct EnableEnergy<TypeTag, TTag::FlashModel>
136{ static constexpr bool value = false; };
137
138template<class TypeTag>
139struct EnableWater<TypeTag, TTag::MultiPhaseBaseModel>
141
142} // namespace Opm::Properties
143
144namespace Opm {
145
188template <class TypeTag>
190 : public MultiPhaseBaseModel<TypeTag>
191{
193
197
199
200 enum { numComponents = getPropValue<TypeTag, Properties::NumComponents>() };
201 enum { enableDiffusion = getPropValue<TypeTag, Properties::EnableDiffusion>() };
202 enum { enableEnergy = getPropValue<TypeTag, Properties::EnableEnergy>() };
203
204
206
207public:
208 explicit FlashModel(Simulator& simulator)
209 : ParentType(simulator)
210 {}
211
215 static void registerParameters()
216 {
218
219 // register runtime parameters of the VTK output modules
222
223 if (enableDiffusion)
225
226 if (enableEnergy)
228
229 Parameters::Register<Parameters::FlashTolerance<Scalar>>
230 ("The maximum tolerance for the flash solver to "
231 "consider the solution converged");
232 Parameters::Register<Parameters::FlashVerbosity>
233 ("Flash solver verbosity level");
234 Parameters::Register<Parameters::FlashTwoPhaseMethod>
235 ("Method for solving vapor-liquid composition. Available options include: "
236 "ssi, newton, ssi+newton");
237
238 Parameters::SetDefault<Parameters::FlashTolerance<Scalar>>(1.e-8);
239 Parameters::SetDefault<Parameters::EnableIntensiveQuantityCache>(true);
240
241 // since thermodynamic hints are basically free if the cache for intensive quantities is
242 // enabled, and this model usually shows quite a performance improvment if they are
243 // enabled, let's enable them by default.
244 Parameters::SetDefault<Parameters::EnableThermodynamicHints>(true);
245 }
246
250 std::string primaryVarName(unsigned pvIdx) const
251 {
252 const std::string& tmp = EnergyModule::primaryVarName(pvIdx);
253 if (!tmp.empty())
254 return tmp;
255
256 std::ostringstream oss;
257 if (Indices::z0Idx <= pvIdx && pvIdx < Indices::z0Idx + numComponents - 1)
258 oss << "z_," << FluidSystem::componentName(/*compIdx=*/pvIdx - Indices::z0Idx);
259 else if (pvIdx==Indices::pressure0Idx)
260 oss << "pressure_" << FluidSystem::phaseName(0);
261 else
262 assert(false);
263
264 return oss.str();
265 }
266
270 std::string eqName(unsigned eqIdx) const
271 {
272 const std::string& tmp = EnergyModule::eqName(eqIdx);
273 if (!tmp.empty())
274 return tmp;
275
276 std::ostringstream oss;
277 if (Indices::conti0EqIdx <= eqIdx && eqIdx < Indices::conti0EqIdx
278 + numComponents) {
279 unsigned compIdx = eqIdx - Indices::conti0EqIdx;
280 oss << "continuity^" << FluidSystem::componentName(compIdx);
281 }
282 else
283 assert(false);
284
285 return oss.str();
286 }
287
288 void registerOutputModules_()
289 {
290 ParentType::registerOutputModules_();
291
292 // add the VTK output modules which are meaningful for the model
293 this->addOutputModule(new Opm::VtkCompositionModule<TypeTag>(this->simulator_));
294 this->addOutputModule(new Opm::VtkPTFlashModule<TypeTag>(this->simulator_));
295 if (enableDiffusion)
296 this->addOutputModule(new Opm::VtkDiffusionModule<TypeTag>(this->simulator_));
297 if (enableEnergy)
298 this->addOutputModule(new Opm::VtkEnergyModule<TypeTag>(this->simulator_));
299 }
300};
301
302} // namespace Opm
303
304#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 flashindices.hh:47
Contains the intensive quantities of the flash-based compositional multi-phase model.
Definition flashintensivequantities.hh:58
Calculates the local residual of the compositional multi-phase model based on flash calculations.
Definition flashlocalresidual.hh:46
A compositional multi-phase model based on flash-calculations.
Definition flashmodel.hh:191
static void registerParameters()
Register all run-time parameters for the immiscible model.
Definition flashmodel.hh:215
std::string primaryVarName(unsigned pvIdx) const
Given an primary variable index, return a human readable name.
Definition flashmodel.hh:250
std::string eqName(unsigned eqIdx) const
Given an equation index, return a human readable name.
Definition flashmodel.hh:270
A Newton solver specific to the PTFlash model.
Definition flashnewtonmethod.hh:53
Represents the primary variables used by the compositional flow model based on flash calculations.
Definition flashprimaryvariables.hh:51
Implements a vector representing rates of conserved quantities.
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
VTK output module for the fluid composition.
Definition vtkcompositionmodule.hpp:57
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition vtkcompositionmodule.hpp:86
VTK output module for quantities which make sense for models which incorperate molecular diffusion.
Definition vtkdiffusionmodule.hpp:58
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition vtkdiffusionmodule.hpp:87
VTK output module for quantities which make sense for models which assume thermal equilibrium.
Definition vtkenergymodule.hpp:58
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition vtkenergymodule.hpp:86
VTK output module for the PT Flash calculation This module deals with the following quantities: K,...
Definition vtkptflashmodule.hpp:53
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition vtkptflashmodule.hpp:82
Contains the classes required to consider energy as a conservation quantity in a multi-phase module.
Implements a boundary vector for the fully implicit compositional multi-phase model which is based on...
This template class contains the data which is required to calculate all fluxes of components over a ...
A Newton solver specific to the PTFlash model.
Declares the properties required by the compositional multi-phase model based on flash calculations.
Implements a vector representing rates of conserved quantities.
A base class for fully-implicit multi-phase porous-media flow models which assume multiple fluid phas...
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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
Defines the primary variable and equation indices for the compositional multi-phase model based on fl...
Contains the intensive quantities of the flash-based compositional multi-phase model.
Calculates the local residual of the compositional multi-phase model based on flash calculations.
Declares the parameters for the compositional multi-phase model based on flash calculations.
Represents the primary variables used by the compositional flow model based on flash calculations.
Definition multiphasebaseproperties.hh:87
The type of the flash constraint solver.
Definition flashproperties.hh:39
Specifies the type of the actual Newton method.
Definition newtonmethodproperties.hh:32
VTK output module for the fluid composition.
VTK output module for quantities which make sense for models which incorperate molecular diffusion.
VTK output module for quantities which make sense for models which assume thermal equilibrium.
VTK output module for the PT Flash calculation This module deals with the following quantities: K,...