28 #ifndef EWOMS_IMMISCIBLE_RATE_VECTOR_HH 29 #define EWOMS_IMMISCIBLE_RATE_VECTOR_HH 31 #include <dune/common/fvector.hh> 33 #include <opm/common/Valgrind.hpp> 34 #include <opm/material/constraintsolvers/NcpFlash.hpp> 47 template <
class TypeTag>
49 :
public Dune::FieldVector<typename GET_PROP_TYPE(TypeTag, Evaluation),
50 GET_PROP_VALUE(TypeTag, NumEq)>
52 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
53 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
54 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
55 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
58 enum { conti0EqIdx = Indices::conti0EqIdx };
62 typedef Dune::FieldVector<Evaluation, numEq> ParentType;
70 { Opm::Valgrind::SetUndefined(*
this); }
100 { ParentType::operator=(value); }
114 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
115 ParentType::operator[](conti0EqIdx + compIdx) =
116 value[conti0EqIdx + compIdx]*FluidSystem::molarMass(compIdx);
126 template <
class RhsEval>
128 { EnergyModule::setEnthalpyRate(*
this, rate); }
145 template <
class Flu
idState,
class RhsEval>
148 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
149 (*
this)[conti0EqIdx + compIdx] =
150 fluidState.density(phaseIdx, compIdx)
151 * fluidState.massFraction(phaseIdx, compIdx)
154 EnergyModule::setEnthalpyRate(*
this, fluidState, phaseIdx, volume);
160 template <
class RhsEval>
163 for (
unsigned i=0; i < this->size(); ++i)
173 for (
unsigned i=0; i < this->size(); ++i)
174 (*
this)[i] = other[i];
Contains the quantities which are are constant within a finite volume for the immiscible multi-phase ...
ImmiscibleRateVector(const ImmiscibleRateVector &value)
Copy constructor.
Definition: immiscibleratevector.hh:86
void setMassRate(const ParentType &value)
Set a mass rate of the conservation quantities.
Definition: immiscibleratevector.hh:99
Definition: baseauxiliarymodule.hh:37
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
Implements a vector representing rates of conserved quantities.
Definition: immiscibleratevector.hh:48
ImmiscibleRateVector(const Evaluation &value)
Constructor with assignment from scalar.
Definition: immiscibleratevector.hh:77
void setVolumetricRate(const FluidState &fluidState, unsigned phaseIdx, const RhsEval &volume)
Set a volumetric rate of a phase.
Definition: immiscibleratevector.hh:146
ImmiscibleRateVector()
Default constructor.
Definition: immiscibleratevector.hh:69
void setEnthalpyRate(const RhsEval &rate)
Set an enthalpy rate [J/As] where .
Definition: immiscibleratevector.hh:127
ImmiscibleRateVector & operator=(const RhsEval &value)
Assignment operator from a scalar or a function evaluation.
Definition: immiscibleratevector.hh:161
ImmiscibleRateVector & operator=(const ImmiscibleRateVector &other)
Assignment operator from another rate vector.
Definition: immiscibleratevector.hh:171
void setMolarRate(const ParentType &value)
Set a molar rate of the conservation quantities.
Definition: immiscibleratevector.hh:111