28 #ifndef EWOMS_PVS_PRIMARY_VARIABLES_HH
29 #define EWOMS_PVS_PRIMARY_VARIABLES_HH
37 #include <opm/material/constraintsolvers/NcpFlash.hpp>
38 #include <opm/material/fluidstates/CompositionalFluidState.hpp>
39 #include <opm/common/Valgrind.hpp>
40 #include <opm/common/ErrorMacros.hpp>
41 #include <opm/common/Exceptions.hpp>
43 #include <dune/common/fvector.hh>
58 template <
class TypeTag>
63 typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) Implementation;
65 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
66 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
67 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
68 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
69 typedef typename GET_PROP_TYPE(TypeTag, MaterialLawParams) MaterialLawParams;
70 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
73 enum { pressure0Idx = Indices::pressure0Idx };
74 enum { switch0Idx = Indices::switch0Idx };
80 typedef typename Opm::MathToolbox<Evaluation> Toolbox;
81 typedef Dune::FieldVector<Scalar, numComponents> ComponentVector;
83 typedef Opm::NcpFlash<Scalar, FluidSystem> NcpFlash;
87 { Opm::Valgrind::SetDefined(*
this); }
94 Opm::Valgrind::CheckDefined(value);
95 Opm::Valgrind::SetDefined(*
this);
106 Opm::Valgrind::SetDefined(*
this);
108 phasePresence_ = value.phasePresence_;
114 template <
class Flu
idState>
116 const MaterialLawParams& matParams,
117 bool isInEquilibrium =
false)
121 for (
unsigned phaseIdx = 1; phaseIdx < numPhases; ++phaseIdx) {
122 assert(std::abs(fluidState.temperature(0) - fluidState.temperature(phaseIdx)) < 1e-30);
128 if (isInEquilibrium) {
135 typename FluidSystem::template ParameterCache<Scalar> paramCache;
136 Opm::CompositionalFluidState<Scalar, FluidSystem> fsFlash;
140 fsFlash.assign(fluidState);
143 paramCache.updateAll(fsFlash);
144 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
145 Scalar rho = FluidSystem::density(fsFlash, paramCache, phaseIdx);
146 fsFlash.setDensity(phaseIdx, rho);
149 ComponentVector globalMolarities(0.0);
150 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
151 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
152 globalMolarities[compIdx] +=
153 fsFlash.saturation(phaseIdx) * fsFlash.molarity(phaseIdx, compIdx);
158 NcpFlash::template solve<MaterialLaw>(fsFlash, matParams, paramCache, globalMolarities);
169 {
return phasePresence_; }
178 { phasePresence_ = value; }
211 {
return phasePresence& (1 << phaseIdx); }
220 {
return phasePresence_& (1 << phaseIdx); }
226 {
return static_cast<unsigned>(ffs(phasePresence_) - 1); }
234 phasePresence_ = value.phasePresence_;
263 unsigned varIdx = switch0Idx + phaseIdx - 1;
265 Toolbox::createConstant((*
this)[varIdx]);
266 return Toolbox::createVariable((*
this)[varIdx], varIdx);
272 template <
class Flu
idState>
275 typedef Opm::MathToolbox<typename FluidState::Scalar> FsToolbox;
279 EnergyModule::setPriVarTemperatures(*
this, fluidState);
282 (*this)[pressure0Idx] = FsToolbox::value(fluidState.pressure(0));
283 Opm::Valgrind::CheckDefined((*
this)[pressure0Idx]);
287 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
291 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
292 a -= FsToolbox::value(fluidState.moleFraction(phaseIdx, compIdx));
294 Scalar b = FsToolbox::value(fluidState.saturation(phaseIdx));
297 phasePresence_ |= (1 << phaseIdx);
301 if (phasePresence_ == 0)
302 OPM_THROW(Opm::NumericalProblem,
303 "Phase state was 0, i.e., no fluid is present");
308 for (
unsigned switchIdx = 0; switchIdx < numPhases - 1; ++switchIdx) {
309 unsigned phaseIdx = switchIdx;
310 unsigned compIdx = switchIdx + 1;
311 if (switchIdx >= lowestPhaseIdx)
315 (*this)[switch0Idx + switchIdx] = FsToolbox::value(fluidState.saturation(phaseIdx));
316 Opm::Valgrind::CheckDefined((*
this)[switch0Idx + switchIdx]);
319 (*this)[switch0Idx + switchIdx] =
320 FsToolbox::value(fluidState.moleFraction(lowestPhaseIdx, compIdx));
321 Opm::Valgrind::CheckDefined((*
this)[switch0Idx + switchIdx]);
327 for (
unsigned compIdx = numPhases - 1; compIdx < numComponents - 1; ++compIdx) {
328 (*this)[switch0Idx + compIdx] =
329 FsToolbox::value(fluidState.moleFraction(lowestPhaseIdx, compIdx + 1));
330 Opm::Valgrind::CheckDefined((*
this)[switch0Idx + compIdx]);
337 void print(std::ostream& os = std::cout)
const
339 os <<
"(p_" << FluidSystem::phaseName(0) <<
" = "
340 << this->operator[](pressure0Idx);
342 for (
unsigned switchIdx = 0; switchIdx < numPhases - 1; ++switchIdx) {
343 unsigned phaseIdx = switchIdx;
344 unsigned compIdx = switchIdx + 1;
345 if (phaseIdx >= lowestPhaseIdx)
350 os <<
", S_" << FluidSystem::phaseName(phaseIdx) <<
" = "
351 << (*this)[switch0Idx + switchIdx];
354 os <<
", x_" << FluidSystem::phaseName(lowestPhaseIdx) <<
"^"
355 << FluidSystem::componentName(compIdx) <<
" = "
356 << (*this)[switch0Idx + switchIdx];
359 for (
unsigned compIdx = numPhases - 1; compIdx < numComponents - 1;
361 os <<
", x_" << FluidSystem::phaseName(lowestPhaseIdx) <<
"^"
362 << FluidSystem::componentName(compIdx + 1) <<
" = "
363 << (*this)[switch0Idx + compIdx];
366 os <<
", phase presence: " <<
static_cast<int>(phasePresence_) << std::flush;
370 short phasePresence_;
Represents the primary variables used in the primary variable switching compositional model...
Definition: pvsprimaryvariables.hh:59
void print(std::ostream &os=std::cout) const
Prints the names of the primary variables and their values.
Definition: pvsprimaryvariables.hh:337
bool phaseIsPresent(unsigned phaseIdx) const
Returns true iff a phase is present for the current phase presence.
Definition: pvsprimaryvariables.hh:219
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
short phasePresence() const
Return the fluid phases which are present in a given control volume.
Definition: pvsprimaryvariables.hh:168
Contains the classes required to consider energy as a conservation quantity in a multi-phase module...
unsigned lowestPresentPhaseIdx() const
Returns the phase with the lowest index that is present.
Definition: pvsprimaryvariables.hh:225
static bool phaseIsPresent(unsigned phaseIdx, short phasePresence)
Returns true iff a phase is present for a given phase presence.
Definition: pvsprimaryvariables.hh:210
The indices for the compositional multi-phase primary variable switching model.
ThisType & operator=(Scalar value)
Assignment operator from a scalar value.
Definition: pvsprimaryvariables.hh:242
void setPhasePresence(short value)
Set which fluid phases are present in a given control volume.
Definition: pvsprimaryvariables.hh:177
PvsPrimaryVariables(const PvsPrimaryVariables &value)
Default constructor.
Definition: pvsprimaryvariables.hh:104
Represents the primary variables used by the a model.
Definition: fvbaseprimaryvariables.hh:48
FvBasePrimaryVariables & operator=(const FvBasePrimaryVariables &value)=default
Assignment from another primary variables object.
void setPhasePresent(unsigned phaseIdx, bool yesno=true)
Set whether a given indivividual phase should be present or not.
Definition: pvsprimaryvariables.hh:187
Represents the primary variables used by the a model.
Declares the properties required for the compositional multi-phase primary variable switching model...
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: pvsprimaryvariables.hh:115
Evaluation explicitSaturationValue(unsigned phaseIdx, unsigned timeIdx) const
Returns an explcitly stored saturation for a given phase.
Definition: pvsprimaryvariables.hh:257
ThisType & operator=(const Implementation &value)
Assignment operator from an other primary variables object.
Definition: pvsprimaryvariables.hh:231
PvsPrimaryVariables(Scalar value)
Constructor with assignment from scalar.
Definition: pvsprimaryvariables.hh:92
unsigned implicitSaturationIdx() const
Returns the index of the phase with's its saturation is determined by the closure condition of satura...
Definition: pvsprimaryvariables.hh:199
void assignNaive(const FluidState &fluidState)
Directly retrieve the primary variables from an arbitrary fluid state.
Definition: pvsprimaryvariables.hh:273