diff --git a/ApplicationCode/ReservoirDataModel/RigFlowDiagSolverInterface.cpp b/ApplicationCode/ReservoirDataModel/RigFlowDiagSolverInterface.cpp index c943688b2e..50d3d9bad8 100644 --- a/ApplicationCode/ReservoirDataModel/RigFlowDiagSolverInterface.cpp +++ b/ApplicationCode/ReservoirDataModel/RigFlowDiagSolverInterface.cpp @@ -682,7 +682,6 @@ std::vector RigFlowDiagSolverInterface if (!useEps) { scaling.enable = static_cast(0); } - scaling.invalid = Opm::SatFunc::EPSEvalInterface::InvalidEndpointBehaviour::IgnorePoint; std::vector graphArr = m_opmFlowDiagStaticData->m_eclSaturationFunc->getSatFuncCurve(satFuncRequests, static_cast(activeCellIndex), scaling); for (size_t i = 0; i < graphArr.size(); i++) { diff --git a/ThirdParty/custom-opm-flowdiag-app/CMakeLists.txt b/ThirdParty/custom-opm-flowdiag-app/CMakeLists.txt index c37bd33e02..d54ae136a9 100644 --- a/ThirdParty/custom-opm-flowdiag-app/CMakeLists.txt +++ b/ThirdParty/custom-opm-flowdiag-app/CMakeLists.txt @@ -5,7 +5,6 @@ project (custom-opm-flowdiag-app) include_directories( ../custom-opm-flowdiagnostics/opm-flowdiagnostics opm-flowdiagnostics-applications - opmCore ) include (opm-flowdiagnostics-applications/CMakeLists_files.cmake) diff --git a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/examples/extractPropCurves.cpp b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/examples/extractPropCurves.cpp index cbfbfed881..acf60e4dd4 100644 --- a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/examples/extractPropCurves.cpp +++ b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/examples/extractPropCurves.cpp @@ -54,13 +54,13 @@ namespace { const auto& x = graph.first; const auto& y = graph.second; - os << name << '{' << k << "} = [\n"; + os << name << '{' << k << "} = extendTab([\n"; for (auto n = x.size(), i = 0*n; i < n; ++i) { os << x[i] << ' ' << y[i] << '\n'; } - os << "];\n\n"; + os << "]);\n\n"; k += 1; } @@ -118,7 +118,7 @@ namespace { const auto graph = sfunc.getSatFuncCurve(func, activeCell, scaling); - printGraph(std::cout, "krg", graph); + printGraph(std::cout, "crv.krg", graph); } void krog(const Opm::ECLSaturationFunc& sfunc, @@ -140,7 +140,7 @@ namespace { const auto graph = sfunc.getSatFuncCurve(func, activeCell, scaling); - printGraph(std::cout, "krog", graph); + printGraph(std::cout, "crv.krog", graph); } void krow(const Opm::ECLSaturationFunc& sfunc, @@ -162,7 +162,7 @@ namespace { const auto graph = sfunc.getSatFuncCurve(func, activeCell, scaling); - printGraph(std::cout, "krow", graph); + printGraph(std::cout, "crv.krow", graph); } void krw(const Opm::ECLSaturationFunc& sfunc, @@ -184,7 +184,7 @@ namespace { const auto graph = sfunc.getSatFuncCurve(func, activeCell, scaling); - printGraph(std::cout, "krw", graph); + printGraph(std::cout, "crv.krw", graph); } // ----------------------------------------------------------------- @@ -209,7 +209,7 @@ namespace { const auto graph = sfunc.getSatFuncCurve(func, activeCell, scaling); - printGraph(std::cout, "pcgo", graph); + printGraph(std::cout, "crv.pcgo", graph); } void pcow(const Opm::ECLSaturationFunc& sfunc, @@ -231,7 +231,7 @@ namespace { const auto graph = sfunc.getSatFuncCurve(func, activeCell, scaling); - printGraph(std::cout, "pcow", graph); + printGraph(std::cout, "crv.pcow", graph); } // ----------------------------------------------------------------- @@ -245,7 +245,7 @@ namespace { const auto graph = pvtCurves .getPvtCurve(RC::FVF, Opm::ECLPhaseIndex::Vapour, activeCell); - printGraph(std::cout, "Bg", graph); + printGraph(std::cout, "crv.Bg", graph); } void mu_g(const Opm::ECLPVT::ECLPvtCurveCollection& pvtCurves, @@ -256,7 +256,7 @@ namespace { const auto graph = pvtCurves .getPvtCurve(RC::Viscosity, Opm::ECLPhaseIndex::Vapour, activeCell); - printGraph(std::cout, "mu_g", graph); + printGraph(std::cout, "crv.mu_g", graph); } void Bo(const Opm::ECLPVT::ECLPvtCurveCollection& pvtCurves, @@ -267,7 +267,7 @@ namespace { const auto graph = pvtCurves .getPvtCurve(RC::FVF, Opm::ECLPhaseIndex::Liquid, activeCell); - printGraph(std::cout, "Bo", graph); + printGraph(std::cout, "crv.Bo", graph); } void mu_o(const Opm::ECLPVT::ECLPvtCurveCollection& pvtCurves, @@ -278,7 +278,7 @@ namespace { const auto graph = pvtCurves .getPvtCurve(RC::Viscosity, Opm::ECLPhaseIndex::Liquid, activeCell); - printGraph(std::cout, "mu_o", graph); + printGraph(std::cout, "crv.mu_o", graph); } // ----------------------------------------------------------------- @@ -293,7 +293,7 @@ namespace { const auto graph = pvtCurves .getPvtCurve(RC::SaturatedState, PI::Vapour, activeCell); - printGraph(std::cout, "rvSat", graph); + printGraph(std::cout, "crv.rvSat", graph); } void rsSat(const Opm::ECLPVT::ECLPvtCurveCollection& pvtCurves, @@ -305,7 +305,7 @@ namespace { const auto graph = pvtCurves .getPvtCurve(RC::SaturatedState, PI::Liquid, activeCell); - printGraph(std::cout, "rsSat", graph); + printGraph(std::cout, "crv.rsSat", graph); } // ----------------------------------------------------------------- @@ -340,37 +340,6 @@ namespace { return Opm::ECLUnits::serialisedUnitConventions(init); } - Opm::SatFunc::EPSEvalInterface::InvalidEndpointBehaviour - handleInvalid(const std::string& behaviour) - { - using IEB = Opm::SatFunc:: - EPSEvalInterface::InvalidEndpointBehaviour; - - if ((behaviour == "ignore") || - (behaviour == "ignore-point") || - (behaviour == "ignore_point") || - (behaviour == "ignorepoint")) - { - return IEB::IgnorePoint; - } - - return IEB::UseUnscaled; - } - - auto handleInvalid(const Opm::ParameterGroup& prm) - -> decltype(handleInvalid("ignore")) - { - for (const auto* param : { "handle_invalid" , - "hInv", "handleInv" }) - { - if (prm.has(param)) { - return handleInvalid(prm.get(param)); - } - } - - return handleInvalid("ignore"); - } - int getActiveCell(const Opm::ECLGraph& G, const Opm::ParameterGroup& prm) { @@ -424,8 +393,7 @@ namespace { prm.getDefault("useEPS", std::string{"off"}); auto scaling = Opm::ECLSaturationFunc::SatFuncScaling{}; - scaling.enable = static_cast(0); - scaling.invalid = handleInvalid(prm); + scaling.enable = static_cast(0); if (std::regex_search(useEPS, horiz)) { scaling.enable |= T::Horizontal; @@ -460,6 +428,8 @@ try { pvtCC.setOutputUnits(std::move(units)); } + std::cout << "function crv = pcurves()\n"; + // ----------------------------------------------------------------- // Relative permeability @@ -470,7 +440,8 @@ try { // ----------------------------------------------------------------- // Capillary pressure - if (prm.getDefault("pcgo", false)) { pcgo(sfunc, cellID, scaling); } + if (prm.getDefault("pcog", false) || // Alias pcog -> pcgo + prm.getDefault("pcgo", false)) { pcgo(sfunc, cellID, scaling); } if (prm.getDefault("pcow", false)) { pcow(sfunc, cellID, scaling); } // ----------------------------------------------------------------- diff --git a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLEndPointScaling.cpp b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLEndPointScaling.cpp index 47b5024db6..6b16d2c88b 100644 --- a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLEndPointScaling.cpp +++ b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLEndPointScaling.cpp @@ -1,4 +1,5 @@ /* + Copyright 2018 Equinor ASA. Copyright 2017 Statoil ASA. This file is part of the Open Porous Media Project (OPM). @@ -34,8 +35,8 @@ #include #include #include -#include #include +#include #include #include #include @@ -217,21 +218,6 @@ namespace { return (std::abs(s) < 1.0e+20) ? s : dflt; } - bool validSaturation(const double s) - { - return (! (s < 0.0)) && (! (s > 1.0)); - } - - bool validSaturations(std::initializer_list sats) - { - return std::accumulate(std::begin(sats), - std::end (sats), true, - [](const bool result, const double s) -> bool - { - return result && validSaturation(s); - }); - } - bool haveScaledRelPermAtCritSat(const ::Opm::ECLGraph& G, const ::Opm::ECLInitFileData& init, @@ -263,9 +249,7 @@ class Opm::SatFunc::TwoPointScaling::Impl public: Impl(std::vector smin, std::vector smax) - : smin_ (std::move(smin)) - , smax_ (std::move(smax)) - , handle_invalid_(InvalidEndpointBehaviour::UseUnscaled) + : smin_(std::move(smin)), smax_(std::move(smax)) { if (this->smin_.size() != this->smax_.size()) { throw std::invalid_argument { @@ -287,11 +271,6 @@ private: std::vector smin_; std::vector smax_; - InvalidEndpointBehaviour handle_invalid_; - - void handleInvalidEndpoint(const SaturationAssoc& sp, - std::vector& effsat) const; - double sMin(const std::vector::size_type cell, const TableEndPoints& tep) const { @@ -327,12 +306,6 @@ Impl::eval(const TableEndPoints& tep, const auto sLO = this->sMin(cell, tep); const auto sHI = this->sMax(cell, tep); - if (! validSaturations({ sLO, sHI })) { - this->handleInvalidEndpoint(eval_pt, effsat); - - continue; - } - effsat.push_back(0.0); auto& s_eff = effsat.back(); @@ -372,12 +345,6 @@ Impl::reverse(const TableEndPoints& tep, const auto sLO = this->sMin(cell, tep); const auto sHI = this->sMax(cell, tep); - if (! validSaturations({ sLO, sHI })) { - this->handleInvalidEndpoint(eval_pt, unscaledsat); - - continue; - } - unscaledsat.push_back(0.0); auto& s_unsc = unscaledsat.back(); @@ -404,26 +371,6 @@ Impl::reverse(const TableEndPoints& tep, return unscaledsat; } -void -Opm::SatFunc::TwoPointScaling::Impl:: -handleInvalidEndpoint(const SaturationAssoc& sp, - std::vector& effsat) const -{ - if (this->handle_invalid_ == InvalidEndpointBehaviour::UseUnscaled) { - // User requests that invalid scaling be treated as unscaled - // saturations. Pick that. - effsat.push_back(sp.sat); - return; - } - - if (this->handle_invalid_ == InvalidEndpointBehaviour::IgnorePoint) { - // User requests that invalid scaling be ignored. Signal invalid - // scaled saturation to caller as NaN. - effsat.push_back(std::nan("")); - return; - } -} - // --------------------------------------------------------------------- // Class Opm::SatFunc::PureVerticalScaling::Impl // --------------------------------------------------------------------- @@ -433,6 +380,9 @@ class Opm::SatFunc::PureVerticalScaling::Impl public: explicit Impl(std::vector fmax) : fmax_(std::move(fmax)) + , inf_ (std::numeric_limits::has_infinity + ? std::numeric_limits::infinity() + : std::numeric_limits::max()) {} std::vector @@ -442,6 +392,16 @@ public: private: std::vector fmax_; + + double inf_; + + std::vector + zeroMaxVal(const SaturationPoints& sp) const; + + std::vector + nonZeroMaxVal(const SaturationPoints& sp, + const double maxVal, + std::vector val) const; }; std::vector @@ -452,10 +412,40 @@ vertScale(const FunctionValues& f, { assert (sp.size() == val.size() && "Internal Error in Vertical Scaling"); - auto ret = std::move(val); - const auto maxVal = f.max.val; + if (! (std::abs(maxVal) > 0.0)) { + return this->zeroMaxVal(sp); + } + + return this->nonZeroMaxVal(sp, maxVal, std::move(val)); +} + +std::vector +Opm::SatFunc::PureVerticalScaling::Impl:: +zeroMaxVal(const SaturationPoints& sp) const +{ + auto ret = std::vector(sp.size(), 0.0); + + for (auto n = sp.size(), i = 0*n; i < n; ++i) { + const auto fmax = this->fmax_[ sp[i].cell ]; + + ret[i] = (std::abs(fmax) > 0.0) + ? (std::signbit(fmax) ? -this->inf_ : this->inf_) + : 0.0; + } + + return ret; +} + +std::vector +Opm::SatFunc::PureVerticalScaling::Impl:: +nonZeroMaxVal(const SaturationPoints& sp, + const double maxVal, + std::vector val) const +{ + auto ret = std::move(val); + for (auto n = sp.size(), i = 0*n; i < n; ++i) { ret[i] *= this->fmax_[ sp[i].cell ] / maxVal; } @@ -513,24 +503,29 @@ vertScale(const FunctionValues& f, const auto c = sp[i].cell; const auto s = sp[i].sat; - const auto sr = this->sdisp_[c]; + const auto sr = std::min(this->sdisp_[c], this->smax_[c]); const auto fr = this->fdisp_[c]; const auto sm = this->smax_ [c]; const auto fm = this->fmax_ [c]; - if (! (s > sr)) { - // s <= sr: Pure vertical scaling in left interval. + if (s < sr) { + // s < sr: Pure vertical scaling in left interval. y *= fr / fdisp; } else if (sepfv) { - // s > sr; normal case: Kr(Smax) > Kr(Sr) + // s \in [sr, sm), sm > sr; normal case: Kr(Smax) > Kr(Sr). + // + // Linear function between (sr,fr) and (sm,fm) in terms of + // function value 'y'. const auto t = (y - fdisp) / (fmax - fdisp); y = fr + t*(fm - fr); } else if (s < sm) { - // s > sr; special case: Kr(Smax) == Kr(Sr) in table. Use - // linear function between (sr,fr) and (sm,fm). + // s \in [sr, sm), sm > sr; special case: Kr(Smax) == Kr(Sr). + // + // Use linear function between (sr,fr) and (sm,fm) in terms of + // saturation value 's'. const auto t = (s - sr) / (sm - sr); y = fr + t*(fm - fr); @@ -554,10 +549,7 @@ public: Impl(std::vector smin , std::vector sdisp, std::vector smax ) - : smin_ (std::move(smin )) - , sdisp_ (std::move(sdisp)) - , smax_ (std::move(smax )) - , handle_invalid_(InvalidEndpointBehaviour::UseUnscaled) + : smin_(std::move(smin)), sdisp_(std::move(sdisp)), smax_(std::move(smax)) { if ((this->sdisp_.size() != this->smin_.size()) || (this->sdisp_.size() != this->smax_.size())) @@ -582,11 +574,6 @@ private: std::vector sdisp_; std::vector smax_; - InvalidEndpointBehaviour handle_invalid_; - - void handleInvalidEndpoint(const SaturationAssoc& sp, - std::vector& effsat) const; - double sMin(const std::vector::size_type cell, const TableEndPoints& tep) const { @@ -630,12 +617,6 @@ Impl::eval(const TableEndPoints& tep, const auto sR = this->sDisp(cell, tep); const auto sHI = this->sMax (cell, tep); - if (! validSaturations({ sLO, sR, sHI })) { - this->handleInvalidEndpoint(eval_pt, effsat); - - continue; - } - effsat.push_back(0.0); auto& s_eff = effsat.back(); @@ -643,22 +624,26 @@ Impl::eval(const TableEndPoints& tep, // s <= sLO s_eff = tep.low; } - else if (! (eval_pt.sat < sHI)) { - // s >= sHI - s_eff = tep.high; - } - else if (eval_pt.sat < sR) { - // s \in (sLO, sR) + else if (eval_pt.sat < std::min(sR, sHI)) { + // s in scaled interval [sLO, sR) + // Map to tabulated saturation in [tep.low, tep.disp) const auto t = (eval_pt.sat - sLO) / (sR - sLO); s_eff = tep.low + t*(tep.disp - tep.low); } - else { - // s \in (sR, sHI) + else if (eval_pt.sat < sHI) { + // s in scaled interval [sR, sHI) + // Map to tabulated saturation in [tep.disp, tep.high) + assert (sHI > sR); + const auto t = (eval_pt.sat - sR) / (sHI - sR); s_eff = tep.disp + t*(tep.high - tep.disp); } + else { + // s >= sHI + s_eff = tep.high; + } } return effsat; @@ -679,12 +664,6 @@ Impl::reverse(const TableEndPoints& tep, const auto sR = this->sDisp(cell, tep); const auto sHI = this->sMax (cell, tep); - if (! validSaturations({ sLO, sR, sHI })) { - this->handleInvalidEndpoint(eval_pt, unscaledsat); - - continue; - } - unscaledsat.push_back(0.0); auto& s_unsc = unscaledsat.back(); @@ -693,54 +672,45 @@ Impl::reverse(const TableEndPoints& tep, // Map to Minimum Input Saturation in cell (sLO). s_unsc = sLO; } - else if (! (eval_pt.sat < tep.high)) { - // s >= maximum tabulated saturation. - // Map to Maximum Input Saturation in cell (sHI). - s_unsc = sHI; - } else if (eval_pt.sat < tep.disp) { - // s in tabulated interval (tep.low, tep.disp) + // s in tabulated interval [tep.low, tep.disp) // Map to Input Saturation in (sLO, sR) const auto t = (eval_pt.sat - tep.low) / (tep.disp - tep.low); - s_unsc = sLO + t*(sR - sLO); + s_unsc = std::min(sLO + t*(sR - sLO), sHI); } - else { - // s in tabulated interval (tep.disp, tep.high) - // Map to Input Saturation in (sR, sHI) + else if (eval_pt.sat < tep.high) { + // s in tabulated interval [tep.disp, tep.high) + // Map to Input Saturation in [sR, sHI) + assert (tep.high > tep.disp); + const auto t = (eval_pt.sat - tep.disp) / (tep.high - tep.disp); - s_unsc = sR + t*(sHI - sR); + s_unsc = std::min(sR + t*std::max(sHI - sR, 0.0), sHI); + } + else { + // s >= maximum tabulated saturation. + // + // Map to Maximum Input Saturation in cell (sHI) if maximum + // tabulated saturation is strictly greater than critical + // displacing saturation--otherwise map to critical displacing + // saturation. + // + // Needed to handle cases in which \code tep.disp==tep.high + // \endcode but scaled versions of these might differ, i.e. when + // sR < sHI, but the corresponding saturation points in the + // underlying input table coincide. + s_unsc = (tep.high > tep.disp) ? sHI : sR; } } return unscaledsat; } -void -Opm::SatFunc::ThreePointScaling::Impl:: -handleInvalidEndpoint(const SaturationAssoc& sp, - std::vector& effsat) const -{ - if (this->handle_invalid_ == InvalidEndpointBehaviour::UseUnscaled) { - // User requests that invalid scaling be treated as unscaled - // saturations. Pick that. - effsat.push_back(sp.sat); - return; - } - - if (this->handle_invalid_ == InvalidEndpointBehaviour::IgnorePoint) { - // User requests that invalid scaling be ignored. Signal invalid - // scaled saturation to caller as NaN. - effsat.push_back(std::nan("")); - return; - } -} - // --------------------------------------------------------------------- // EPS factory functions for two-point and three-point scaling options // --------------------------------------------------------------------- @@ -749,7 +719,6 @@ namespace Create { using EPSOpt = ::Opm::SatFunc::CreateEPS::EPSOptions; using RTEP = ::Opm::SatFunc::CreateEPS::RawTableEndPoints; using TEP = ::Opm::SatFunc::EPSEvalInterface::TableEndPoints; - using InvBeh = ::Opm::SatFunc::EPSEvalInterface::InvalidEndpointBehaviour; namespace TwoPoint { using EPS = ::Opm::SatFunc::TwoPointScaling; @@ -1742,7 +1711,7 @@ KrGO::twoPointMethod(const ::Opm::ECLGraph& G, auto t = std::vector(sogcr.size(), 0.0); { const auto& sgco = tep.conn.gas; - const auto& swco = tep.conn.gas; + const auto& swco = tep.conn.water; const auto& sgcr = tep.crit.gas; for (auto n = sgcr.size(), i = 0*n; i < n; ++i) { @@ -1790,7 +1759,7 @@ KrOW::twoPointMethod(const ::Opm::ECLGraph& G, auto t = std::vector(sowcr.size(), 0.0); { const auto& sgco = tep.conn.gas; - const auto& swco = tep.conn.gas; + const auto& swco = tep.conn.water; const auto& swcr = tep.crit.water; for (auto n = swcr.size(), i = 0*n; i < n; ++i) { diff --git a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLEndPointScaling.hpp b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLEndPointScaling.hpp index 384ce66c97..4fbc84aba2 100644 --- a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLEndPointScaling.hpp +++ b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLEndPointScaling.hpp @@ -67,19 +67,6 @@ namespace Opm { namespace SatFunc { double sat; }; - /// Policy for how to handle an invalid end-point scaling (e.g., if - /// lower and/or upper scaled saturations have nonsensical values - /// like -1.0E+20). - enum class InvalidEndpointBehaviour { - /// Use the unscaled value for this point. - UseUnscaled, - - /// Ignore this scaling request (e.g., produce no graphical - /// output for a scaled saturation function when the scaled - /// end-points are meaningless). - IgnorePoint, - }; - /// Convenience type alias. using SaturationPoints = std::vector; @@ -574,13 +561,6 @@ namespace Opm { namespace SatFunc { /// auto eps = CreateEPS::fromECLOutput(G, init, opt); /// \endcode ::Opm::ECLPhaseIndex thisPh; - - /// How to handle an invalid end-point scaling (e.g., if lower - /// and/or upper scaled saturations have nonsensical values like - /// -1.0E+20). - EPSEvalInterface::InvalidEndpointBehaviour handle_invalid { - EPSEvalInterface::InvalidEndpointBehaviour::UseUnscaled - }; }; /// Collection of raw saturation table end points. diff --git a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLPvtGas.cpp b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLPvtGas.cpp index fef7dd5140..c712476ee3 100644 --- a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLPvtGas.cpp +++ b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLPvtGas.cpp @@ -90,6 +90,8 @@ namespace { class PVxGBase { public: + virtual ~PVxGBase() {} + virtual std::vector formationVolumeFactor(const std::vector& rv, const std::vector& pg) const = 0; diff --git a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLPvtOil.cpp b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLPvtOil.cpp index 362ae167f8..7ace98f6c7 100644 --- a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLPvtOil.cpp +++ b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLPvtOil.cpp @@ -129,6 +129,8 @@ namespace { class PVxOBase { public: + virtual ~PVxOBase() {} + virtual std::vector formationVolumeFactor(const std::vector& rs, const std::vector& po) const = 0; diff --git a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLResultData.cpp b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLResultData.cpp index d2e977cd5f..d0effa9510 100644 --- a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLResultData.cpp +++ b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLResultData.cpp @@ -473,13 +473,13 @@ namespace { explicit InitFileSections(const ecl_file_type* init); struct Section { - Section(const ecl_file_view_type* blk) + Section(ecl_file_view_type* blk) : block (blk) , first_kw(Details::firstBlockKeyword(block)) {} - const ecl_file_view_type* block; - std::string first_kw; + ecl_file_view_type* block; + std::string first_kw; }; const std::vector
& sections() const @@ -502,7 +502,7 @@ namespace { } private: - const ecl_file_view_type* init_; + ecl_file_view_type* init_; std::vector
sect_; }; @@ -796,7 +796,7 @@ ECLImpl::InitFileSections::InitFileSections(const ecl_file_type* init) ? sectID - 1 : 0*sectID; // Main section 'sectID': [ start_kw, LGRSGONE ] - const auto* sect = + auto* sect = ecl_file_view_add_blockview2(this->init_, start_kw, end_kw, start_kw_occurrence); @@ -808,12 +808,12 @@ ECLImpl::InitFileSections::InitFileSections(const ecl_file_type* init) const auto occurrence = 0; // Main grid sub-section of 'sectID': [ start_kw, LGR ] - const auto* main_grid_sect = + auto* main_grid_sect = ecl_file_view_add_blockview2(sect, firstkw.c_str(), LGR_KW, occurrence); // LGR sub-section of 'sectID': [ LGR, LGRSGONE ] - const auto* lgr_sect = + auto* lgr_sect = ecl_file_view_add_blockview2(sect, LGR_KW, end_kw, occurrence); @@ -965,7 +965,7 @@ private: /// Saved original active view from host (prior to restricting view /// to single grid ID). - const ecl_file_view_type* save_; + ecl_file_view_type* save_; }; /// Casename prefix. Mostly to implement copy ctor. @@ -986,7 +986,7 @@ private: std::unique_ptr gridIDCache_; /// Current active result-set view. - mutable const ecl_file_view_type* activeBlock_{ nullptr }; + mutable ecl_file_view_type* activeBlock_{ nullptr }; /// Support for passing \code *this \endcode to ERT functions that /// require an \c ecl_file_type, particularly the function that selects @@ -1276,7 +1276,7 @@ private: /// Sections of the INIT result set. ECLImpl::InitFileSections sections_; - mutable const ecl_file_view_type* activeBlock_{ nullptr }; + mutable ecl_file_view_type* activeBlock_{ nullptr }; /// Negative look-up cache for haveKeywordData() queries. mutable MissingKW missing_kw_; @@ -1415,6 +1415,10 @@ lookup(const std::string& vector, const std::string& gridName) const } } + // Status of 'vector' unknown for 'gridName'. Actually look for the + // vector in gridName's sections (main grid if gridName.empty(), + // otherwise named local grid). + if (gridName.empty()) { return this->lookupMainGrid(key); } diff --git a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLSaturationFunc.cpp b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLSaturationFunc.cpp index d7c8c3a72b..0465752cc2 100644 --- a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLSaturationFunc.cpp +++ b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLSaturationFunc.cpp @@ -758,8 +758,7 @@ private: const ActPh& active) { auto opt = Create::EPSOptions{}; - opt.use3PtScaling = use3PtScaling; - opt.handle_invalid = SatFuncScaling::IEB::UseUnscaled; + opt.use3PtScaling = use3PtScaling; if (active.oil) { this->create_oil_eps(host, G, init, ep, active, opt); diff --git a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLSaturationFunc.hpp b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLSaturationFunc.hpp index ff9926d0aa..7b94f868c6 100644 --- a/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLSaturationFunc.hpp +++ b/ThirdParty/custom-opm-flowdiag-app/opm-flowdiagnostics-applications/opm/utility/ECLSaturationFunc.hpp @@ -85,19 +85,12 @@ namespace Opm { Vertical = 1 << 1u, }; - using IEB = SatFunc::EPSEvalInterface::InvalidEndpointBehaviour; - SatFuncScaling() - : enable (Type::Horizontal | Type::Vertical) - , invalid(IEB::UseUnscaled) + : enable(Type::Horizontal | Type::Vertical) {} // Default: Use both Horizontal and Vertical if specified. unsigned char enable; - - // Default: Use unscaled values in case of invalid scaled - // saturations occurring in the result set. - IEB invalid; }; /// Constructor diff --git a/ThirdParty/custom-opm-flowdiag-app/opmCore/opm/parser/eclipse/Units/Units.hpp b/ThirdParty/custom-opm-flowdiag-app/opmCore/opm/parser/eclipse/Units/Units.hpp deleted file mode 100644 index cc1dde7e66..0000000000 --- a/ThirdParty/custom-opm-flowdiag-app/opmCore/opm/parser/eclipse/Units/Units.hpp +++ /dev/null @@ -1,326 +0,0 @@ -//=========================================================================== -// -// File: Units.hpp -// -// Created: Thu Jul 2 09:19:08 2009 -// -// Author(s): Halvor M Nilsen -// -// $Date$ -// -// $Revision$ -// -//=========================================================================== - -/* -Copyright 2009, 2010, 2011, 2012 SINTEF ICT, Applied Mathematics. -Copyright 2009, 2010, 2011, 2012 Statoil ASA. - -This file is part of the Open Porous Media project (OPM). - -OPM is free software: you can redistribute it and/or modify -it under the terms of the GNU General Public License as published by -the Free Software Foundation, either version 3 of the License, or -(at your option) any later version. - -OPM is distributed in the hope that it will be useful, -but WITHOUT ANY WARRANTY; without even the implied warranty of -MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -GNU General Public License for more details. - -You should have received a copy of the GNU General Public License -along with OPM. If not, see . -*/ - -#ifndef OPM_UNITS_HEADER -#define OPM_UNITS_HEADER - -/** -The unit sets emplyed in ECLIPSE, in particular the FIELD units, -are quite inconsistent. Ideally one should choose units for a set -of base quantities like Mass,Time and Length and then derive the -units for e.g. pressure and flowrate in a consisten manner. However -that is not the case; for instance in the metric system we have: - -[Length] = meters -[time] = days -[mass] = kg - -This should give: - -[Pressure] = [mass] / ([length] * [time]^2) = kg / (m * days * days) - -Instead pressure is given in Bars. When it comes to FIELD units the -number of such examples is long. -*/ -namespace Opm { -namespace prefix - /// Conversion prefix for units. -{ -constexpr const double micro = 1.0e-6; /**< Unit prefix [\f$\mu\f$] */ -constexpr const double milli = 1.0e-3; /**< Unit prefix [m] */ -constexpr const double centi = 1.0e-2; /**< Non-standard unit prefix [c] */ -constexpr const double deci = 1.0e-1; /**< Non-standard unit prefix [d] */ -constexpr const double kilo = 1.0e3; /**< Unit prefix [k] */ -constexpr const double mega = 1.0e6; /**< Unit prefix [M] */ -constexpr const double giga = 1.0e9; /**< Unit prefix [G] */ -} // namespace prefix - -namespace unit - /// Definition of various units. - /// All the units are defined in terms of international standard - /// units (SI). Example of use: We define a variable \c k which - /// gives a permeability. We want to set \c k to \f$1\,mD\f$. - /// \code - /// using namespace Opm::unit - /// double k = 0.001*darcy; - /// \endcode - /// We can also use one of the prefixes defined in Opm::prefix - /// \code - /// using namespace Opm::unit - /// using namespace Opm::prefix - /// double k = 1.0*milli*darcy; - /// \endcode -{ -///\name Common powers -/// @{ -constexpr double square(double v) { return v * v; } -constexpr double cubic (double v) { return v * v * v; } -/// @} - -// -------------------------------------------------------------- -// Basic (fundamental) units and conversions -// -------------------------------------------------------------- - -/// \name Length -/// @{ -constexpr const double meter = 1; -constexpr const double inch = 2.54 * prefix::centi*meter; -constexpr const double feet = 12 * inch; -/// @} - -/// \name Time -/// @{ -constexpr const double second = 1; -constexpr const double minute = 60 * second; -constexpr const double hour = 60 * minute; -constexpr const double day = 24 * hour; -constexpr const double year = 365 * day; -/// @} - -/// \name Volume -/// @{ -constexpr const double gallon = 231 * cubic(inch); -constexpr const double stb = 42 * gallon; -constexpr const double liter = 1 * cubic(prefix::deci*meter); -/// @} - -/// \name Mass -/// @{ -constexpr const double kilogram = 1; -constexpr const double gram = 1.0e-3 * kilogram; -// http://en.wikipedia.org/wiki/Pound_(mass)#Avoirdupois_pound -constexpr const double pound = 0.45359237 * kilogram; -/// @} - -// -------------------------------------------------------------- -// Standardised constants -// -------------------------------------------------------------- - -/// \name Standardised constant -/// @{ -constexpr const double gravity = 9.80665 * meter/square(second); -/// @} - -// -------------------------------------------------------------- -// Derived units and conversions -// -------------------------------------------------------------- - -/// \name Force -/// @{ -constexpr const double Newton = kilogram*meter / square(second); // == 1 -constexpr const double dyne = 1e-5*Newton; -constexpr const double lbf = pound * gravity; // Pound-force - /// @} - - /// \name Pressure - /// @{ -constexpr const double Pascal = Newton / square(meter); // == 1 -constexpr const double barsa = 100000 * Pascal; -constexpr const double atm = 101325 * Pascal; -constexpr const double psia = lbf / square(inch); -/// @} - -/// \name Temperature. This one is more complicated -/// because the unit systems used by Eclipse (i.e. degrees -/// Celsius and degrees Fahrenheit require to add or -/// subtract an offset for the conversion between from/to -/// Kelvin -/// @{ -constexpr const double degCelsius = 1.0; // scaling factor °C -> K -constexpr const double degCelsiusOffset = 273.15; // offset for the °C -> K conversion - -constexpr const double degFahrenheit = 5.0/9; // scaling factor °F -> K -constexpr const double degFahrenheitOffset = 255.37; // offset for the °C -> K conversion - /// @} - - /// \name Viscosity - /// @{ -constexpr const double Pas = Pascal * second; // == 1 -constexpr const double Poise = prefix::deci*Pas; -/// @} - -namespace perm_details { -constexpr const double p_grad = atm / (prefix::centi*meter); -constexpr const double area = square(prefix::centi*meter); -constexpr const double flux = cubic (prefix::centi*meter) / second; -constexpr const double velocity = flux / area; -constexpr const double visc = prefix::centi*Poise; -constexpr const double darcy = (velocity * visc) / p_grad; -// == 1e-7 [m^2] / 101325 -// == 9.869232667160130e-13 [m^2] -} -/// \name Permeability -/// @{ -/// -/// A porous medium with a permeability of 1 darcy permits a flow (flux) -/// of \f$1\,\mathit{cm}^3/s\f$ of a fluid with viscosity -/// \f$1\,\mathit{cP}\f$ (\f$1\,mPa\cdot s\f$) under a pressure gradient -/// of \f$1\,\mathit{atm}/\mathit{cm}\f$ acting across an area of -/// \f$1\,\mathit{cm}^2\f$. -/// -constexpr const double darcy = perm_details::darcy; -/// @} - -/** -* Unit conversion routines. -*/ -namespace convert { -/** -* Convert from external units of measurements to equivalent -* internal units of measurements. Note: The internal units of -* measurements are *ALWAYS*, and exclusively, SI. -* -* Example: Convert a double @c kx, containing a permeability value -* in units of milli-darcy (mD) to the equivalent value in SI units -* (i.e., \f$m^2\f$). -* \code -* using namespace Opm::unit; -* using namespace Opm::prefix; -* convert::from(kx, milli*darcy); -* \endcode -* -* @param[in] q Physical quantity. -* @param[in] unit Physical unit of measurement. -* @return Value of @c q in equivalent SI units of measurements. -*/ -constexpr double from(const double q, const double unit) -{ - return q * unit; -} - -/** -* Convert from internal units of measurements to equivalent -* external units of measurements. Note: The internal units of -* measurements are *ALWAYS*, and exclusively, SI. -* -* Example: Convert a std::vector p, containing -* pressure values in the SI unit Pascal (i.e., unit::Pascal) to the -* equivalent values in Psi (unit::psia). -* \code -* using namespace Opm::unit; -* std::transform(p.begin(), p.end(), p.begin(), -* boost::bind(convert::to, _1, psia)); -* \endcode -* -* @param[in] q Physical quantity, measured in SI units. -* @param[in] unit Physical unit of measurement. -* @return Value of @c q in unit unit. -*/ -constexpr double to(const double q, const double unit) -{ - return q / unit; -} -} // namespace convert -} - -namespace Metric { -using namespace prefix; -using namespace unit; -constexpr const double Pressure = barsa; -constexpr const double Temperature = degCelsius; -constexpr const double TemperatureOffset = degCelsiusOffset; -constexpr const double AbsoluteTemperature = degCelsius; // actually [K], but the these two are identical -constexpr const double Length = meter; -constexpr const double Time = day; -constexpr const double Mass = kilogram; -constexpr const double Permeability = milli*darcy; -constexpr const double Transmissibility = centi*Poise*cubic(meter)/(day*barsa); -constexpr const double LiquidSurfaceVolume = cubic(meter); -constexpr const double GasSurfaceVolume = cubic(meter); -constexpr const double ReservoirVolume = cubic(meter); -constexpr const double GasDissolutionFactor = GasSurfaceVolume/LiquidSurfaceVolume; -constexpr const double OilDissolutionFactor = LiquidSurfaceVolume/GasSurfaceVolume; -constexpr const double Density = kilogram/cubic(meter); -constexpr const double PolymerDensity = kilogram/cubic(meter); -constexpr const double Salinity = kilogram/cubic(meter); -constexpr const double Viscosity = centi*Poise; -constexpr const double Timestep = day; -constexpr const double SurfaceTension = dyne/(centi*meter); -} - - -namespace Field { -using namespace prefix; -using namespace unit; -constexpr const double Pressure = psia; -constexpr const double Temperature = degFahrenheit; -constexpr const double TemperatureOffset = degFahrenheitOffset; -constexpr const double AbsoluteTemperature = degFahrenheit; // actually [°R], but the these two are identical -constexpr const double Length = feet; -constexpr const double Time = day; -constexpr const double Mass = pound; -constexpr const double Permeability = milli*darcy; -constexpr const double Transmissibility = centi*Poise*stb/(day*psia); -constexpr const double LiquidSurfaceVolume = stb; -constexpr const double GasSurfaceVolume = 1000*cubic(feet); -constexpr const double ReservoirVolume = stb; -constexpr const double GasDissolutionFactor = GasSurfaceVolume/LiquidSurfaceVolume; -constexpr const double OilDissolutionFactor = LiquidSurfaceVolume/GasSurfaceVolume; -constexpr const double Density = pound/cubic(feet); -constexpr const double PolymerDensity = pound/stb; -constexpr const double Salinity = pound/stb; -constexpr const double Viscosity = centi*Poise; -constexpr const double Timestep = day; -constexpr const double SurfaceTension = dyne/(centi*meter); -} - - -namespace Lab { -using namespace prefix; -using namespace unit; -constexpr const double Pressure = atm; -constexpr const double Temperature = degCelsius; -constexpr const double TemperatureOffset = degCelsiusOffset; -constexpr const double AbsoluteTemperature = degCelsius; // actually [K], but the these two are identical -constexpr const double Length = centi*meter; -constexpr const double Time = hour; -constexpr const double Mass = gram; -constexpr const double Permeability = milli*darcy; -constexpr const double Transmissibility = centi*Poise*cubic(centi*meter)/(hour*atm); -constexpr const double LiquidSurfaceVolume = cubic(centi*meter); -constexpr const double GasSurfaceVolume = cubic(centi*meter); -constexpr const double ReservoirVolume = cubic(centi*meter); -constexpr const double GasDissolutionFactor = GasSurfaceVolume/LiquidSurfaceVolume; -constexpr const double OilDissolutionFactor = LiquidSurfaceVolume/GasSurfaceVolume; -constexpr const double Density = gram/cubic(centi*meter); -constexpr const double PolymerDensity = gram/cubic(centi*meter); -constexpr const double Salinity = gram/cubic(centi*meter); -constexpr const double Viscosity = centi*Poise; -constexpr const double Timestep = hour; -constexpr const double SurfaceTension = dyne/(centi*meter); -} - -} - -#endif // OPM_UNITS_HEADER