From 4aa0eaff6755af898dbc956bde73f91a34a50584 Mon Sep 17 00:00:00 2001 From: Tor Harald Sandve Date: Thu, 28 Nov 2013 11:27:25 +0100 Subject: [PATCH 1/2] Whether the fluid is saturated or not is explicitly passed to the pvts The criteria for whether the fluid is saturated or not is moved from the within the pvt calculations to the solver, and passed to the pvt calculations as a array of boolean values. --- opm/autodiff/BlackoilPropsAd.cpp | 10 +++- opm/autodiff/BlackoilPropsAd.hpp | 8 ++++ opm/autodiff/BlackoilPropsAdFromDeck.cpp | 15 ++++-- opm/autodiff/BlackoilPropsAdFromDeck.hpp | 8 ++++ opm/autodiff/BlackoilPropsAdInterface.hpp | 8 ++++ opm/autodiff/FullyImplicitBlackoilSolver.cpp | 50 +++++++++++++++----- opm/autodiff/FullyImplicitBlackoilSolver.hpp | 6 +++ opm/autodiff/ImpesTPFAAD.cpp | 12 +++-- 8 files changed, 96 insertions(+), 21 deletions(-) diff --git a/opm/autodiff/BlackoilPropsAd.cpp b/opm/autodiff/BlackoilPropsAd.cpp index a4ac552da..9c4827e25 100644 --- a/opm/autodiff/BlackoilPropsAd.cpp +++ b/opm/autodiff/BlackoilPropsAd.cpp @@ -121,10 +121,12 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. V BlackoilPropsAd::muOil(const V& po, const V& rs, + const bool* /*isSat*/, const Cells& cells) const { if (!pu_.phase_used[Oil]) { @@ -197,14 +199,16 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. ADB BlackoilPropsAd::muOil(const ADB& po, const ADB& rs, + const bool* isSat, const Cells& cells) const { #if 1 - return ADB::constant(muOil(po.value(), rs.value(), cells), po.blockPattern()); + return ADB::constant(muOil(po.value(), rs.value(), isSat,cells), po.blockPattern()); #else if (!pu_.phase_used[Oil]) { OPM_THROW(std::runtime_error, "Cannot call muOil(): oil phase not present."); @@ -306,10 +310,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. V BlackoilPropsAd::bOil(const V& po, const V& rs, + const bool* /*isSat*/, const Cells& cells) const { if (!pu_.phase_used[Oil]) { @@ -382,10 +388,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. ADB BlackoilPropsAd::bOil(const ADB& po, const ADB& rs, + const bool* /*isSat*/, const Cells& cells) const { if (!pu_.phase_used[Oil]) { diff --git a/opm/autodiff/BlackoilPropsAd.hpp b/opm/autodiff/BlackoilPropsAd.hpp index 70545d4f6..07e4559b7 100644 --- a/opm/autodiff/BlackoilPropsAd.hpp +++ b/opm/autodiff/BlackoilPropsAd.hpp @@ -109,10 +109,12 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. V muOil(const V& po, const V& rs, + const bool* isSat, const Cells& cells) const; /// Gas viscosity. @@ -132,10 +134,12 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. ADB muOil(const ADB& po, const ADB& rs, + const bool* isSat, const Cells& cells) const; /// Gas viscosity. @@ -158,10 +162,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. V bOil(const V& po, const V& rs, + const bool* isSat, const Cells& cells) const; /// Gas formation volume factor. @@ -181,10 +187,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. ADB bOil(const ADB& po, const ADB& rs, + const bool* isSat, const Cells& cells) const; /// Gas formation volume factor. diff --git a/opm/autodiff/BlackoilPropsAdFromDeck.cpp b/opm/autodiff/BlackoilPropsAdFromDeck.cpp index f3b8c0259..f86d958eb 100644 --- a/opm/autodiff/BlackoilPropsAdFromDeck.cpp +++ b/opm/autodiff/BlackoilPropsAdFromDeck.cpp @@ -212,10 +212,12 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. V BlackoilPropsAdFromDeck::muOil(const V& po, const V& rs, + const bool* isSat, const Cells& cells) const { if (!phase_usage_.phase_used[Oil]) { @@ -227,7 +229,7 @@ namespace Opm V dmudp(n); V dmudr(n); - props_[phase_usage_.phase_pos[Oil]]->mu(n, po.data(), rs.data(), + props_[phase_usage_.phase_pos[Oil]]->mu(n, po.data(), rs.data(),isSat, mu.data(), dmudp.data(), dmudr.data()); return mu; } @@ -285,10 +287,12 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. ADB BlackoilPropsAdFromDeck::muOil(const ADB& po, const ADB& rs, + const bool* isSat, const Cells& cells) const { if (!phase_usage_.phase_used[Oil]) { @@ -300,7 +304,7 @@ namespace Opm V dmudp(n); V dmudr(n); - props_[phase_usage_.phase_pos[Oil]]->mu(n, po.value().data(), rs.value().data(), + props_[phase_usage_.phase_pos[Oil]]->mu(n, po.value().data(), rs.value().data(), isSat, mu.data(), dmudp.data(), dmudr.data()); ADB::M dmudp_diag = spdiag(dmudp); @@ -391,6 +395,7 @@ namespace Opm /// \return Array of n formation volume factor values. V BlackoilPropsAdFromDeck::bOil(const V& po, const V& rs, + const bool* isSat, const Cells& cells) const { if (!phase_usage_.phase_used[Oil]) { @@ -403,7 +408,7 @@ namespace Opm V dbdp(n); V dbdr(n); - props_[phase_usage_.phase_pos[Oil]]->b(n, po.data(), rs.data(), + props_[phase_usage_.phase_pos[Oil]]->b(n, po.data(), rs.data(),isSat, b.data(), dbdp.data(), dbdr.data()); return b; @@ -466,10 +471,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. ADB BlackoilPropsAdFromDeck::bOil(const ADB& po, const ADB& rs, + const bool* isSat, const Cells& cells) const { if (!phase_usage_.phase_used[Oil]) { @@ -482,7 +489,7 @@ namespace Opm V dbdp(n); V dbdr(n); - props_[phase_usage_.phase_pos[Oil]]->b(n, po.value().data(), rs.value().data(), + props_[phase_usage_.phase_pos[Oil]]->b(n, po.value().data(), rs.value().data(),isSat, b.data(), dbdp.data(), dbdr.data()); ADB::M dbdp_diag = spdiag(dbdp); diff --git a/opm/autodiff/BlackoilPropsAdFromDeck.hpp b/opm/autodiff/BlackoilPropsAdFromDeck.hpp index f8cf5f7ad..3f4a8f4e5 100644 --- a/opm/autodiff/BlackoilPropsAdFromDeck.hpp +++ b/opm/autodiff/BlackoilPropsAdFromDeck.hpp @@ -110,10 +110,12 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. V muOil(const V& po, const V& rs, + const bool* isSat, const Cells& cells) const; /// Gas viscosity. @@ -133,10 +135,12 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. ADB muOil(const ADB& po, const ADB& rs, + const bool* isSat, const Cells& cells) const; /// Gas viscosity. @@ -159,10 +163,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. V bOil(const V& po, const V& rs, + const bool* isSat, const Cells& cells) const; /// Gas formation volume factor. @@ -182,10 +188,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. ADB bOil(const ADB& po, const ADB& rs, + const bool* isSat, const Cells& cells) const; /// Gas formation volume factor. diff --git a/opm/autodiff/BlackoilPropsAdInterface.hpp b/opm/autodiff/BlackoilPropsAdInterface.hpp index 98e63fa77..6b4928c48 100644 --- a/opm/autodiff/BlackoilPropsAdInterface.hpp +++ b/opm/autodiff/BlackoilPropsAdInterface.hpp @@ -100,11 +100,13 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. virtual V muOil(const V& po, const V& rs, + const bool* isSat, const Cells& cells) const = 0; /// Gas viscosity. @@ -126,11 +128,13 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. virtual ADB muOil(const ADB& po, const ADB& rs, + const bool* isSat, const Cells& cells) const = 0; /// Gas viscosity. @@ -155,11 +159,13 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. virtual V bOil(const V& po, const V& rs, + const bool* isSat, const Cells& cells) const = 0; /// Gas formation volume factor. @@ -181,11 +187,13 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. virtual ADB bOil(const ADB& po, const ADB& rs, + const bool* isSat, const Cells& cells) const = 0; /// Gas formation volume factor. diff --git a/opm/autodiff/FullyImplicitBlackoilSolver.cpp b/opm/autodiff/FullyImplicitBlackoilSolver.cpp index 096befcf1..2db07478c 100644 --- a/opm/autodiff/FullyImplicitBlackoilSolver.cpp +++ b/opm/autodiff/FullyImplicitBlackoilSolver.cpp @@ -71,8 +71,6 @@ namespace { return all_cells; } - - template AutoDiffBlock::M gravityOperator(const UnstructuredGrid& grid, @@ -530,13 +528,16 @@ namespace { const std::vector& sat = state.saturation; const ADB& rs = state.rs; + bool isSat[rs.size()]; + getSaturatedCells(state,&isSat[0]); + const ADB pv_mult = poroMult(press); const int maxnp = Opm::BlackoilPhases::MaxNumPhases; for (int phase = 0; phase < maxnp; ++phase) { if (active_[ phase ]) { const int pos = pu.phase_pos[ phase ]; - rq_[pos].b = fluidReciprocFVF(phase, press, rs, cells_); + rq_[pos].b = fluidReciprocFVF(phase, press, rs, &isSat[0], cells_); rq_[pos].accum[aix] = pv_mult * rq_[pos].b * sat[pos]; // DUMP(rq_[pos].b); // DUMP(rq_[pos].accum[aix]); @@ -594,6 +595,8 @@ namespace { // DUMP(residual_.mass_balance[phase]); } + bool isSat[grid_.number_of_cells]; + getSaturatedCells(state,isSat); // -------- Extra (optional) sg or rs equation, and rs contributions to the mass balance equations -------- // Add the extra (flux) terms to the gas mass balance equations @@ -613,8 +616,9 @@ namespace { const ADB sg_eq = state.saturation[pg]; const ADB rs_max = fluidRsMax(state.pressure, cells_); const ADB rs_eq = state.rs - rs_max; - Selector use_rs_eq(rs_eq.value()); - residual_.rs_or_sg_eq = use_rs_eq.select(rs_eq, sg_eq); + // Consider the fluid to be saturated if sg >= 1e-14 (a small number) + Selector use_sat_eq(sg_eq.value()-1e-14); + residual_.rs_or_sg_eq = use_sat_eq.select(rs_eq, sg_eq); // DUMP(residual_.rs_or_sg_eq); } @@ -654,7 +658,7 @@ namespace { for (int phase = 0; phase < 3; ++phase) { if (active_[phase]) { const int pos = pu.phase_pos[phase]; - const ADB cell_rho = fluidDensity(phase, state.pressure, state.rs, cells_); + const ADB cell_rho = fluidDensity(phase, state.pressure, state.rs, &isSat[0],cells_); cell_rho_total += state.saturation[pos] * cell_rho; } } @@ -664,7 +668,7 @@ namespace { for (int phase = 0; phase < 3; ++phase) { if (active_[phase]) { const int pos = pu.phase_pos[phase]; - const ADB cell_rho = fluidDensity(phase, state.pressure, state.rs, cells_); + const ADB cell_rho = fluidDensity(phase, state.pressure, state.rs, &isSat[0],cells_); const V fraction = compi.col(pos); inj_rho_total += (wops_.w2p * fraction.matrix()).array() * subset(cell_rho, well_cells); } @@ -885,6 +889,7 @@ namespace { // rs(sg > 0) = rs_sat(sg > 0); // rs(rs > rs_sat*rs_adjust) = rs_sat(rs > rs_sat*rs_adjust); for (int c = 0; c < nc; ++c) { + if (zerosg[c]) { sg[c] = 0.0; } @@ -1002,12 +1007,14 @@ namespace { const SolutionState& state ) { const int phase = canph_[ actph ]; - const ADB mu = fluidViscosity(phase, state.pressure, state.rs, cells_); + bool isSat[grid_.number_of_cells]; + getSaturatedCells(state,&isSat[0]); + const ADB mu = fluidViscosity(phase, state.pressure, state.rs, &isSat[0],cells_); const ADB tr_mult = transMult(state.pressure); rq_[ actph ].mob = tr_mult * kr[ phase ] / mu; - const ADB rho = fluidDensity(phase, state.pressure, state.rs, cells_); + const ADB rho = fluidDensity(phase, state.pressure, state.rs, &isSat[0],cells_); ADB& head = rq_[ actph ].head; @@ -1060,13 +1067,14 @@ namespace { FullyImplicitBlackoilSolver::fluidViscosity(const int phase, const ADB& p , const ADB& rs , + const bool* isSat, const std::vector& cells) const { switch (phase) { case Water: return fluid_.muWat(p, cells); case Oil: { - return fluid_.muOil(p, rs, cells); + return fluid_.muOil(p, rs, isSat,cells); } case Gas: return fluid_.muGas(p, cells); @@ -1083,13 +1091,14 @@ namespace { FullyImplicitBlackoilSolver::fluidReciprocFVF(const int phase, const ADB& p , const ADB& rs , + const bool* isSat, const std::vector& cells) const { switch (phase) { case Water: return fluid_.bWat(p, cells); case Oil: { - return fluid_.bOil(p, rs, cells); + return fluid_.bOil(p, rs, isSat, cells); } case Gas: return fluid_.bGas(p, cells); @@ -1106,10 +1115,11 @@ namespace { FullyImplicitBlackoilSolver::fluidDensity(const int phase, const ADB& p , const ADB& rs , + const bool* isSat, const std::vector& cells) const { const double* rhos = fluid_.surfaceDensity(); - ADB b = fluidReciprocFVF(phase, p, rs, cells); + ADB b = fluidReciprocFVF(phase, p, rs, isSat, cells); ADB rho = V::Constant(p.size(), 1, rhos[phase]) * b; if (phase == Oil && active_[Gas]) { // It is correct to index into rhos with canonical phase indices. @@ -1195,4 +1205,20 @@ namespace { } + void + FullyImplicitBlackoilSolver::getSaturatedCells(const SolutionState& state, bool* isSat) const + { + const int pg = fluid_.phaseUsage().phase_pos[ Gas ]; + const V sg = state.saturation[pg].value(); + for (int c=0; c < sg.size(); ++ c) { + if (sg[c]>0){ + isSat[c] = true; + } + else{ + isSat[c] = false; + } + } + } + + } // namespace Opm diff --git a/opm/autodiff/FullyImplicitBlackoilSolver.hpp b/opm/autodiff/FullyImplicitBlackoilSolver.hpp index 399fb4051..f440c44b0 100644 --- a/opm/autodiff/FullyImplicitBlackoilSolver.hpp +++ b/opm/autodiff/FullyImplicitBlackoilSolver.hpp @@ -189,18 +189,21 @@ namespace Opm { fluidViscosity(const int phase, const ADB& p , const ADB& rs , + const bool* isSat, const std::vector& cells) const; ADB fluidReciprocFVF(const int phase, const ADB& p , const ADB& rs , + const bool* isSat, const std::vector& cells) const; ADB fluidDensity(const int phase, const ADB& p , const ADB& rs , + const bool* isSat, const std::vector& cells) const; V @@ -216,6 +219,9 @@ namespace Opm { ADB transMult(const ADB& p) const; + + void + getSaturatedCells(const SolutionState& state, bool* isSat) const; }; } // namespace Opm diff --git a/opm/autodiff/ImpesTPFAAD.cpp b/opm/autodiff/ImpesTPFAAD.cpp index e7a70292f..abe6367f8 100644 --- a/opm/autodiff/ImpesTPFAAD.cpp +++ b/opm/autodiff/ImpesTPFAAD.cpp @@ -533,7 +533,8 @@ namespace { return fluid_.muWat(p, cells); case Oil: { V dummy_rs = V::Zero(p.size(), 1) * p; - return fluid_.muOil(p, dummy_rs, cells); + bool dummy_isSat[p.size()]; + return fluid_.muOil(p, dummy_rs, dummy_isSat, cells); } case Gas: return fluid_.muGas(p, cells); @@ -553,7 +554,8 @@ namespace { return fluid_.muWat(p, cells); case Oil: { ADB dummy_rs = V::Zero(p.size(), 1) * p; - return fluid_.muOil(p, dummy_rs, cells); + bool dummy_isSat[p.size()]; + return fluid_.muOil(p, dummy_rs, dummy_isSat, cells); } case Gas: return fluid_.muGas(p, cells); @@ -573,7 +575,8 @@ namespace { return fluid_.bWat(p, cells); case Oil: { V dummy_rs = V::Zero(p.size(), 1) * p; - return fluid_.bOil(p, dummy_rs, cells); + bool dummy_isSat[p.size()]; + return fluid_.bOil(p, dummy_rs, dummy_isSat,cells); } case Gas: return fluid_.bGas(p, cells); @@ -593,7 +596,8 @@ namespace { return fluid_.bWat(p, cells); case Oil: { ADB dummy_rs = V::Zero(p.size(), 1) * p; - return fluid_.bOil(p, dummy_rs, cells); + bool dummy_isSat[p.size()]; + return fluid_.bOil(p, dummy_rs, dummy_isSat,cells); } case Gas: return fluid_.bGas(p, cells); From cb483e92cc0521b92d1346cf718c30f1f9a4276a Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Tue, 3 Dec 2013 18:12:54 +0100 Subject: [PATCH 2/2] Switch condition interface to phase presence facility Commit 4aa0eaf introduced density and viscosity evaluators into the BlackoilPropsAdInterface that accepted an externally assignable condition to distinguish saturated from unsaturated cases. As a result of a few low-level technical problems with that approach, this commit changes those affected interfaces to use the black-oil specific 'PhasePresence' facility of opm-core's commit a033329. Update callers accordingly. --- opm/autodiff/BlackoilPropsAd.cpp | 18 ++--- opm/autodiff/BlackoilPropsAd.hpp | 16 ++-- opm/autodiff/BlackoilPropsAdFromDeck.cpp | 27 +++---- opm/autodiff/BlackoilPropsAdFromDeck.hpp | 16 ++-- opm/autodiff/BlackoilPropsAdInterface.hpp | 16 ++-- opm/autodiff/FullyImplicitBlackoilSolver.cpp | 79 +++++++++++++------- opm/autodiff/FullyImplicitBlackoilSolver.hpp | 9 ++- opm/autodiff/ImpesTPFAAD.cpp | 20 +++-- 8 files changed, 118 insertions(+), 83 deletions(-) diff --git a/opm/autodiff/BlackoilPropsAd.cpp b/opm/autodiff/BlackoilPropsAd.cpp index 9c4827e25..0234f54b0 100644 --- a/opm/autodiff/BlackoilPropsAd.cpp +++ b/opm/autodiff/BlackoilPropsAd.cpp @@ -121,12 +121,12 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. V BlackoilPropsAd::muOil(const V& po, const V& rs, - const bool* /*isSat*/, + const std::vector& /*cond*/, const Cells& cells) const { if (!pu_.phase_used[Oil]) { @@ -199,16 +199,16 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. ADB BlackoilPropsAd::muOil(const ADB& po, const ADB& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const { #if 1 - return ADB::constant(muOil(po.value(), rs.value(), isSat,cells), po.blockPattern()); + return ADB::constant(muOil(po.value(), rs.value(), cond, cells), po.blockPattern()); #else if (!pu_.phase_used[Oil]) { OPM_THROW(std::runtime_error, "Cannot call muOil(): oil phase not present."); @@ -310,12 +310,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. V BlackoilPropsAd::bOil(const V& po, const V& rs, - const bool* /*isSat*/, + const std::vector& /*cond*/, const Cells& cells) const { if (!pu_.phase_used[Oil]) { @@ -388,12 +388,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. ADB BlackoilPropsAd::bOil(const ADB& po, const ADB& rs, - const bool* /*isSat*/, + const std::vector& /*cond*/, const Cells& cells) const { if (!pu_.phase_used[Oil]) { diff --git a/opm/autodiff/BlackoilPropsAd.hpp b/opm/autodiff/BlackoilPropsAd.hpp index 07e4559b7..2ed939509 100644 --- a/opm/autodiff/BlackoilPropsAd.hpp +++ b/opm/autodiff/BlackoilPropsAd.hpp @@ -109,12 +109,12 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. V muOil(const V& po, const V& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const; /// Gas viscosity. @@ -134,12 +134,12 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. ADB muOil(const ADB& po, const ADB& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const; /// Gas viscosity. @@ -162,12 +162,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. V bOil(const V& po, const V& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const; /// Gas formation volume factor. @@ -187,12 +187,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. ADB bOil(const ADB& po, const ADB& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const; /// Gas formation volume factor. diff --git a/opm/autodiff/BlackoilPropsAdFromDeck.cpp b/opm/autodiff/BlackoilPropsAdFromDeck.cpp index f86d958eb..9b4a57ad4 100644 --- a/opm/autodiff/BlackoilPropsAdFromDeck.cpp +++ b/opm/autodiff/BlackoilPropsAdFromDeck.cpp @@ -212,12 +212,12 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. V BlackoilPropsAdFromDeck::muOil(const V& po, const V& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const { if (!phase_usage_.phase_used[Oil]) { @@ -229,7 +229,7 @@ namespace Opm V dmudp(n); V dmudr(n); - props_[phase_usage_.phase_pos[Oil]]->mu(n, po.data(), rs.data(),isSat, + props_[phase_usage_.phase_pos[Oil]]->mu(n, po.data(), rs.data(), &cond[0], mu.data(), dmudp.data(), dmudr.data()); return mu; } @@ -287,12 +287,12 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. ADB BlackoilPropsAdFromDeck::muOil(const ADB& po, const ADB& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const { if (!phase_usage_.phase_used[Oil]) { @@ -304,8 +304,8 @@ namespace Opm V dmudp(n); V dmudr(n); - props_[phase_usage_.phase_pos[Oil]]->mu(n, po.value().data(), rs.value().data(), isSat, - mu.data(), dmudp.data(), dmudr.data()); + props_[phase_usage_.phase_pos[Oil]]->mu(n, po.value().data(), rs.value().data(), + &cond[0], mu.data(), dmudp.data(), dmudr.data()); ADB::M dmudp_diag = spdiag(dmudp); ADB::M dmudr_diag = spdiag(dmudr); @@ -391,11 +391,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. V BlackoilPropsAdFromDeck::bOil(const V& po, const V& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const { if (!phase_usage_.phase_used[Oil]) { @@ -408,7 +409,7 @@ namespace Opm V dbdp(n); V dbdr(n); - props_[phase_usage_.phase_pos[Oil]]->b(n, po.data(), rs.data(),isSat, + props_[phase_usage_.phase_pos[Oil]]->b(n, po.data(), rs.data(), &cond[0], b.data(), dbdp.data(), dbdr.data()); return b; @@ -471,12 +472,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. ADB BlackoilPropsAdFromDeck::bOil(const ADB& po, const ADB& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const { if (!phase_usage_.phase_used[Oil]) { @@ -489,8 +490,8 @@ namespace Opm V dbdp(n); V dbdr(n); - props_[phase_usage_.phase_pos[Oil]]->b(n, po.value().data(), rs.value().data(),isSat, - b.data(), dbdp.data(), dbdr.data()); + props_[phase_usage_.phase_pos[Oil]]->b(n, po.value().data(), rs.value().data(), + &cond[0], b.data(), dbdp.data(), dbdr.data()); ADB::M dbdp_diag = spdiag(dbdp); ADB::M dbdr_diag = spdiag(dbdr); diff --git a/opm/autodiff/BlackoilPropsAdFromDeck.hpp b/opm/autodiff/BlackoilPropsAdFromDeck.hpp index 3f4a8f4e5..9331a8e11 100644 --- a/opm/autodiff/BlackoilPropsAdFromDeck.hpp +++ b/opm/autodiff/BlackoilPropsAdFromDeck.hpp @@ -110,12 +110,12 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. V muOil(const V& po, const V& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const; /// Gas viscosity. @@ -135,12 +135,12 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. ADB muOil(const ADB& po, const ADB& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const; /// Gas viscosity. @@ -163,12 +163,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. V bOil(const V& po, const V& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const; /// Gas formation volume factor. @@ -188,12 +188,12 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. ADB bOil(const ADB& po, const ADB& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const; /// Gas formation volume factor. diff --git a/opm/autodiff/BlackoilPropsAdInterface.hpp b/opm/autodiff/BlackoilPropsAdInterface.hpp index 6b4928c48..abfe39da8 100644 --- a/opm/autodiff/BlackoilPropsAdInterface.hpp +++ b/opm/autodiff/BlackoilPropsAdInterface.hpp @@ -100,13 +100,13 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. virtual V muOil(const V& po, const V& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const = 0; /// Gas viscosity. @@ -128,13 +128,13 @@ namespace Opm /// Oil viscosity. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n viscosity values. virtual ADB muOil(const ADB& po, const ADB& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const = 0; /// Gas viscosity. @@ -159,13 +159,13 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. virtual V bOil(const V& po, const V& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const = 0; /// Gas formation volume factor. @@ -187,13 +187,13 @@ namespace Opm /// Oil formation volume factor. /// \param[in] po Array of n oil pressure values. /// \param[in] rs Array of n gas solution factor values. - /// \param[in] isSat Array of n booleans telling whether the fluid is saturated or not. + /// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \return Array of n formation volume factor values. virtual ADB bOil(const ADB& po, const ADB& rs, - const bool* isSat, + const std::vector& cond, const Cells& cells) const = 0; /// Gas formation volume factor. diff --git a/opm/autodiff/FullyImplicitBlackoilSolver.cpp b/opm/autodiff/FullyImplicitBlackoilSolver.cpp index 2db07478c..260100bac 100644 --- a/opm/autodiff/FullyImplicitBlackoilSolver.cpp +++ b/opm/autodiff/FullyImplicitBlackoilSolver.cpp @@ -528,8 +528,8 @@ namespace { const std::vector& sat = state.saturation; const ADB& rs = state.rs; - bool isSat[rs.size()]; - getSaturatedCells(state,&isSat[0]); + std::vector cond; + classifyCondition(state, cond); const ADB pv_mult = poroMult(press); @@ -537,7 +537,7 @@ namespace { for (int phase = 0; phase < maxnp; ++phase) { if (active_[ phase ]) { const int pos = pu.phase_pos[ phase ]; - rq_[pos].b = fluidReciprocFVF(phase, press, rs, &isSat[0], cells_); + rq_[pos].b = fluidReciprocFVF(phase, press, rs, cond, cells_); rq_[pos].accum[aix] = pv_mult * rq_[pos].b * sat[pos]; // DUMP(rq_[pos].b); // DUMP(rq_[pos].accum[aix]); @@ -595,8 +595,6 @@ namespace { // DUMP(residual_.mass_balance[phase]); } - bool isSat[grid_.number_of_cells]; - getSaturatedCells(state,isSat); // -------- Extra (optional) sg or rs equation, and rs contributions to the mass balance equations -------- // Add the extra (flux) terms to the gas mass balance equations @@ -654,11 +652,15 @@ namespace { assert(g[dd] == 0.0); } } + + std::vector cond; + classifyCondition(state, cond); + ADB cell_rho_total = ADB::constant(V::Zero(nc), state.pressure.blockPattern()); for (int phase = 0; phase < 3; ++phase) { if (active_[phase]) { const int pos = pu.phase_pos[phase]; - const ADB cell_rho = fluidDensity(phase, state.pressure, state.rs, &isSat[0],cells_); + const ADB cell_rho = fluidDensity(phase, state.pressure, state.rs, cond, cells_); cell_rho_total += state.saturation[pos] * cell_rho; } } @@ -668,7 +670,7 @@ namespace { for (int phase = 0; phase < 3; ++phase) { if (active_[phase]) { const int pos = pu.phase_pos[phase]; - const ADB cell_rho = fluidDensity(phase, state.pressure, state.rs, &isSat[0],cells_); + const ADB cell_rho = fluidDensity(phase, state.pressure, state.rs, cond, cells_); const V fraction = compi.col(pos); inj_rho_total += (wops_.w2p * fraction.matrix()).array() * subset(cell_rho, well_cells); } @@ -1007,14 +1009,16 @@ namespace { const SolutionState& state ) { const int phase = canph_[ actph ]; - bool isSat[grid_.number_of_cells]; - getSaturatedCells(state,&isSat[0]); - const ADB mu = fluidViscosity(phase, state.pressure, state.rs, &isSat[0],cells_); + + std::vector cond; + classifyCondition(state, cond); + + const ADB mu = fluidViscosity(phase, state.pressure, state.rs, cond, cells_); const ADB tr_mult = transMult(state.pressure); rq_[ actph ].mob = tr_mult * kr[ phase ] / mu; - const ADB rho = fluidDensity(phase, state.pressure, state.rs, &isSat[0],cells_); + const ADB rho = fluidDensity(phase, state.pressure, state.rs, cond, cells_); ADB& head = rq_[ actph ].head; @@ -1067,14 +1071,14 @@ namespace { FullyImplicitBlackoilSolver::fluidViscosity(const int phase, const ADB& p , const ADB& rs , - const bool* isSat, + const std::vector& cond, const std::vector& cells) const { switch (phase) { case Water: return fluid_.muWat(p, cells); case Oil: { - return fluid_.muOil(p, rs, isSat,cells); + return fluid_.muOil(p, rs, cond, cells); } case Gas: return fluid_.muGas(p, cells); @@ -1091,14 +1095,14 @@ namespace { FullyImplicitBlackoilSolver::fluidReciprocFVF(const int phase, const ADB& p , const ADB& rs , - const bool* isSat, + const std::vector& cond, const std::vector& cells) const { switch (phase) { case Water: return fluid_.bWat(p, cells); case Oil: { - return fluid_.bOil(p, rs, isSat, cells); + return fluid_.bOil(p, rs, cond, cells); } case Gas: return fluid_.bGas(p, cells); @@ -1115,11 +1119,11 @@ namespace { FullyImplicitBlackoilSolver::fluidDensity(const int phase, const ADB& p , const ADB& rs , - const bool* isSat, + const std::vector& cond, const std::vector& cells) const { const double* rhos = fluid_.surfaceDensity(); - ADB b = fluidReciprocFVF(phase, p, rs, isSat, cells); + ADB b = fluidReciprocFVF(phase, p, rs, cond, cells); ADB rho = V::Constant(p.size(), 1, rhos[phase]) * b; if (phase == Oil && active_[Gas]) { // It is correct to index into rhos with canonical phase indices. @@ -1206,16 +1210,41 @@ namespace { void - FullyImplicitBlackoilSolver::getSaturatedCells(const SolutionState& state, bool* isSat) const + FullyImplicitBlackoilSolver:: + classifyCondition(const SolutionState& state, + std::vector& cond ) const { - const int pg = fluid_.phaseUsage().phase_pos[ Gas ]; - const V sg = state.saturation[pg].value(); - for (int c=0; c < sg.size(); ++ c) { - if (sg[c]>0){ - isSat[c] = true; + const PhaseUsage& pu = fluid_.phaseUsage(); + + if (active_[ Gas ]) { + // Oil/Gas or Water/Oil/Gas system + const int po = pu.phase_pos[ Oil ]; + const int pg = pu.phase_pos[ Gas ]; + + const V& so = state.saturation[ po ].value(); + const V& sg = state.saturation[ pg ].value(); + + cond.resize(sg.size()); + + for (V::Index c = 0, e = sg.size(); c != e; ++c) { + if (so[c] > 0) { cond[c].setFreeOil (); } + if (sg[c] > 0) { cond[c].setFreeGas (); } + if (active_[ Water ]) { cond[c].setFreeWater(); } } - else{ - isSat[c] = false; + } + else { + // Water/Oil system + assert (active_[ Water ]); + + const int po = pu.phase_pos[ Oil ]; + const V& so = state.saturation[ po ].value(); + + cond.resize(so.size()); + + for (V::Index c = 0, e = so.size(); c != e; ++c) { + cond[c].setFreeWater(); + + if (so[c] > 0) { cond[c].setFreeOil(); } } } } diff --git a/opm/autodiff/FullyImplicitBlackoilSolver.hpp b/opm/autodiff/FullyImplicitBlackoilSolver.hpp index f440c44b0..815c8a4fd 100644 --- a/opm/autodiff/FullyImplicitBlackoilSolver.hpp +++ b/opm/autodiff/FullyImplicitBlackoilSolver.hpp @@ -189,21 +189,21 @@ namespace Opm { fluidViscosity(const int phase, const ADB& p , const ADB& rs , - const bool* isSat, + const std::vector& cond, const std::vector& cells) const; ADB fluidReciprocFVF(const int phase, const ADB& p , const ADB& rs , - const bool* isSat, + const std::vector& cond, const std::vector& cells) const; ADB fluidDensity(const int phase, const ADB& p , const ADB& rs , - const bool* isSat, + const std::vector& cond, const std::vector& cells) const; V @@ -221,7 +221,8 @@ namespace Opm { transMult(const ADB& p) const; void - getSaturatedCells(const SolutionState& state, bool* isSat) const; + classifyCondition(const SolutionState& state, + std::vector& cond ) const; }; } // namespace Opm diff --git a/opm/autodiff/ImpesTPFAAD.cpp b/opm/autodiff/ImpesTPFAAD.cpp index abe6367f8..292031345 100644 --- a/opm/autodiff/ImpesTPFAAD.cpp +++ b/opm/autodiff/ImpesTPFAAD.cpp @@ -533,8 +533,9 @@ namespace { return fluid_.muWat(p, cells); case Oil: { V dummy_rs = V::Zero(p.size(), 1) * p; - bool dummy_isSat[p.size()]; - return fluid_.muOil(p, dummy_rs, dummy_isSat, cells); + std::vector cond(dummy_rs.size()); + + return fluid_.muOil(p, dummy_rs, cond, cells); } case Gas: return fluid_.muGas(p, cells); @@ -554,8 +555,9 @@ namespace { return fluid_.muWat(p, cells); case Oil: { ADB dummy_rs = V::Zero(p.size(), 1) * p; - bool dummy_isSat[p.size()]; - return fluid_.muOil(p, dummy_rs, dummy_isSat, cells); + std::vector cond(dummy_rs.size()); + + return fluid_.muOil(p, dummy_rs, cond, cells); } case Gas: return fluid_.muGas(p, cells); @@ -575,8 +577,9 @@ namespace { return fluid_.bWat(p, cells); case Oil: { V dummy_rs = V::Zero(p.size(), 1) * p; - bool dummy_isSat[p.size()]; - return fluid_.bOil(p, dummy_rs, dummy_isSat,cells); + std::vector cond(dummy_rs.size()); + + return fluid_.bOil(p, dummy_rs, cond, cells); } case Gas: return fluid_.bGas(p, cells); @@ -596,8 +599,9 @@ namespace { return fluid_.bWat(p, cells); case Oil: { ADB dummy_rs = V::Zero(p.size(), 1) * p; - bool dummy_isSat[p.size()]; - return fluid_.bOil(p, dummy_rs, dummy_isSat,cells); + std::vector cond(dummy_rs.size()); + + return fluid_.bOil(p, dummy_rs, cond, cells); } case Gas: return fluid_.bGas(p, cells);