27 #ifndef OPM_COMPUTE_FROM_REFERENCE_PHASE_HPP
28 #define OPM_COMPUTE_FROM_REFERENCE_PHASE_HPP
34 #include <dune/common/fvector.hh>
63 template <
class Scalar,
class Flu
idSystem,
class Evaluation = Scalar>
66 enum { numPhases = FluidSystem::numPhases };
67 enum { numComponents = FluidSystem::numComponents };
69 typedef Dune::FieldVector<Evaluation, numComponents> ComponentVector;
106 template <
class Flu
idState>
107 static void solve(FluidState& fluidState,
108 typename FluidSystem::template ParameterCache<typename FluidState::Scalar>& paramCache,
109 unsigned refPhaseIdx,
115 paramCache.updatePhase(fluidState, refPhaseIdx);
116 fluidState.setDensity(refPhaseIdx,
117 FluidSystem::density(fluidState,
122 fluidState.setEnthalpy(refPhaseIdx,
123 FluidSystem::enthalpy(fluidState,
128 fluidState.setViscosity(refPhaseIdx,
129 FluidSystem::viscosity(fluidState,
134 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
135 fluidState.setFugacityCoefficient(refPhaseIdx,
137 FluidSystem::fugacityCoefficient(fluidState,
144 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
145 if (phaseIdx == refPhaseIdx)
148 ComponentVector fugVec;
149 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
150 const auto& fug = fluidState.fugacity(refPhaseIdx, compIdx);
151 fugVec[compIdx] = decay<Evaluation>(fug);
157 fluidState.setViscosity(phaseIdx,
158 FluidSystem::viscosity(fluidState,
163 fluidState.setEnthalpy(phaseIdx,
164 FluidSystem::enthalpy(fluidState,
Calculates the chemical equilibrium from the component fugacities in a phase.
Some templates to wrap the valgrind client request macros.
Calculates the chemical equilibrium from the component fugacities in a phase.
Definition: CompositionFromFugacities.hpp:47
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:80
Computes all quantities of a generic fluid state if a reference phase has been specified.
Definition: ComputeFromReferencePhase.hpp:65
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:107