diff --git a/CMakeLists_files.cmake b/CMakeLists_files.cmake index 327f93c2aee..8fca374eea7 100644 --- a/CMakeLists_files.cmake +++ b/CMakeLists_files.cmake @@ -1050,6 +1050,8 @@ if(ENABLE_ECL_INPUT) opm/input/eclipse/EclipseState/Tables/PvdoTable.hpp opm/input/eclipse/EclipseState/Tables/OilvisctTable.hpp opm/input/eclipse/EclipseState/Tables/SgfnTable.hpp + opm/input/eclipse/EclipseState/Tables/WsfTable.hpp + opm/input/eclipse/EclipseState/Tables/GsfTable.hpp opm/input/eclipse/EclipseState/Tables/MiscTable.hpp opm/input/eclipse/EclipseState/Tables/SgwfnTable.hpp opm/input/eclipse/EclipseState/Tables/PvdsTable.hpp diff --git a/opm/input/eclipse/EclipseState/Runspec.hpp b/opm/input/eclipse/EclipseState/Runspec.hpp index e0452850b3f..c8ac14a78cf 100644 --- a/opm/input/eclipse/EclipseState/Runspec.hpp +++ b/opm/input/eclipse/EclipseState/Runspec.hpp @@ -330,6 +330,7 @@ class SatFuncControls { enum class KeywordFamily { Family_I, // SGOF, SWOF, SLGOF Family_II, // SGFN, SOF{2,3}, SWFN + Family_III, // GSF, WSF Undefined, }; diff --git a/opm/input/eclipse/EclipseState/Tables/GsfTable.hpp b/opm/input/eclipse/EclipseState/Tables/GsfTable.hpp new file mode 100644 index 00000000000..e5608680b51 --- /dev/null +++ b/opm/input/eclipse/EclipseState/Tables/GsfTable.hpp @@ -0,0 +1,41 @@ +/* + Copyright (C) 2023 Equinor AS + + 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_PARSER_GSF_TABLE_HPP +#define OPM_PARSER_GSF_TABLE_HPP + +#include "SimpleTable.hpp" + +namespace Opm { + + class DeckItem; + + class GsfTable : public SimpleTable { + + public: + GsfTable( const DeckItem& item, const int tableID ); + + const TableColumn& getSgColumn() const; + const TableColumn& getKrgColumn() const; + const TableColumn& getPcgwColumn() const; + + }; +} + +#endif + diff --git a/opm/input/eclipse/EclipseState/Tables/TableManager.hpp b/opm/input/eclipse/EclipseState/Tables/TableManager.hpp index 819a2b3a42c..f88d81fe5d0 100644 --- a/opm/input/eclipse/EclipseState/Tables/TableManager.hpp +++ b/opm/input/eclipse/EclipseState/Tables/TableManager.hpp @@ -132,6 +132,9 @@ namespace Opm { const TableContainer& getMsfnTables() const; const TableContainer& getTlpmixpaTables() const; + const TableContainer& getWsfTables() const; + const TableContainer& getGsfTables() const; + const JFunc& getJFunc() const; const std::vector& getPvtgTables() const; diff --git a/opm/input/eclipse/EclipseState/Tables/WsfTable.hpp b/opm/input/eclipse/EclipseState/Tables/WsfTable.hpp new file mode 100644 index 00000000000..f9d2d3b6881 --- /dev/null +++ b/opm/input/eclipse/EclipseState/Tables/WsfTable.hpp @@ -0,0 +1,39 @@ +/* + Copyright (C) 2023 Equinor AS + + 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_PARSER_WSF_TABLE_HPP +#define OPM_PARSER_WSF_TABLE_HPP + +#include "SimpleTable.hpp" + +namespace Opm { + + class DeckItem; + + class WsfTable : public SimpleTable { + + public: + WsfTable( const DeckItem& item, const int tableID ); + + const TableColumn& getSwColumn() const; + const TableColumn& getKrwColumn() const; + }; +} + +#endif + diff --git a/opm/output/eclipse/Tables.hpp b/opm/output/eclipse/Tables.hpp index ee4fde576bb..5322162570e 100644 --- a/opm/output/eclipse/Tables.hpp +++ b/opm/output/eclipse/Tables.hpp @@ -105,6 +105,15 @@ namespace Opm { const bool oil, const bool wat); + /// Add saturation function tables corresponding to family III (WSF, + /// GSF) to the tabular data (TABDIMS and TAB vectors). + /// + /// \param[in] es Valid \c EclipseState object with accurate table + /// dimensions ("TABDIMS" keyword) and an initialised \c + /// TableManager sub-object. + /// + void addSatFunc_FamilyThree(const EclipseState& es); + /// Add gas PVT tables (keywords PVDG and PVTG) to the tabular data /// (TABDIMS and TAB vectors). /// diff --git a/src/opm/input/eclipse/EclipseState/Grid/SatfuncPropertyInitializers.cpp b/src/opm/input/eclipse/EclipseState/Grid/SatfuncPropertyInitializers.cpp index f32608347fb..463dfc7ecda 100644 --- a/src/opm/input/eclipse/EclipseState/Grid/SatfuncPropertyInitializers.cpp +++ b/src/opm/input/eclipse/EclipseState/Grid/SatfuncPropertyInitializers.cpp @@ -27,6 +27,8 @@ #include #include #include +#include +#include #include #include #include @@ -72,7 +74,7 @@ namespace { * If SWFN, SGFN and SOF3 are specified in the deck it return II * If keywords are missing or mixed, an error is given. */ - enum class SatfuncFamily { none = 0, I = 1, II = 2 }; + enum class SatfuncFamily { none = 0, I = 1, II = 2, III = 3 }; SatfuncFamily getSaturationFunctionFamily(const Opm::TableManager& tm, @@ -101,6 +103,10 @@ namespace { (twoP && tm.hasTables("SOF2")))) || (wat && tm.hasTables("SWFN")); + const auto family3 = //WSF, GSF gas-water CO2STORE case + tm.hasTables("GSF") && + tm.hasTables("WSF"); + if (gas && tm.hasTables("SGOF") && tm.hasTables("SLGOF")) { throw std::invalid_argument("Both SGOF and SLGOF have been specified but these tables are mutually exclusive!"); } @@ -110,12 +116,19 @@ namespace { "Use either SGOF (or SLGOF) and/or SWOF or SGFN/SWFN and SOF2/SOF3"); } + if (family3 && !(gas && wat)) { + throw std::invalid_argument("Saturation family 3 should only be used for GAS and WATER case"); + } + if (family1) return SatfuncFamily::I; if (family2) return SatfuncFamily::II; + if (family3) + return SatfuncFamily::III; + throw std::invalid_argument { "Saturation functions must be specified using " "either family 1 or family 2 keywords\n" @@ -135,6 +148,7 @@ namespace { const auto& swofTables = tm.getSwofTables(); const auto& swofLetTables = tm.getSwofletTable(); const auto& swfnTables = tm.getSwfnTables(); + const auto& wsfTables = tm.getWsfTables(); const auto famI = [&swofTables]( int i ) { return swofTables.getTable( i ).getSwColumn().front(); @@ -147,7 +161,9 @@ namespace { const auto famII = [&swfnTables]( int i ) { return swfnTables.getTable( i ).getSwColumn().front(); }; - + const auto famIII = [&wsfTables]( int i ) { + return wsfTables.getTable( i ).getSwColumn().front(); + }; switch( getSaturationFunctionFamily( tm, ph ) ) { case SatfuncFamily::I: if( !swofTables.empty() ) @@ -158,6 +174,7 @@ namespace { throw std::domain_error("Either SWOF or SWOFLET tables must be provided"); case SatfuncFamily::II: return map( famII, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: return map( famIII, Opm::fun::iota( num_tables ) ); default: throw std::domain_error("No valid saturation keyword family specified"); } @@ -175,6 +192,7 @@ namespace { const auto& swofTables = tm.getSwofTables(); const auto& swofLetTables = tm.getSwofletTable(); const auto& swfnTables = tm.getSwfnTables(); + const auto& wsfTables = tm.getWsfTables(); const auto famI = [&swofTables]( int i ) { return swofTables.getTable( i ).getSwColumn().back(); @@ -187,7 +205,9 @@ namespace { const auto famII = [&swfnTables]( int i ) { return swfnTables.getTable( i ).getSwColumn().back(); }; - + const auto famIII = [&wsfTables]( int i ) { + return wsfTables.getTable( i ).getSwColumn().back(); + }; switch( getSaturationFunctionFamily( tm, ph ) ) { case SatfuncFamily::I: if( !swofTables.empty() ) @@ -198,6 +218,8 @@ namespace { throw std::domain_error("Either SWOF or SWOFLET tables must be provided"); case SatfuncFamily::II: return map( famII, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: return map( famIII, Opm::fun::iota( num_tables ) ); + default: throw std::domain_error("No valid saturation keyword family specified"); } @@ -216,6 +238,7 @@ namespace { const auto& sgofLetTables = tm.getSgofletTable(); const auto& slgofTables = tm.getSlgofTables(); const auto& sgfnTables = tm.getSgfnTables(); + const auto& gsfTables = tm.getGsfTables(); const auto famI_sgof = [&sgofTables]( int i ) { return sgofTables.getTable( i ).getSgColumn().front(); @@ -232,7 +255,9 @@ namespace { const auto famII = [&sgfnTables]( int i ) { return sgfnTables.getTable( i ).getSgColumn().front(); }; - + const auto famIII = [&gsfTables]( int i ) { + return gsfTables.getTable( i ).getSgColumn().front(); + }; switch( getSaturationFunctionFamily( tm, ph ) ) { case SatfuncFamily::I: if( sgofTables.empty() && sgofLetTables.empty() && slgofTables.empty() ) @@ -247,6 +272,8 @@ namespace { case SatfuncFamily::II: return Opm::fun::map( famII, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: + return Opm::fun::map( famIII, Opm::fun::iota( num_tables ) ); default: throw std::domain_error("No valid saturation keyword family specified"); @@ -267,6 +294,7 @@ namespace { const auto& swofLetTables = tm.getSwofletTable(); const auto& slgofTables = tm.getSlgofTables(); const auto& sgfnTables = tm.getSgfnTables(); + const auto& gsfTables = tm.getGsfTables(); const auto famI_sgof = [&sgofTables]( int i ) { return sgofTables.getTable( i ).getSgColumn().back(); @@ -284,6 +312,9 @@ namespace { const auto famII = [&sgfnTables]( int i ) { return sgfnTables.getTable( i ).getSgColumn().back(); }; + const auto famIII = [&gsfTables]( int i ) { + return gsfTables.getTable( i ).getSgColumn().back(); + }; switch( getSaturationFunctionFamily( tm, ph ) ) { case SatfuncFamily::I: @@ -299,6 +330,8 @@ namespace { case SatfuncFamily::II: return Opm::fun::map( famII, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: + return Opm::fun::map( famIII, Opm::fun::iota( num_tables ) ); default: throw std::domain_error("No valid saturation keyword family specified"); @@ -375,6 +408,7 @@ namespace { const auto& swofTables = tm.getSwofTables(); const auto& swofLetTables = tm.getSwofletTable(); const auto& swfnTables = tm.getSwfnTables(); + const auto& wsfTables = tm.getWsfTables(); const auto famI = [&swofTables, tolcrit](const int i) -> double { @@ -389,7 +423,10 @@ namespace { { return critical_water(swfnTables.getTable(i), tolcrit); }; - + const auto famIII = [&wsfTables, tolcrit](const int i) -> double + { + return critical_water(wsfTables.getTable(i), tolcrit); + }; switch( getSaturationFunctionFamily( tm, ph ) ) { case SatfuncFamily::I: if( !swofTables.empty() ) @@ -400,6 +437,8 @@ namespace { throw std::domain_error("Either SWOF or SWOFLET tables must be provided"); case SatfuncFamily::II: return Opm::fun::map( famII, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: return Opm::fun::map( famIII, Opm::fun::iota( num_tables ) ); + default: throw std::domain_error("No valid saturation keyword family specified"); } } @@ -448,6 +487,7 @@ namespace { const auto& sgofTables = tm.getSgofTables(); const auto& sgofLetTables = tm.getSgofletTable(); const auto& slgofTables = tm.getSlgofTables(); + const auto& gsfTables = tm.getGsfTables(); const auto famI_sgof = [&sgofTables, tolcrit](const int i) -> double { @@ -468,6 +508,11 @@ namespace { return critical_gas(sgfnTables.getTable(i), tolcrit); }; + const auto famIII = [&gsfTables, tolcrit](const int i) -> double + { + return critical_gas(gsfTables.getTable(i), tolcrit); + }; + switch( getSaturationFunctionFamily( tm, ph ) ) { case SatfuncFamily::I: if( sgofTables.empty() && sgofLetTables.empty() && slgofTables.empty() ) @@ -482,6 +527,8 @@ namespace { case SatfuncFamily::II: return Opm::fun::map( famII, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: + return Opm::fun::map( famIII, Opm::fun::iota( num_tables ) ); default: throw std::domain_error("No valid saturation keyword family specified"); @@ -580,7 +627,8 @@ namespace { return ph.active(::Opm::Phase::GAS) ? Opm::fun::map( famII_3p, Opm::fun::iota( num_tables ) ) : Opm::fun::map( famII_2p, Opm::fun::iota( num_tables ) ); - + case SatfuncFamily::III: + throw std::domain_error("Saturation keyword family III is not applicable for a oil system"); default: throw std::domain_error("No valid saturation keyword family specified"); } } @@ -674,7 +722,8 @@ namespace { return ph.active(::Opm::Phase::WATER) ? Opm::fun::map( famII_3p, Opm::fun::iota( num_tables ) ) : Opm::fun::map( famII_2p, Opm::fun::iota( num_tables ) ); - + case SatfuncFamily::III: + throw std::domain_error("Saturation keyword family III is not applicable for a oil system"); default: throw std::domain_error("No valid saturation keyword family specified"); } @@ -693,6 +742,7 @@ namespace { const auto& sgofLetTables = tm.getSgofletTable(); const auto& slgofTables = tm.getSlgofTables(); const auto& sgfnTables = tm.getSgfnTables(); + const auto& gsfTables = tm.getGsfTables(); const auto& famI_sgof = [&sgofTables]( int i ) { return sgofTables.getTable( i ).getKrgColumn().back(); @@ -710,6 +760,10 @@ namespace { return sgfnTables.getTable( i ).getKrgColumn().back(); }; + const auto& famIII = [&gsfTables]( int i ) { + return gsfTables.getTable( i ).getKrgColumn().back(); + }; + switch( getSaturationFunctionFamily( tm, ph ) ) { case SatfuncFamily::I: if( sgofTables.empty() && sgofLetTables.empty() && slgofTables.empty() ) @@ -722,6 +776,8 @@ namespace { return Opm::fun::map( famI_slgof, Opm::fun::iota( num_tables ) ); case SatfuncFamily::II: return Opm::fun::map( famII, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: + return Opm::fun::map( famIII, Opm::fun::iota( num_tables ) ); default: throw std::domain_error("No valid saturation keyword family specified"); } @@ -741,6 +797,7 @@ namespace { const auto& sgofLetTables = tm.getSgofletTable(); const auto& slgofTables = tm.getSlgofTables(); const auto& sgfnTables = tm.getSgfnTables(); + const auto& gsfTables = tm.getGsfTables(); auto sr = std::vector(num_tables, 0.0); if (ph.active(Opm::Phase::OIL)) { @@ -786,6 +843,14 @@ namespace { return sgfn.getKrgColumn().eval(ix); }; + const auto famIII = [&gsfTables, &sr](const int i) -> double + { + const auto& gsf = gsfTables.getTable(i); + const auto ix = gsf.getSgColumn().lookup(sr[i]); + + return gsf.getKrgColumn().eval(ix); + }; + switch( getSaturationFunctionFamily( tm, ph ) ) { case SatfuncFamily::I: if( sgofTables.empty() && sgofLetTables.empty() && slgofTables.empty() ) @@ -798,6 +863,8 @@ namespace { return Opm::fun::map( famI_slgof, Opm::fun::iota( num_tables ) ); case SatfuncFamily::II: return Opm::fun::map( famII, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: + return Opm::fun::map( famIII, Opm::fun::iota( num_tables ) ); default: throw std::domain_error("No valid saturation keyword family specified"); } @@ -816,6 +883,7 @@ namespace { const auto& swofTables = tm.getSwofTables(); const auto& swofLetTables = tm.getSwofletTable(); const auto& swfnTables = tm.getSwfnTables(); + const auto& wsfTables = tm.getWsfTables(); auto sr = std::vector(num_tables, 0.0); if (ph.active(Opm::Phase::OIL)) { @@ -853,6 +921,14 @@ namespace { return swfn.getKrwColumn().eval(ix); }; + const auto& famIII = [&wsfTables, &sr](const int i) -> double + { + const auto& wsf = wsfTables.getTable(i); + const auto ix = wsf.getSwColumn().lookup(sr[i]); + + return wsf.getKrwColumn().eval(ix); + }; + switch( getSaturationFunctionFamily( tm, ph ) ) { case SatfuncFamily::I: if( !swofTables.empty() ) @@ -864,6 +940,8 @@ namespace { case SatfuncFamily::II: return Opm::fun::map( famII, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: + return Opm::fun::map( famIII, Opm::fun::iota( num_tables ) ); default: throw std::domain_error("No valid saturation keyword family specified"); } @@ -930,6 +1008,8 @@ namespace { return ph.active(::Opm::Phase::GAS) ? Opm::fun::map( famII_3p, Opm::fun::iota( num_tables ) ) : Opm::fun::map( famII_2p, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: + throw std::domain_error("Saturation keyword family III is not applicable for a oil system"); default: throw std::domain_error("No valid saturation keyword family specified"); } @@ -1006,6 +1086,8 @@ namespace { return ph.active(::Opm::Phase::WATER) ? Opm::fun::map( famII_3p, Opm::fun::iota( num_tables ) ) : Opm::fun::map( famII_2p, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: + throw std::domain_error("Saturation keyword family III is not applicable for a oil system"); default: throw std::domain_error("No valid saturation keyword family specified"); } @@ -1063,6 +1145,8 @@ namespace { return Opm::fun::map( famI_slgof, Opm::fun::iota( num_tables ) ); case SatfuncFamily::II: return Opm::fun::map( famII, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: + throw std::domain_error("Saturation keyword family III is not applicable for a oil system"); default: throw std::domain_error("No valid saturation keyword family specified"); } @@ -1105,6 +1189,8 @@ namespace { case SatfuncFamily::II: return Opm::fun::map( famII, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: + throw std::domain_error("Saturation keyword family III is not applicable for a oil system"); default: throw std::domain_error("No valid saturation keyword family specified"); } @@ -1171,6 +1257,8 @@ namespace { return ph.active(::Opm::Phase::GAS) && ph.active(::Opm::Phase::WATER) ? Opm::fun::map( famII_3p, Opm::fun::iota( num_tables ) ) : Opm::fun::map( famII_2p, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: + throw std::domain_error("Saturation keyword family III is not applicable for a oil system"); default: throw std::domain_error("No valid saturation keyword family specified"); } @@ -1188,6 +1276,7 @@ namespace { const auto& swofTables = tm.getSwofTables(); const auto& swofLetTables = tm.getSwofletTable(); const auto& swfnTables = tm.getSwfnTables(); + const auto& wsfTables = tm.getWsfTables(); const auto& famI = [&swofTables]( int i ) { return swofTables.getTable( i ).getKrwColumn().back(); @@ -1201,6 +1290,10 @@ namespace { return swfnTables.getTable( i ).getKrwColumn().back(); }; + const auto& famIII = [&wsfTables]( int i ) { + return wsfTables.getTable( i ).getKrwColumn().back(); + }; + switch( getSaturationFunctionFamily( tm, ph ) ) { case SatfuncFamily::I: if( !swofTables.empty() ) @@ -1213,6 +1306,8 @@ namespace { case SatfuncFamily::II: return Opm::fun::map( famII, Opm::fun::iota( num_tables ) ); + case SatfuncFamily::III: + return Opm::fun::map( famIII, Opm::fun::iota( num_tables ) ); default: throw std::domain_error("No valid saturation keyword family specified"); } diff --git a/src/opm/input/eclipse/EclipseState/Runspec.cpp b/src/opm/input/eclipse/EclipseState/Runspec.cpp index 6f241c2342f..80de0f07e97 100644 --- a/src/opm/input/eclipse/EclipseState/Runspec.cpp +++ b/src/opm/input/eclipse/EclipseState/Runspec.cpp @@ -103,12 +103,19 @@ namespace { (twoP && deck.hasKeyword()))) || (wat && deck.hasKeyword()); + const auto family3 = //WSF, GSF gas-water CO2STORE case + deck.hasKeyword() && + deck.hasKeyword(); + if (family1) return Opm::SatFuncControls::KeywordFamily::Family_I; if (family2) return Opm::SatFuncControls::KeywordFamily::Family_II; + if (family3) + return Opm::SatFuncControls::KeywordFamily::Family_III; + return Opm::SatFuncControls::KeywordFamily::Undefined; } diff --git a/src/opm/input/eclipse/EclipseState/Tables/TableManager.cpp b/src/opm/input/eclipse/EclipseState/Tables/TableManager.cpp index ea1b9cd9795..bfa5af013a9 100644 --- a/src/opm/input/eclipse/EclipseState/Tables/TableManager.cpp +++ b/src/opm/input/eclipse/EclipseState/Tables/TableManager.cpp @@ -101,6 +101,8 @@ #include #include #include +#include +#include #include #include #include @@ -421,6 +423,9 @@ std::optional make_jfunc(const Deck& deck) { addTables( "SSFN", m_tabdims.getNumSatTables() ); addTables( "MSFN", m_tabdims.getNumSatTables() ); + addTables( "GSF", m_tabdims.getNumSatTables() ); + addTables( "WSF", m_tabdims.getNumSatTables() ); + addTables( "PLYADS", m_tabdims.getNumSatTables() ); addTables( "PLYROCK", m_tabdims.getNumSatTables()); addTables( "PLYVISC", m_tabdims.getNumPVTTables()); @@ -508,6 +513,9 @@ std::optional make_jfunc(const Deck& deck) { initSimpleTableContainer(deck, "SSFN" , m_tabdims.getNumSatTables()); initSimpleTableContainer(deck, "MSFN" , m_tabdims.getNumSatTables()); + initSimpleTableContainer(deck, "WSF" , m_tabdims.getNumSatTables()); + initSimpleTableContainer(deck, "GSF" , m_tabdims.getNumSatTables()); + initSimpleTableContainer(deck, "RSVD" , m_eqldims.getNumEquilRegions()); initSimpleTableContainer(deck, "RVVD" , m_eqldims.getNumEquilRegions()); initSimpleTableContainer(deck, "RVWVD" , m_eqldims.getNumEquilRegions()); @@ -917,6 +925,14 @@ std::optional make_jfunc(const Deck& deck) { return getTables("SSFN"); } + const TableContainer& TableManager::getWsfTables() const { + return getTables("WSF"); + } + + const TableContainer& TableManager::getGsfTables() const { + return getTables("GSF"); + } + const TableContainer& TableManager::getRsvdTables() const { return getTables("RSVD"); } diff --git a/src/opm/input/eclipse/EclipseState/Tables/Tables.cpp b/src/opm/input/eclipse/EclipseState/Tables/Tables.cpp index ffca0553f4f..7c26b13de33 100644 --- a/src/opm/input/eclipse/EclipseState/Tables/Tables.cpp +++ b/src/opm/input/eclipse/EclipseState/Tables/Tables.cpp @@ -90,6 +90,8 @@ #include #include #include +#include +#include #include @@ -730,6 +732,55 @@ SgfnTable::getJFuncColumn() const return SimpleTable::getColumn(2); } +GsfTable::GsfTable(const DeckItem& item, const int tableID) +{ + m_schema.addColumn(ColumnSchema("SG", Table::STRICTLY_INCREASING, Table::DEFAULT_NONE)); + m_schema.addColumn(ColumnSchema("KRG", Table::INCREASING, Table::DEFAULT_LINEAR)); + m_schema.addColumn(ColumnSchema("PCGW", Table::INCREASING, Table::DEFAULT_LINEAR)); + + SimpleTable::init("GSF", item, tableID); +} + + +const TableColumn& +GsfTable::getSgColumn() const +{ + return SimpleTable::getColumn(0); +} + +const TableColumn& +GsfTable::getKrgColumn() const +{ + return SimpleTable::getColumn(1); +} + +const TableColumn& +GsfTable::getPcgwColumn() const +{ + return SimpleTable::getColumn(2); +} + +WsfTable::WsfTable(const DeckItem& item, const int tableID) +{ + m_schema.addColumn(ColumnSchema("SW", Table::STRICTLY_INCREASING, Table::DEFAULT_NONE)); + m_schema.addColumn(ColumnSchema("KRW", Table::INCREASING, Table::DEFAULT_LINEAR)); + + SimpleTable::init("WSF", item, tableID); +} + + +const TableColumn& +WsfTable::getSwColumn() const +{ + return SimpleTable::getColumn(0); +} + +const TableColumn& +WsfTable::getKrwColumn() const +{ + return SimpleTable::getColumn(1); +} + SsfnTable::SsfnTable(const DeckItem& item, const int tableID) { m_schema.addColumn(ColumnSchema("SolventFraction", Table::STRICTLY_INCREASING, Table::DEFAULT_NONE)); diff --git a/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SGFN b/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SGFN index 224b06394b0..a83b21a1a7c 100644 --- a/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SGFN +++ b/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SGFN @@ -8,6 +8,7 @@ "item": "NTSFUN" }, "requires": ["GAS"], + "prohibits": ["GSF"], "items": [ { "name": "DATA", diff --git a/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SGOF b/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SGOF index 3868aff9cbe..7fb3c97ce9a 100644 --- a/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SGOF +++ b/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SGOF @@ -8,7 +8,7 @@ "item": "NTSFUN" }, "requires": ["GAS", "OIL"], - "prohibits": ["SLGOF"], + "prohibits": ["SLGOF", "GSF"], "items": [ { "name": "DATA", diff --git a/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SLGOF b/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SLGOF index e8312b356e6..d95c55f66ed 100644 --- a/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SLGOF +++ b/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SLGOF @@ -8,7 +8,7 @@ "item": "NTSFUN" }, "requires": ["GAS"], - "prohibits": ["SGOF"], + "prohibits": ["SGOF", "GSF"], "items": [ { "name": "DATA", diff --git a/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SWFN b/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SWFN index 3497208c34e..051cf0d6c7f 100644 --- a/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SWFN +++ b/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SWFN @@ -8,6 +8,7 @@ "item": "NTSFUN" }, "requires": ["WATER"], + "prohibits": ["WSF"], "items": [ { "name": "DATA", diff --git a/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SWOF b/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SWOF index aac5c41b02f..b30ae0000ce 100644 --- a/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SWOF +++ b/src/opm/input/eclipse/share/keywords/000_Eclipse100/S/SWOF @@ -8,6 +8,7 @@ "item": "NTSFUN" }, "requires": ["OIL", "WATER"], + "prohibits": ["WSF"], "items": [ { "name": "DATA", diff --git a/src/opm/input/eclipse/share/keywords/001_Eclipse300/G/GSF b/src/opm/input/eclipse/share/keywords/001_Eclipse300/G/GSF new file mode 100644 index 00000000000..21a951a88c4 --- /dev/null +++ b/src/opm/input/eclipse/share/keywords/001_Eclipse300/G/GSF @@ -0,0 +1,24 @@ +{ + "name": "GSF", + "sections": [ + "PROPS" + ], + "size": { + "keyword": "TABDIMS", + "item": "NTSFUN" + }, + "requires": ["GAS"], + "prohibits": ["SGOF", "SLGOF", "SGFN"], + "items": [ + { + "name": "DATA", + "value_type": "DOUBLE", + "size_type": "ALL", + "dimension": [ + "1", + "1", + "Pressure" + ] + } + ] +} diff --git a/src/opm/input/eclipse/share/keywords/001_Eclipse300/W/WSF b/src/opm/input/eclipse/share/keywords/001_Eclipse300/W/WSF new file mode 100644 index 00000000000..2df42de067a --- /dev/null +++ b/src/opm/input/eclipse/share/keywords/001_Eclipse300/W/WSF @@ -0,0 +1,23 @@ +{ + "name": "WSF", + "sections": [ + "PROPS" + ], + "size": { + "keyword": "TABDIMS", + "item": "NTSFUN" + }, + "requires": ["WATER"], + "prohibits": ["SWOF", "SWFN"], + "items": [ + { + "name": "DATA", + "value_type": "DOUBLE", + "size_type": "ALL", + "dimension": [ + "1", + "1" + ] + } + ] +} diff --git a/src/opm/input/eclipse/share/keywords/keyword_list.cmake b/src/opm/input/eclipse/share/keywords/keyword_list.cmake index 15073c22569..2d7e93a1d42 100644 --- a/src/opm/input/eclipse/share/keywords/keyword_list.cmake +++ b/src/opm/input/eclipse/share/keywords/keyword_list.cmake @@ -1056,6 +1056,7 @@ set( keywords 001_Eclipse300/D/DZV 001_Eclipse300/G/GASVISCT 001_Eclipse300/G/GCONPROD + 001_Eclipse300/G/GSF 001_Eclipse300/H/HEATCR 001_Eclipse300/H/HEATCRT 001_Eclipse300/L/LIVEOIL @@ -1082,6 +1083,7 @@ set( keywords 001_Eclipse300/T/TREFS 001_Eclipse300/W/WINJTEMP 001_Eclipse300/W/WATDENT + 001_Eclipse300/W/WSF 001_Eclipse300/Z/ZFACT1 001_Eclipse300/Z/ZFACT1S 001_Eclipse300/Z/ZFACTOR diff --git a/src/opm/material/fluidmatrixinteractions/EclMaterialLawManagerReadEffectiveParams.cpp b/src/opm/material/fluidmatrixinteractions/EclMaterialLawManagerReadEffectiveParams.cpp index 32a3f769f24..a9d1f2c249f 100644 --- a/src/opm/material/fluidmatrixinteractions/EclMaterialLawManagerReadEffectiveParams.cpp +++ b/src/opm/material/fluidmatrixinteractions/EclMaterialLawManagerReadEffectiveParams.cpp @@ -26,6 +26,9 @@ #include #include #include +#include +#include + #include namespace Opm { @@ -164,6 +167,11 @@ readGasOilParameters_(GasOilEffectiveParamVector& dest, unsigned satRegionIdx) break; } + case SatFuncControls::KeywordFamily::Family_III: + { + throw std::domain_error("Saturation keyword family III is not applicable for a gas-oil system"); + } + case SatFuncControls::KeywordFamily::Undefined: throw std::domain_error("No valid saturation keyword family specified"); } @@ -290,6 +298,29 @@ readGasWaterParameters_(GasWaterEffectiveParamVector& dest, unsigned satRegionId break; } + case SatFuncControls::KeywordFamily::Family_III: + { + const GsfTable& gsfTable = tableManager.getGsfTables().getTable( satRegionIdx ); + const WsfTable& wsfTable = tableManager.getWsfTables().getTable( satRegionIdx ); + + effParams.setApproach(SatCurveMultiplexerApproach::PiecewiseLinear); + auto& realParams = effParams.template getRealParams(); + + std::vector SwColumn = wsfTable.getColumn("SW").vectorCopy(); + + realParams.setKrwSamples(SwColumn, normalizeKrValues_(tolcrit, wsfTable.getColumn("KRW"))); + std::vector SwSamples(gsfTable.numRows()); + for (size_t sampleIdx = 0; sampleIdx < gsfTable.numRows(); ++ sampleIdx) + SwSamples[sampleIdx] = 1 - gsfTable.get("SG", sampleIdx); + realParams.setKrnSamples(SwSamples, normalizeKrValues_(tolcrit, gsfTable.getColumn("KRG"))); + //Capillary pressure is read from GSF. + // TODO need to check if gas/water sg + std::vector SgColumn = gsfTable.getColumn("SG").vectorCopy(); + realParams.setPcnwSamples(SgColumn, gsfTable.getColumn("PCGW").vectorCopy()); + realParams.finalize(); + + break; + } case SatFuncControls::KeywordFamily::Undefined: throw std::domain_error("No valid saturation keyword family specified"); } @@ -401,6 +432,11 @@ readOilWaterParameters_(OilWaterEffectiveParamVector& dest, unsigned satRegionId break; } + case SatFuncControls::KeywordFamily::Family_III: + { + throw std::domain_error("Saturation keyword family III is not applicable for a oil-water system"); + } + case SatFuncControls::KeywordFamily::Undefined: throw std::domain_error("No valid saturation keyword family specified"); } diff --git a/src/opm/output/eclipse/Tables.cpp b/src/opm/output/eclipse/Tables.cpp index b74b2d2025c..64f10eb57a1 100644 --- a/src/opm/output/eclipse/Tables.cpp +++ b/src/opm/output/eclipse/Tables.cpp @@ -34,6 +34,8 @@ #include #include #include +#include +#include #include #include #include @@ -420,6 +422,84 @@ namespace { namespace SatFunc { return numActRows; }); } + + + /// Create linearised and padded 'TAB' vector entries of normalised + /// GSF tables for all saturation function regions from Family Two + /// table data (GSF 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] gsf Collection of GSF tables for all saturation + /// regions. + /// + /// \return Linearised and padded 'TAB' vector values for output + /// GSF tables. A unit-converted copy of the input table \p + /// gsf with added derivatives. + std::vector + fromGSF(const std::size_t numRows, + const double tolcrit, + const Opm::UnitSystem& units, + const Opm::TableContainer& gsf) + { + using Gsf = ::Opm::GsfTable; + + const auto numTab = gsf.size(); + const auto numDep = std::size_t{2}; // Krg, Pcgw + + return detail::createSatfuncTable(numTab, numRows, numDep, + [tolcrit, &units, &gsf](const std::size_t tableID, + const std::size_t primID, + Opm::LinearisedOutputTable& linTable) + -> std::size_t + { + const auto& t = gsf.getTable(tableID); + + auto numActRows = std::size_t{0}; + + // Sg + { + const auto& Sg = t.getSgColumn(); + + numActRows = Sg.size(); + std::copy(std::begin(Sg), std::end(Sg), + linTable.column(tableID, primID, 0)); + } + + // Krg(Sg) + detail::outputRelperm(t.getKrgColumn(), tolcrit, + linTable.column(tableID, primID, 1)); + + // Pcgw(Sg) + { + const auto uPress = ::Opm::UnitSystem::measure::pressure; + + const auto& pc = t.getPcgwColumn(); + std::transform(std::begin(pc), std::end(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; + }); + } } // Gas /// Functions to create linearised, padded, and normalised SOFN output @@ -1164,6 +1244,82 @@ namespace { namespace SatFunc { return numActRows; }); } + + /// Create linearised and padded 'TAB' vector entries of normalised + /// WSF tables for all saturation function regions from Family Three + /// table data (WSF keyword). + /// + /// \param[in] numRows Number of rows to allocate for each table in + /// the output vector. 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] swfn Collection of WSF tables for all saturation + /// regions. + /// + /// \return Linearised and padded 'TAB' vector values for output + /// WSF tables. A unit-converted copy of the input table \p + /// wsf with added derivatives. + std::vector + fromWsf(const std::size_t numRows, + const double tolcrit, + const Opm::TableContainer& wsf) + { + using WSF = ::Opm::WsfTable; + + const auto numTab = wsf.size(); + const auto numDep = std::size_t{2}; // Krw, {zero pc} + + return detail::createSatfuncTable(numTab, numRows, numDep, + [tolcrit, &wsf] + (const std::size_t tableID, + const std::size_t primID, + Opm::LinearisedOutputTable& linTable) + -> std::size_t + { + const auto& t = wsf.getTable(tableID); + + auto numActRows = std::size_t{0}; + + // Sw + { + const auto& Sw = t.getSwColumn(); + + numActRows = Sw.size(); + std::copy(std::begin(Sw), std::end(Sw), + linTable.column(tableID, primID, 0)); + } + + // Krw(Sw) + detail::outputRelperm(t.getKrwColumn(), tolcrit, + linTable.column(tableID, primID, 1)); + + // Pcow = zero + { + const auto& krw = t.getKrwColumn(); + std::transform(std::begin(krw), std::end(krw), + linTable.column(tableID, primID, 2), + [](const double) -> double + { + return 0.0; + }); + } + + // Inform createSatfuncTable() of number of active rows in + // this table. Needed to compute slopes of piecewise linear + // interpolants. + return numActRows; + }); + } + } // Water @@ -2306,17 +2462,26 @@ namespace Opm { tabMgr.hasTables("SOF2"))) || (wat && tabMgr.hasTables("SWFN")); - if ((famI + famII) != 1) { - // Both Fam I && Fam II or neither of them. Can't have that. + const auto famIII = // GSF WSF + tabMgr.hasTables("GSF") && + tabMgr.hasTables("WSF"); + + if ((famI + famII + famIII) != 1) { + // No known saturation function family or at least two of the families I, II, or II. Can't have that. return; // Logging here? } else if (famI) { this->addSatFunc_FamilyOne(es, gas, oil, wat); } - else { + else if (famII){ // Family II this->addSatFunc_FamilyTwo(es, gas, oil, wat); } + else { + // Family II + assert(gas && wat && !oil); + this->addSatFunc_FamilyThree(es); + } } const std::vector& Tables::tabdims() const @@ -2502,6 +2667,38 @@ namespace Opm { } } + void Tables::addSatFunc_FamilyThree(const EclipseState& es) + { + const auto& tabMgr = es.getTableManager(); + const auto& tabd = es.runspec().tabdims(); + const auto nssfun = tabd.getNumSatNodes(); + const auto tolcrit = es.runspec() + .saturationFunctionControls() + .minimumRelpermMobilityThreshold(); + + { + const auto& tables = tabMgr.getGsfTables(); + const auto gsf = + SatFunc::Gas::fromGSF(nssfun, tolcrit, + this->units, tables); + + this->addData(Ix::SgfnTableStart, gsf); + this->m_tabdims[Ix::SgfnNumSatNodes] = nssfun; + this->m_tabdims[Ix::SgfnNumTables] = tables.size(); + } + + { + const auto& tables = tabMgr.getWsfTables(); + + const auto wsf = + SatFunc::Water::fromWsf(nssfun, tolcrit, tables); + + this->addData(Ix::SwfnTableStart, wsf); + this->m_tabdims[Ix::SwfnNumSatNodes] = nssfun; + this->m_tabdims[Ix::SwfnNumTables] = tables.size(); + } + } + void Tables::addGasPVTTables(const EclipseState& es) { const auto& tabMgr = es.getTableManager(); diff --git a/tests/material/test_eclmateriallawmanager.cpp b/tests/material/test_eclmateriallawmanager.cpp index 00cd05a28af..83b4921f2eb 100644 --- a/tests/material/test_eclmateriallawmanager.cpp +++ b/tests/material/test_eclmateriallawmanager.cpp @@ -452,6 +452,119 @@ static const char* letDeckString = " 0.0 0.03 1.8 1.9 1.0 0.95 0.0 0.01 3.5 4.0 1.1 1.0 1.0 1.0 1.0 0.2 0.01 /\n" "\n"; +static const char* fam3DeckStringGasWater = + "RUNSPEC\n" + "\n" + "DIMENS\n" + " 10 10 3 /\n" + "\n" + "TABDIMS\n" + "/\n" + "\n" + "WATER\n" + "GAS\n" + "\n" + "DISGAS\n" + "\n" + "FIELD\n" + "\n" + "GRID\n" + "\n" + "DX\n" + " 300*1000 /\n" + "DY\n" + " 300*1000 /\n" + "DZ\n" + " 100*20 100*30 100*50 /\n" + "\n" + "TOPS\n" + " 100*8325 /\n" + "\n" + "\n" + "PORO\n" + " 300*0.15 /\n" + "PROPS\n" + "\n" + "\n" + "GSF\n" + "0 0 0\n" + "0.001 0 0\n" + "0.02 0 0\n" + "0.05 0.005 0\n" + "0.12 0.025 0\n" + "0.2 0.075 0\n" + "0.25 0.125 0\n" + "0.3 0.190 0\n" + "0.4 0.410 0\n" + "0.45 0.60 0\n" + "0.5 0.72 0\n" + "0.6 0.87 0\n" + "0.7 0.94 0\n" + "0.85 0.98 0\n" + "0.88 0.984 0 /\n" + "\n" + "WSF\n" + "0.12 0 \n" + "0.22 0 \n" + "0.55 0.005 \n" + "0.88 0.984 /\n"; + +static const char* fam2DeckStringGasWater = + "RUNSPEC\n" + "\n" + "DIMENS\n" + " 10 10 3 /\n" + "\n" + "TABDIMS\n" + "/\n" + "\n" + "WATER\n" + "GAS\n" + "\n" + "DISGAS\n" + "\n" + "FIELD\n" + "\n" + "GRID\n" + "\n" + "DX\n" + " 300*1000 /\n" + "DY\n" + " 300*1000 /\n" + "DZ\n" + " 100*20 100*30 100*50 /\n" + "\n" + "TOPS\n" + " 100*8325 /\n" + "\n" + "\n" + "PORO\n" + " 300*0.15 /\n" + "PROPS\n" + "\n" + "\n" + "SGFN\n" + "0 0 0\n" + "0.001 0 0\n" + "0.02 0 0\n" + "0.05 0.005 0\n" + "0.12 0.025 0\n" + "0.2 0.075 0\n" + "0.25 0.125 0\n" + "0.3 0.190 0\n" + "0.4 0.410 0\n" + "0.45 0.60 0\n" + "0.5 0.72 0\n" + "0.6 0.87 0\n" + "0.7 0.94 0\n" + "0.85 0.98 0\n" + "0.88 0.984 0 /\n" + "\n" + "SWFN\n" + "0.12 0 0\n" + "0.22 0 0\n" + "0.55 0.005 0\n" + "0.88 0.984 0 /\n"; template inline Scalar computeLetCurve(const Scalar S, const Scalar L, const Scalar E, const Scalar T) @@ -676,6 +789,60 @@ inline void testAll() } } + // Water gas + { + const auto fam2Deck = parser.parseString(fam2DeckStringGasWater); + const Opm::EclipseState fam2EclState(fam2Deck); + + Opm::EclMaterialLawManager fam2materialLawManager; + fam2materialLawManager.initFromState(fam2EclState); + fam2materialLawManager.initParamsForElements(fam2EclState, n); + + const auto fam3Deck = parser.parseString(fam3DeckStringGasWater); + const Opm::EclipseState fam3EclState(fam3Deck); + + MaterialLawManager fam3materialLawManager; + fam3materialLawManager.initFromState(fam3EclState); + fam3materialLawManager.initParamsForElements(fam3EclState, n); + + for (unsigned elemIdx = 0; elemIdx < n; ++ elemIdx) { + for (int i = 0; i < 100; ++ i) { + Scalar Sw = 0; + Scalar So = Scalar(i)/100; + Scalar Sg = 1 - Sw - So; + FluidState fs; + fs.setSaturation(waterPhaseIdx, Sw); + fs.setSaturation(oilPhaseIdx, So); + fs.setSaturation(gasPhaseIdx, Sg); + + Scalar pcFam2[numPhases] = { 0.0, 0.0 }; + Scalar pcFam3[numPhases] = { 0.0, 0.0 }; + MaterialLaw::capillaryPressures(pcFam2, + fam2materialLawManager.materialLawParams(elemIdx), + fs); + MaterialLaw::capillaryPressures(pcFam3, + fam3materialLawManager.materialLawParams(elemIdx), + fs); + + Scalar krFam2[numPhases] = { 0.0, 0.0 }; + Scalar krFam3[numPhases] = { 0.0, 0.0 }; + MaterialLaw::relativePermeabilities(krFam2, + fam2materialLawManager.materialLawParams(elemIdx), + fs); + MaterialLaw::relativePermeabilities(krFam3, + fam3materialLawManager.materialLawParams(elemIdx), + fs); + + for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) { + if (std::abs(pcFam2[phaseIdx] - pcFam3[phaseIdx]) > 1e-5) + throw std::logic_error("Discrepancy between capillary pressure of family 2 and family 3 keywords"); + if (std::abs(krFam2[phaseIdx] - krFam3[phaseIdx]) > 1e-1) + throw std::logic_error("Discrepancy between relative permeabilities of family 2 and family 3 keywords"); + } + } + } + } + // LET: Curves from input parameters interpreted via MaterialLawManager versus directly computed curves. { const auto letDeck = parser.parseString(letDeckString);