From 662d6e28cdbe7795d0dc8a3d8d31024a38d9a01b Mon Sep 17 00:00:00 2001 From: Andreas Lauser Date: Wed, 28 Sep 2016 16:38:33 +0200 Subject: [PATCH] threshold pressures: properly consider the dissolved components for the density even if a phase is saturated, the dissolved component must be considered when calculating the fluid density. --- opm/simulators/thresholdPressures.hpp | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/opm/simulators/thresholdPressures.hpp b/opm/simulators/thresholdPressures.hpp index f3007c69cc..066bfa5c78 100644 --- a/opm/simulators/thresholdPressures.hpp +++ b/opm/simulators/thresholdPressures.hpp @@ -208,6 +208,10 @@ void computeMaxDp(std::map, double>& maxDp, if (Rs >= RsSat) { double b = pvto.saturatedInverseFormationVolumeFactor(pvtRegionIdx, T, p); rho[opos][cellIdx] = surfaceDensity[pvtRegionIdx][opos]*b; + if (pu.phase_used[BlackoilPhases::Vapour]) { + int gpos = pu.phase_pos[BlackoilPhases::Vapour]; + rho[opos][cellIdx] += surfaceDensity[pvtRegionIdx][gpos]*RsSat*b; + } } else { double b = pvto.inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rs); @@ -234,11 +238,14 @@ void computeMaxDp(std::map, double>& maxDp, if (Rv >= RvSat) { double b = pvtg.saturatedInverseFormationVolumeFactor(pvtRegionIdx, T, p); rho[gpos][cellIdx] = surfaceDensity[pvtRegionIdx][gpos]*b; + if (pu.phase_used[BlackoilPhases::Liquid]) { + int opos = pu.phase_pos[BlackoilPhases::Liquid]; + rho[gpos][cellIdx] += surfaceDensity[pvtRegionIdx][opos]*RvSat*b; + } } else { double b = pvtg.inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rv); rho[gpos][cellIdx] = surfaceDensity[pvtRegionIdx][gpos]*b; - if (pu.phase_used[BlackoilPhases::Liquid]) { int opos = pu.phase_pos[BlackoilPhases::Liquid]; rho[gpos][cellIdx] += surfaceDensity[pvtRegionIdx][opos]*Rv*b;