28#ifndef EWOMS_FLASH_MODEL_HH
29#define EWOMS_FLASH_MODEL_HH
31#include <opm/material/constraintsolvers/NcpFlash.hpp>
33#include <opm/material/densead/Math.hpp>
35#include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
36#include <opm/material/fluidmatrixinteractions/NullMaterial.hpp>
59template <
class TypeTag>
63namespace Opm::Properties {
67struct FlashModel {
using InheritsFrom = std::tuple<MultiPhaseBaseModel>; };
71template<
class TypeTag>
75template<
class TypeTag>
76struct FlashSolver<TypeTag, TTag::FlashModel>
77{
using type = Opm::NcpFlash<GetPropType<TypeTag, Properties::Scalar>,
78 GetPropType<TypeTag, Properties::FluidSystem>>; };
81template<
class TypeTag>
85template<
class TypeTag>
89template<
class TypeTag>
93template<
class TypeTag>
97template<
class TypeTag>
101template<
class TypeTag>
105template<
class TypeTag>
106struct Indices<TypeTag, TTag::FlashModel> {
using type =
Opm::FlashIndices<TypeTag, 0>; };
109template<
class TypeTag>
110struct EnableDiffusion<TypeTag, TTag::FlashModel> {
static constexpr bool value =
false; };
113template<
class TypeTag>
114struct EnableEnergy<TypeTag, TTag::FlashModel> {
static constexpr bool value =
false; };
179template <
class TypeTag>
181 :
public MultiPhaseBaseModel<TypeTag>
183 using ParentType = MultiPhaseBaseModel<TypeTag>;
185 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
186 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
187 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
189 using Indices = GetPropType<TypeTag, Properties::Indices>;
199 FlashModel(Simulator& simulator)
200 : ParentType(simulator)
219 Parameters::Register<Parameters::FlashTolerance<Scalar>>
220 (
"The maximum tolerance for the flash solver to "
221 "consider the solution converged");
225 Parameters::SetDefault<Parameters::EnableIntensiveQuantityCache>(
true);
230 Parameters::SetDefault<Parameters::EnableThermodynamicHints>(
true);
244 const std::string& tmp = EnergyModule::primaryVarName(pvIdx);
248 std::ostringstream oss;
249 if (Indices::cTot0Idx <= pvIdx && pvIdx < Indices::cTot0Idx
251 oss <<
"c_tot," << FluidSystem::componentName(pvIdx
252 - Indices::cTot0Idx);
264 const std::string& tmp = EnergyModule::eqName(eqIdx);
268 std::ostringstream oss;
269 if (Indices::conti0EqIdx <= eqIdx && eqIdx < Indices::conti0EqIdx
271 unsigned compIdx = eqIdx - Indices::conti0EqIdx;
272 oss <<
"continuity^" << FluidSystem::componentName(compIdx);
285 Scalar tmp = EnergyModule::primaryVarWeight(*
this, globalDofIdx, pvIdx);
289 unsigned compIdx = pvIdx - Indices::cTot0Idx;
296 return FluidSystem::molarMass(compIdx) / 100.0;
302 Scalar
eqWeight(
unsigned globalDofIdx,
unsigned eqIdx)
const
304 Scalar tmp = EnergyModule::eqWeight(*
this, globalDofIdx, eqIdx);
308 unsigned compIdx = eqIdx - Indices::conti0EqIdx;
311 return FluidSystem::molarMass(compIdx);
314 void registerOutputModules_()
316 ParentType::registerOutputModules_();
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:48
Contains the intensive quantities of the flash-based compositional multi-phase model.
Definition flashintensivequantities.hh:60
Calculates the local residual of the compositional multi-phase model based on flash calculations.
Definition flashlocalresidual.hh:45
A compositional multi-phase model based on flash-calculations.
Definition flashmodel.hh:187
static void registerParameters()
Register all run-time parameters for the immiscible model.
Definition flashmodel.hh:206
Scalar primaryVarWeight(unsigned globalDofIdx, unsigned pvIdx) const
Returns the relative weight of a primary variable for calculating relative errors.
Definition flashmodel.hh:283
static std::string name()
Definition flashmodel.hh:236
Scalar eqWeight(unsigned globalDofIdx, unsigned eqIdx) const
Returns the relative weight of an equation.
Definition flashmodel.hh:302
std::string primaryVarName(unsigned pvIdx) const
Given an primary variable index, return a human readable name.
Definition flashmodel.hh:242
std::string eqName(unsigned eqIdx) const
Given an equation index, return a human readable name.
Definition flashmodel.hh:262
Represents the primary variables used by the compositional flow model based on flash calculations.
Definition flashprimaryvariables.hh:57
Implements a vector representing rates of conserved quantities.
Definition flashratevector.hh:48
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
Contains the classes required to consider energy as a conservation quantity in a multi-phase module.
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.
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 ...
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
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.