28 #ifndef EWOMS_IMMISCIBLE_PRIMARY_VARIABLES_HH 29 #define EWOMS_IMMISCIBLE_PRIMARY_VARIABLES_HH 36 #include <opm/material/constraintsolvers/ImmiscibleFlash.hpp> 37 #include <opm/material/fluidstates/ImmiscibleFluidState.hpp> 38 #include <opm/common/Valgrind.hpp> 40 #include <dune/common/fvector.hh> 53 template <
class TypeTag>
58 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
59 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
60 typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) Implementation;
61 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
62 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
63 typedef typename GET_PROP_TYPE(TypeTag, MaterialLawParams) MaterialLawParams;
65 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
68 enum { pressure0Idx = Indices::pressure0Idx };
69 enum { saturation0Idx = Indices::saturation0Idx };
74 typedef typename Opm::MathToolbox<Evaluation> Toolbox;
75 typedef Dune::FieldVector<Scalar, numComponents> ComponentVector;
76 typedef Opm::ImmiscibleFlash<Scalar, FluidSystem> ImmiscibleFlash;
84 { Opm::Valgrind::SetUndefined(*
this); }
128 template <
class Flu
idState>
130 const MaterialLawParams& matParams,
131 bool isInEquilibrium =
false)
135 for (
unsigned phaseIdx = 1; phaseIdx < numPhases; ++phaseIdx) {
136 assert(std::abs(fluidState.temperature(0) - fluidState.temperature(phaseIdx)) < 1e-30);
142 if (isInEquilibrium) {
149 typename FluidSystem::template ParameterCache<Scalar> paramCache;
150 Opm::ImmiscibleFluidState<Scalar, FluidSystem> fsFlash;
154 fsFlash.assign(fluidState);
157 paramCache.updateAll(fsFlash);
158 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
159 Scalar rho = FluidSystem::density(fsFlash, paramCache, phaseIdx);
160 fsFlash.setDensity(phaseIdx, rho);
164 ComponentVector globalMolarities(0.0);
165 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
166 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
167 globalMolarities[compIdx] +=
168 fsFlash.saturation(phaseIdx) * fsFlash.molarity(phaseIdx, compIdx);
173 ImmiscibleFlash::template solve<MaterialLaw>(fsFlash, matParams, paramCache, globalMolarities);
195 template <
class Flu
idState>
200 EnergyModule::setPriVarTemperatures(asImp_(), fluidState);
202 (*this)[pressure0Idx] = fluidState.pressure(0);
203 for (
unsigned phaseIdx = 0; phaseIdx < numPhases - 1; ++phaseIdx)
204 (*
this)[saturation0Idx + phaseIdx] = fluidState.saturation(phaseIdx);
208 Implementation& asImp_()
209 {
return *
static_cast<Implementation *
>(
this); }
void assignNaive(const FluidState &fluidState)
Directly retrieve the primary variables from an arbitrary fluid state.
Definition: immiscibleprimaryvariables.hh:196
Definition: baseauxiliarymodule.hh:37
ImmisciblePrimaryVariables(Scalar value)
Constructor with assignment from scalar.
Definition: immiscibleprimaryvariables.hh:91
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...
ImmisciblePrimaryVariables()
Default constructor.
Definition: immiscibleprimaryvariables.hh:83
ImmisciblePrimaryVariables & operator=(const ImmisciblePrimaryVariables &value)=default
Assignment operator.
void assignMassConservative(const FluidState &fluidState, const MaterialLawParams &matParams, bool isInEquilibrium=false)
Set the primary variables from an arbitrary fluid state in a mass conservative way.
Definition: immiscibleprimaryvariables.hh:129
Represents the primary variables used by the a model.
Definition: fvbaseprimaryvariables.hh:48
Represents the primary variables used by the a model.
Defines the properties required for the immiscible multi-phase model.
Represents the primary variables used by the immiscible multi-phase, model.
Definition: immiscibleprimaryvariables.hh:54