28 #ifndef EWOMS_FLASH_MODEL_HH
29 #define EWOMS_FLASH_MODEL_HH
31 #include <opm/material/densead/Math.hpp>
47 #include <opm/material/fluidmatrixinteractions/NullMaterial.hpp>
48 #include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
49 #include <opm/material/constraintsolvers/NcpFlash.hpp>
55 template <
class TypeTag>
60 namespace Properties {
174 template <
class TypeTag>
176 :
public MultiPhaseBaseModel<TypeTag>
178 typedef MultiPhaseBaseModel<TypeTag> ParentType;
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;
184 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
187 enum { enableDiffusion =
GET_PROP_VALUE(TypeTag, EnableDiffusion) };
194 FlashModel(Simulator& simulator)
195 : ParentType(simulator)
215 "The maximum tolerance for the flash solver to "
216 "consider the solution converged");
230 const std::string& tmp = EnergyModule::primaryVarName(pvIdx);
234 std::ostringstream oss;
235 if (Indices::cTot0Idx <= pvIdx && pvIdx < Indices::cTot0Idx
237 oss <<
"c_tot," << FluidSystem::componentName(pvIdx
238 - Indices::cTot0Idx);
250 const std::string& tmp = EnergyModule::eqName(eqIdx);
254 std::ostringstream oss;
255 if (Indices::conti0EqIdx <= eqIdx && eqIdx < Indices::conti0EqIdx
257 unsigned compIdx = eqIdx - Indices::conti0EqIdx;
258 oss <<
"continuity^" << FluidSystem::componentName(compIdx);
271 Scalar tmp = EnergyModule::primaryVarWeight(*
this, globalDofIdx, pvIdx);
275 unsigned compIdx = pvIdx - Indices::cTot0Idx;
282 return FluidSystem::molarMass(compIdx) / 100.0;
288 Scalar
eqWeight(
unsigned globalDofIdx,
unsigned eqIdx)
const
290 Scalar tmp = EnergyModule::eqWeight(*
this, globalDofIdx, eqIdx);
294 unsigned compIdx = eqIdx - Indices::conti0EqIdx;
297 return FluidSystem::molarMass(compIdx);
300 void registerOutputModules_()
302 ParentType::registerOutputModules_();
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