00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00027 #ifndef OPM_COMPUTE_FROM_REFERENCE_PHASE_HPP
00028 #define OPM_COMPUTE_FROM_REFERENCE_PHASE_HPP
00029
00030 #include <opm/material/constraintsolvers/CompositionFromFugacities.hpp>
00031
00032 #include <opm/common/Exceptions.hpp>
00033 #include <opm/common/ErrorMacros.hpp>
00034 #include <opm/common/Valgrind.hpp>
00035
00036 #include <dune/common/fvector.hh>
00037
00038 namespace Opm {
00039
00065 template <class Scalar, class FluidSystem, class Evaluation = Scalar>
00066 class ComputeFromReferencePhase
00067 {
00068 enum { numPhases = FluidSystem::numPhases };
00069 enum { numComponents = FluidSystem::numComponents };
00070 typedef Opm::CompositionFromFugacities<Scalar, FluidSystem, Evaluation> CompositionFromFugacities;
00071 typedef Dune::FieldVector<Evaluation, numComponents> ComponentVector;
00072
00073 public:
00108 template <class FluidState>
00109 static void solve(FluidState& fluidState,
00110 typename FluidSystem::template ParameterCache<typename FluidState::Scalar>& paramCache,
00111 unsigned refPhaseIdx,
00112 bool setViscosity,
00113 bool setEnthalpy)
00114 {
00115
00116
00117 paramCache.updatePhase(fluidState, refPhaseIdx);
00118 fluidState.setDensity(refPhaseIdx,
00119 FluidSystem::density(fluidState,
00120 paramCache,
00121 refPhaseIdx));
00122
00123 if (setEnthalpy)
00124 fluidState.setEnthalpy(refPhaseIdx,
00125 FluidSystem::enthalpy(fluidState,
00126 paramCache,
00127 refPhaseIdx));
00128
00129 if (setViscosity)
00130 fluidState.setViscosity(refPhaseIdx,
00131 FluidSystem::viscosity(fluidState,
00132 paramCache,
00133 refPhaseIdx));
00134
00135
00136 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
00137 fluidState.setFugacityCoefficient(refPhaseIdx,
00138 compIdx,
00139 FluidSystem::fugacityCoefficient(fluidState,
00140 paramCache,
00141 refPhaseIdx,
00142 compIdx));
00143 }
00144
00145
00146 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
00147 if (phaseIdx == refPhaseIdx)
00148 continue;
00149
00150 ComponentVector fugVec;
00151 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
00152 const auto& fug = fluidState.fugacity(refPhaseIdx, compIdx);
00153 fugVec[compIdx] = Opm::decay<Evaluation>(fug);
00154 }
00155
00156 CompositionFromFugacities::solve(fluidState, paramCache, phaseIdx, fugVec);
00157
00158 if (setViscosity)
00159 fluidState.setViscosity(phaseIdx,
00160 FluidSystem::viscosity(fluidState,
00161 paramCache,
00162 phaseIdx));
00163
00164 if (setEnthalpy)
00165 fluidState.setEnthalpy(phaseIdx,
00166 FluidSystem::enthalpy(fluidState,
00167 paramCache,
00168 phaseIdx));
00169 }
00170 }
00171 };
00172
00173 }
00174
00175 #endif