27 #ifndef OPM_COMPUTE_FROM_REFERENCE_PHASE_HPP
28 #define OPM_COMPUTE_FROM_REFERENCE_PHASE_HPP
32 #include <opm/common/Exceptions.hpp>
33 #include <opm/common/ErrorMacros.hpp>
34 #include <opm/common/Valgrind.hpp>
36 #include <dune/common/fvector.hh>
65 template <
class Scalar,
class Flu
idSystem,
class Evaluation = Scalar>
68 enum { numPhases = FluidSystem::numPhases };
69 enum { numComponents = FluidSystem::numComponents };
71 typedef Dune::FieldVector<Evaluation, numComponents> ComponentVector;
108 template <
class Flu
idState>
109 static void solve(FluidState& fluidState,
110 typename FluidSystem::template ParameterCache<typename FluidState::Scalar>& paramCache,
111 unsigned refPhaseIdx,
117 paramCache.updatePhase(fluidState, refPhaseIdx);
118 fluidState.setDensity(refPhaseIdx,
119 FluidSystem::density(fluidState,
124 fluidState.setEnthalpy(refPhaseIdx,
125 FluidSystem::enthalpy(fluidState,
130 fluidState.setViscosity(refPhaseIdx,
131 FluidSystem::viscosity(fluidState,
136 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
137 fluidState.setFugacityCoefficient(refPhaseIdx,
139 FluidSystem::fugacityCoefficient(fluidState,
146 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
147 if (phaseIdx == refPhaseIdx)
150 ComponentVector fugVec;
151 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
152 const auto& fug = fluidState.fugacity(refPhaseIdx, compIdx);
153 fugVec[compIdx] = Opm::decay<Evaluation>(fug);
159 fluidState.setViscosity(phaseIdx,
160 FluidSystem::viscosity(fluidState,
165 fluidState.setEnthalpy(phaseIdx,
166 FluidSystem::enthalpy(fluidState,
static void solve(FluidState &fluidState, typename FluidSystem::template ParameterCache< typename FluidState::Scalar > ¶mCache, unsigned phaseIdx, const ComponentVector &targetFug)
Calculates the chemical equilibrium from the component fugacities in a phase.
Definition: CompositionFromFugacities.hpp:88
Computes all quantities of a generic fluid state if a reference phase has been specified.
Definition: ComputeFromReferencePhase.hpp:66
static void solve(FluidState &fluidState, typename FluidSystem::template ParameterCache< typename FluidState::Scalar > ¶mCache, unsigned refPhaseIdx, bool setViscosity, bool setEnthalpy)
Computes all quantities of a generic fluid state if a reference phase has been specified.
Definition: ComputeFromReferencePhase.hpp:109
Calculates the chemical equilibrium from the component fugacities in a phase.
Calculates the chemical equilibrium from the component fugacities in a phase.
Definition: CompositionFromFugacities.hpp:54