From 4ec2e65db33a396a85344285782a0778491cb71e Mon Sep 17 00:00:00 2001 From: Vegard Kippe Date: Mon, 4 Nov 2024 20:18:06 +0100 Subject: [PATCH] Add support for output of relative permeabilities to the .INIT file for SLGOF input. --- .../EclipseState/Tables/TableColumn.cpp | 8 + .../EclipseState/Tables/TableColumn.hpp | 2 + opm/output/eclipse/Tables.cpp | 384 +++++++++++++++++- 3 files changed, 382 insertions(+), 12 deletions(-) diff --git a/opm/input/eclipse/EclipseState/Tables/TableColumn.cpp b/opm/input/eclipse/EclipseState/Tables/TableColumn.cpp index d31f70ad2b5..cfd756cece6 100644 --- a/opm/input/eclipse/EclipseState/Tables/TableColumn.cpp +++ b/opm/input/eclipse/EclipseState/Tables/TableColumn.cpp @@ -253,6 +253,14 @@ namespace Opm { return m_values.end(); } + std::vector::const_reverse_iterator TableColumn::rbegin() const { + return m_values.rbegin(); + } + + std::vector::const_reverse_iterator TableColumn::rend() const { + return m_values.rend(); + } + bool TableColumn::hasDefault( ) const { if (m_defaultCount > 0) diff --git a/opm/input/eclipse/EclipseState/Tables/TableColumn.hpp b/opm/input/eclipse/EclipseState/Tables/TableColumn.hpp index 1e9c6e7e624..38f8aa662de 100644 --- a/opm/input/eclipse/EclipseState/Tables/TableColumn.hpp +++ b/opm/input/eclipse/EclipseState/Tables/TableColumn.hpp @@ -67,6 +67,8 @@ namespace Opm { std::vector vectorCopy() const; std::vector::const_iterator begin() const; std::vector::const_iterator end() const; + std::vector::const_reverse_iterator rbegin() const; + std::vector::const_reverse_iterator rend() const; bool operator==(const TableColumn& data) const; diff --git a/opm/output/eclipse/Tables.cpp b/opm/output/eclipse/Tables.cpp index a807c4106c9..64109085ae5 100644 --- a/opm/output/eclipse/Tables.cpp +++ b/opm/output/eclipse/Tables.cpp @@ -32,6 +32,7 @@ #include #include #include +#include #include #include #include @@ -1014,11 +1015,100 @@ namespace { namespace SatFunc { }; } + /// Create linearised and padded 'TAB' vector entries of normalised + /// SGFN tables for all saturation function regions from Family One + /// table data (SLGOF keyword). + /// + /// \param[in] numRows Number of rows to allocate in the output + /// vector for each table. Expected to be equal to the number of + /// declared saturation nodes in the simulation run's TABDIMS + /// keyword (Item 3). + /// + /// \param tolcrit Minimum relative permeability threshold value for + /// phase to be considered mobile. Values less than this threshold + /// are output as zero. + /// + /// \param[in] units Active unit system. Needed to convert SI + /// convention capillary pressure values (Pascal) to declared + /// conventions of the run specification. + /// + /// \param[in] slgof Collection of SLGOF tables for all saturation + /// regions. + /// + /// \return Descriptor and linearised/normalised/padded 'TAB' vector + /// values for output SGFN tables. Corresponds to unit-converted + /// copies of columns 1, 2, and 4--with added derivatives--of the + /// input SLGOF tables. + std::tuple> + fromSLGOF(const std::size_t numRows, + const double tolcrit, + const Opm::UnitSystem& units, + const Opm::TableContainer& slgof) + { + using SLGOF = ::Opm::SlgofTable; + + const auto tableSize = detail::allocatedRows(numRows, slgof); + + const auto numTab = slgof.size(); + const auto numDep = std::size_t{2}; // Krg, Pcgo + + return { + static_cast(tableSize), + static_cast(numTab), + detail::createSatfuncTable(numTab, tableSize, numDep, + [tolcrit, &units, &slgof](const std::size_t tableID, + const std::size_t primID, + Opm::LinearisedOutputTable& linTable) + { + const auto& t = slgof.getTable(tableID); + + auto numActRows = std::size_t{0}; + + /// Reverse order in all columns, since Sl is increasing (Sg descreasing) + // Sg = 1.0 - Sl + { + const auto& Sl = t.getSlColumn(); + + numActRows = Sl.size(); + std::vector Sg(std::rbegin(Sl), std::rend(Sl)); + std::transform(Sg.begin(), Sg.end(), Sg.begin(), + [](double x){ return (1.0 - x); }); + std::copy(std::begin(Sg), std::end(Sg), + linTable.column(tableID, primID, 0)); + } + + // Krg(Sg) + const auto& Krg_rev = t.getKrgColumn(); + std::vector Krg(Krg_rev.rbegin(), Krg_rev.rend()); + detail::outputRelperm(Krg, tolcrit, + linTable.column(tableID, primID, 1)); + + // Pcgo(Sg) + { + constexpr auto uPress = ::Opm::UnitSystem::measure::pressure; + + const auto& pc = t.getPcogColumn(); + std::transform(std::rbegin(pc), std::rend(pc), + linTable.column(tableID, primID, 2), + [&units](const double Pc) -> double + { + return units.from_si(uPress, Pc); + }); + } + + // Inform createSatfuncTable() of number of active rows + // in this table. Needed to compute slopes of piecewise + // linear interpolants. + return numActRows; + }) + }; + } + /// Create linearised and padded 'TAB' vector entries of normalised /// SGFN tables for all saturation function regions from Family One /// table data. /// - /// Distinguishes between keywords SGOF and SGOFLET. + /// Distinguishes between keywords SGOF, SLGOF and SGOFLET. /// /// \param[in] numRows Number of rows to allocate in the output /// vector for each table. Expected to be equal to the number of @@ -1047,7 +1137,15 @@ namespace { namespace SatFunc { if (letTables.empty()) { // Typical case: Not LET-based saturation functions. - return fromSGOF(numRows, tolcrit, units, tabMgr.getSgofTables()); + const auto& sgof = tabMgr.getSgofTables(); + if (sgof.size() > 0) { + return fromSGOF(numRows, tolcrit, units, sgof); + } + const auto& slgof = tabMgr.getSlgofTables(); + if (slgof.size() > 0) { + return fromSLGOF(numRows, tolcrit, units, slgof); + } + // Warn here? } // LET-based saturation functions for gas. @@ -1353,6 +1451,94 @@ namespace { namespace SatFunc { }; } + + /// Create linearised and padded 'TAB' vector entries of + /// normalised two-phase SOFN tables for all saturation function + /// regions from Family One table data (SLGOF keyword--G/O System). + /// + /// \param[in] numRows Number of rows to allocate in the output + /// vector for each table. Expected to be equal to the + /// number of declared saturation nodes in the simulation + /// run's TABDIMS keyword (Item 3). + /// + /// \param tolcrit Minimum relative permeability threshold value + /// for phase to be considered mobile. Values less than this + /// threshold are output as zero. + /// + /// \param[in] slgof Collection of SLGOF tables for all saturation + /// regions. + /// + /// \return Descriptor and linearised/normalised/padded 'TAB' + /// vector values for output SOFN tables. Corresponds to + /// translated (1-Sg), reverse saturation column (column 1) + /// and reverse column of relative permeability for oil + /// (column 3) from the input SLGOF table--with added + /// derivatives. + std::tuple> + fromSLGOF(const std::size_t numRows, + const double tolcrit, + const Opm::TableContainer& slgof) + { + using SLGOF = ::Opm::SlgofTable; + + const auto tableSize = detail::allocatedRows(numRows, slgof); + + const auto numTab = slgof.size(); + const auto numDep = std::size_t{1}; // Kro + + return { + static_cast(tableSize), + static_cast(numTab), + detail::createSatfuncTable(numTab, tableSize, numDep, + [tolcrit, &slgof](const std::size_t tableID, + const std::size_t primID, + Opm::LinearisedOutputTable& linTable) + { + const auto& t = slgof.getTable(tableID); + + auto numActRows = std::size_t{0}; + + // So = Sl (two-phase) + { + const auto& Sl = t.getSlColumn(); + numActRows = Sl.size(); + std::copy(Sl.begin(), Sl.end(), linTable.column(tableID, primID, 0)); + + // auto So = std::vector{}; + // So.reserve(numActRows); + + // // Two-phase system => So = 1-Sg + // std::transform(std::begin(Sg), std::end(Sg), + // std::back_inserter(So), + // [](const double sg) -> double + // { + // return 1.0 - sg; + // }); + + // std::copy(So.rbegin(), So.rend(), + // linTable.column(tableID, primID, 0)); + } + + // Kro(So) + { + const auto& kr = t.getKrogColumn(); + + const auto krog = std::vector { + std::begin(kr), std::end(kr) + }; + + detail::outputRelperm(krog.begin(), krog.end(), tolcrit, + linTable.column(tableID, primID, 1)); + } + + // Inform createSatfuncTable() of number of active + // rows in this table. Needed to compute slopes of + // piecewise linear interpolants. + return numActRows; + }) + }; + } + /// Create linearised and padded 'TAB' vector entries of /// normalised SOFN tables for all saturation function regions /// from Family One table data in the case of an oil/gas @@ -1379,7 +1565,15 @@ namespace { namespace SatFunc { const auto& gasLET = tabMgr.getSgofletTable(); if (gasLET.empty()) { - return fromSGOF(numRows, tolcrit, tabMgr.getSgofTables()); + const auto& sgof = tabMgr.getSgofTables(); + if (sgof.size() > 0) { + return fromSGOF(numRows, tolcrit, sgof); + } + const auto& slgof = tabMgr.getSlgofTables(); + if (sgof.size() > 0) { + return fromSLGOF(numRows, tolcrit, slgof); + } + // Warn here? } return { @@ -1697,7 +1891,7 @@ namespace { namespace SatFunc { // tables having sorted phase saturation values (required by ECL // format). std::vector - mergeTables(const std::vector& t) + mergeTables(const std::vector& t, const double tolerance = 1.0e-12) { auto ret = std::vector{}; @@ -1709,10 +1903,11 @@ namespace { namespace SatFunc { std::set_union(std::begin(t0), std::end(t0), std::begin(t1), std::end(t1), std::back_inserter(ret), - [&t](const TableElement& e1, const TableElement& e2) + [&t,tolerance](const TableElement& e1, const TableElement& e2) { - return t[e1.function].So(e1.index) - < t[e2.function].So(e2.index); + const double val1 = t[e1.function].So(e1.index); + const double val2 = t[e2.function].So(e2.index); + return ( (val1 + tolerance) < val2); }); return ret; @@ -1780,6 +1975,74 @@ namespace { namespace SatFunc { return ret; } + // Create collection of individual columns of single SOF3 table + // through joining input SLGOF and SWOF tables on increasing oil + // saturation and appropriate KrOX columns. + std::array, 3> + makeSOF3Table(const Opm::SlgofTable& slgof, + const Opm::SwofTable& swof) + { + auto ret = std::array, 3>{}; + + auto tbl = std::vector{}; + tbl.reserve(2); + + // Note: Order between Krow(So) and Krog(So) matters + // here. This order must match the expected column + // order in SOF3--i.e. [ So, Krow, Krog ]. + + // 1) Krow(So) + { + // So = 1.0 - Sw + const auto& Sw = swof.getSwColumn(); + const auto& Krow = swof.getKrowColumn(); + const auto& So_off = 1.0; + + tbl.emplace_back(Sw .vectorCopy(), + Krow.vectorCopy(), So_off); + } + + // 2) Krog(So) + { + // Map Krog(Sl) [decreasing Sg] to Krog(Sg) [increasing Sg] + const auto swco = swof.getSwColumn()[0]; + const auto& Sl = slgof.getSlColumn(); + const auto& Krog = slgof.getKrogColumn(); + const auto So_off = + 1.0 - swco; + + // Sg = 1.0 - Sl; + std::vector sg( Sl.rbegin(), Sl.rend() ); + std::transform(sg.begin(), sg.end(), sg.begin(), + [](double x){ return (1.0 - x); }); + + std::vector krog ( Krog.rbegin(), Krog.rend() ); + tbl.emplace_back(sg, krog, So_off); + } + + const auto mrg = mergeTables(tbl); + + for (auto& col : ret) { col.reserve(mrg.size()); } + + for (const auto& row : mrg) { + const auto self = row.function; + const auto other = 1 - row.function; + + // 1) Assign So + ret[0].push_back(tbl[self].So(row.index)); + + // 2) Assign Kro for "self" column (the one that got + // picked for this row). + ret[1 + self].push_back(tbl[self].Kro(row.index)); + + // 3) Assign Kro for "other" column (the one that + // did not get picked for this row). + ret[1 + other].push_back(tbl[other].Kro(ret[0].back())); + } + + return ret; + } + /// Create linearised and padded 'TAB' vector entries of /// normalised three-phase SOFN tables for all saturation /// function regions from Family One table data. @@ -1865,12 +2128,99 @@ namespace { namespace SatFunc { }; } + + /// Create linearised and padded 'TAB' vector entries of + /// normalised three-phase SOFN tables for all saturation + /// function regions from Family One table data. + /// + /// \param[in] numRows Number of rows to allocate in the output + /// vector for each table. Expected to be twice the number + /// of declared saturation nodes in the simulation run's + /// TABDIMS keyword (Item 3). + /// + /// \param tolcrit Minimum relative permeability threshold value + /// for phase to be considered mobile. Values less than this + /// threshold are output as zero. + /// + /// \param[in] slgof Collection of SLGOF tables for all saturation + /// regions. + /// + /// \param[in] swof Collection of SWOF tables for all saturation + /// regions. + /// + /// \return Linearised and padded 'TAB' vector values for + /// three-phase SOFN tables. Corresponds to column 1 from + /// both of the input SLGOF and SWOF tables, as well as column + /// 3 from the input SWOF table and column 3 from the input + /// SGOF table--expanded so as to have values for all oil + /// saturation nodes. Derivatives added in columns 4 and 5. + std::tuple> + fromSLGOFandSWOF(const std::size_t numRows, + const double tolcrit, + const Opm::TableContainer& slgof, + const Opm::TableContainer& swof) + { + using SLGOF = ::Opm::SlgofTable; + using SWOF = ::Opm::SwofTable; + + const auto tableSizeGas = detail::allocatedRows(numRows, slgof); + const auto tableSizeWat = detail::allocatedRows(numRows, swof); + + // Merged table has 2*NSSFUN rows to accommodate both the + // gas and the water table. + const auto tableSize = + 2 * std::max(tableSizeGas, tableSizeWat); + + const auto numTab = slgof.size(); + const auto numDep = std::size_t{2}; // Krow, Krog + + return { + static_cast(tableSize), + static_cast(numTab), + detail::createSatfuncTable(numTab, tableSize, numDep, + [tolcrit, &slgof, &swof] + (const std::size_t tableID, + const std::size_t primID, + Opm::LinearisedOutputTable& linTable) + { + const auto sof3 = + makeSOF3Table(slgof.getTable(tableID), + swof.getTable(tableID)); + + auto numActRows = std::size_t{0}; + + // So + { + const auto& So = sof3[0]; + + numActRows = So.size(); + std::copy(std::begin(So), std::end(So), + linTable.column(tableID, primID, 0)); + } + + // Krow(So) + detail::outputRelperm(sof3[1], tolcrit, + linTable.column(tableID, primID, 1)); + + // Krog(So) + detail::outputRelperm(sof3[2], tolcrit, + linTable.column(tableID, primID, 2)); + + // Inform createSatfuncTable() of number of active + // rows in this table. Needed to compute slopes of + // piecewise linear interpolants. + return numActRows; + }) + }; + } + + /// Create linearised and padded 'TAB' vector entries of /// normalised SOFN tables for all saturation function regions /// from Family One table data for a three-phase oil/gas/water /// system. /// - /// Distinguishes between the SGOF/SWOF and SGOFLET/SWOFLET + /// Distinguishes between the SGOF/SWOF, SLGOF/SWOF and SGOFLET/SWOFLET /// family of keywords. /// /// \param[in] numRows Number of rows to allocate in the output @@ -1896,10 +2246,20 @@ namespace { namespace SatFunc { if (gasLET.empty() || watLET.empty()) { // Common case: Saturation functions for oil entered in - // S[GW]OF tables. - return fromSGOFandSWOF(numRows, tolcrit, - tabMgr.getSgofTables(), - tabMgr.getSwofTables()); + // S[L][GW]OF tables. + const auto& sgof = tabMgr.getSgofTables(); + if (sgof.size() > 0) { + return fromSGOFandSWOF(numRows, tolcrit, + sgof, + tabMgr.getSwofTables()); + } + const auto& slgof = tabMgr.getSlgofTables(); + if (slgof.size() > 0) { + return fromSLGOFandSWOF(numRows, tolcrit, + slgof, + tabMgr.getSwofTables()); + } + // Warn here? } // Saturation functions for oil entered in S[GW]OFLET tables.