diff --git a/opm/input/eclipse/EclipseState/SimulationConfig/RockConfig.hpp b/opm/input/eclipse/EclipseState/SimulationConfig/RockConfig.hpp index 44565379e..51f4aeeeb 100644 --- a/opm/input/eclipse/EclipseState/SimulationConfig/RockConfig.hpp +++ b/opm/input/eclipse/EclipseState/SimulationConfig/RockConfig.hpp @@ -20,6 +20,7 @@ #ifndef OPM_ROCK_CONFIG_HPP #define OPM_ROCK_CONFIG_HPP +#include #include #include @@ -28,36 +29,36 @@ namespace Opm { class Deck; class FieldPropsManager; -class RockConfig { +class RockConfig +{ public: - -enum class Hysteresis { - REVERS = 1, - IRREVERS = 2, - HYSTER = 3, - BOBERG = 4, - REVLIMIT = 5, - PALM_MAN = 6, - NONE = 7 -}; - - -struct RockComp { - double pref; - double compressibility; - - RockComp() = default; - RockComp(double pref_arg, double comp_arg); - bool operator==(const RockComp& other) const; - - template - void serializeOp(Serializer& serializer) + enum class Hysteresis { - serializer(pref); - serializer(compressibility); - } -}; + REVERS = 1, + IRREVERS = 2, + HYSTER = 3, + BOBERG = 4, + REVLIMIT = 5, + PALM_MAN = 6, + NONE = 7, + }; + struct RockComp + { + double pref{}; + double compressibility{}; + + RockComp() = default; + RockComp(double pref_arg, double comp_arg); + bool operator==(const RockComp& other) const; + + template + void serializeOp(Serializer& serializer) + { + serializer(pref); + serializer(compressibility); + } + }; RockConfig(); RockConfig(const Deck& deck, const FieldPropsManager& fp); @@ -95,4 +96,4 @@ private: } //namespace Opm -#endif +#endif // OPM_ROCK_CONFIG_HPP diff --git a/opm/input/eclipse/EclipseState/Tables/FlatTable.hpp b/opm/input/eclipse/EclipseState/Tables/FlatTable.hpp index 8fab74733..804982521 100644 --- a/opm/input/eclipse/EclipseState/Tables/FlatTable.hpp +++ b/opm/input/eclipse/EclipseState/Tables/FlatTable.hpp @@ -313,13 +313,22 @@ struct ROCKRecord { } }; -struct RockTable : public FlatTable< ROCKRecord > { - using FlatTable< ROCKRecord >::FlatTable; +struct RockTable : public FlatTableWithCopy +{ + RockTable() = default; + explicit RockTable(const DeckKeyword& kw); + explicit RockTable(std::initializer_list records); static RockTable serializationTestObject() { return RockTable({{1.0, 2.0}}); } + + template + void serializeOp(Serializer& serializer) + { + FlatTableWithCopy::serializeOp(serializer); + } }; struct PVCDORecord { diff --git a/src/opm/input/eclipse/EclipseState/SimulationConfig/RockConfig.cpp b/src/opm/input/eclipse/EclipseState/SimulationConfig/RockConfig.cpp index d18ab389b..9a9e1b657 100644 --- a/src/opm/input/eclipse/EclipseState/SimulationConfig/RockConfig.cpp +++ b/src/opm/input/eclipse/EclipseState/SimulationConfig/RockConfig.cpp @@ -15,73 +15,119 @@ You should have received a copy of the GNU General Public License along with OPM. If not, see . - */ +*/ -#include +#include + +#include +#include + +#include #include -#include -#include -#include -namespace Opm { +#include +#include +#include namespace { - - -std::string num_prop(const std::string& prop_name) { +std::string num_prop(const std::string& prop_name) +{ static const std::unordered_set prop_names = {"PVTNUM", "SATNUM", "ROCKNUM"}; - if (prop_names.count(prop_name) == 1) + if (prop_names.count(prop_name) == 1) { return prop_name; + } throw std::invalid_argument("The rocknum propertype: " + prop_name + " is not valid"); } - -RockConfig::Hysteresis hysteresis(const std::string& s) { +Opm::RockConfig::Hysteresis hysteresis(const std::string& s) +{ if (s == "REVERS") - return RockConfig::Hysteresis::REVERS; + return Opm::RockConfig::Hysteresis::REVERS; if (s == "IRREVERS") - return RockConfig::Hysteresis::IRREVERS; + return Opm::RockConfig::Hysteresis::IRREVERS; if (s == "HYSTER") - return RockConfig::Hysteresis::HYSTER; + return Opm::RockConfig::Hysteresis::HYSTER; if (s == "BOBERG") - return RockConfig::Hysteresis::BOBERG; + return Opm::RockConfig::Hysteresis::BOBERG; if (s == "REVLIMIT") - return RockConfig::Hysteresis::REVLIMIT; + return Opm::RockConfig::Hysteresis::REVLIMIT; if (s == "PALM-MAN") - return RockConfig::Hysteresis::PALM_MAN; + return Opm::RockConfig::Hysteresis::PALM_MAN; if (s == "NONE") - return RockConfig::Hysteresis::NONE; + return Opm::RockConfig::Hysteresis::NONE; throw std::invalid_argument("Not recognized hysteresis option: " + s); } -} +} // Anonymous namespace -RockConfig::RockConfig() : - num_property(ParserKeywords::ROCKOPTS::TABLE_TYPE::defaultValue), - num_tables(ParserKeywords::ROCKCOMP::NTROCC::defaultValue) -{ -} +// =========================================================================== -bool RockConfig::RockComp::operator==(const RockConfig::RockComp& other) const { - return this->pref == other.pref && - this->compressibility == other.compressibility; -} +namespace Opm { -RockConfig::RockComp::RockComp(double pref_arg, double comp_arg) : - pref(pref_arg), - compressibility(comp_arg) +RockConfig::RockComp::RockComp(double pref_arg, double comp_arg) + : pref(pref_arg) + , compressibility(comp_arg) {} +bool RockConfig::RockComp::operator==(const RockConfig::RockComp& other) const +{ + return (this->pref == other.pref) + && (this->compressibility == other.compressibility); +} + +// --------------------------------------------------------------------------- + +RockConfig::RockConfig() + : num_property(ParserKeywords::ROCKOPTS::TABLE_TYPE::defaultValue) + , num_tables (ParserKeywords::ROCKCOMP::NTROCC::defaultValue) +{} + +RockConfig::RockConfig(const Deck& deck, const FieldPropsManager& fp) + : RockConfig {} +{ + using rock = ParserKeywords::ROCK; + using rockopts = ParserKeywords::ROCKOPTS; + using rockcomp = ParserKeywords::ROCKCOMP; + + if (deck.hasKeyword()) { + for (const auto& table : RockTable { deck.get().back() }) { + this->m_comp.emplace_back(table.reference_pressure, + table.compressibility); + } + } + + if (deck.hasKeyword()) { + const auto& record = deck.get().back().getRecord(0); + this->num_property = num_prop( record.getItem().getTrimmedString(0) ); + } + + if (deck.hasKeyword()) { + const auto& record = deck.get().back().getRecord(0); + if (fp.has_int("ROCKNUM")) { + this->num_property = "ROCKNUM"; + } + + this->num_tables = record.getItem().get(0); + this->hyst_mode = hysteresis(record.getItem().getTrimmedString(0)); + this->m_water_compaction = DeckItem::to_bool(record.getItem().getTrimmedString(0)); + + this->m_active = true; + if ((this->hyst_mode == Hysteresis::NONE) && !this->m_water_compaction) { + this->m_active = false; + } + } +} + RockConfig RockConfig::serializationTestObject() { RockConfig result; @@ -95,71 +141,44 @@ RockConfig RockConfig::serializationTestObject() return result; } -bool RockConfig::water_compaction() const { - return this->m_water_compaction; -} - -const std::vector& RockConfig::comp() const { - return this->m_comp; -} - -std::size_t RockConfig::num_rock_tables() const { - return this->num_tables; -} - -const std::string& RockConfig::rocknum_property() const { - return this->num_property; -} - -RockConfig::Hysteresis RockConfig::hysteresis_mode() const { - return this->hyst_mode; -} - -bool RockConfig::active() const { +bool RockConfig::active() const +{ return this->m_active; } -RockConfig::RockConfig(const Deck& deck, const FieldPropsManager& fp) - : RockConfig {} +const std::vector& RockConfig::comp() const { - using rock = ParserKeywords::ROCK; - using rockopts = ParserKeywords::ROCKOPTS; - using rockcomp = ParserKeywords::ROCKCOMP; - if (deck.hasKeyword()) { - const auto& rock_kw = deck.get().back(); - for (const auto& record : rock_kw) - this->m_comp.emplace_back( record.getItem().getSIDouble(0), - record.getItem().getSIDouble(0)); - } - - if (deck.hasKeyword()) { - const auto& record = deck.get().back().getRecord(0); - this->num_property = num_prop( record.getItem().getTrimmedString(0) ); - } - - if (deck.hasKeyword()) { - const auto& record = deck.get().back().getRecord(0); - if (fp.has_int("ROCKNUM")) - this->num_property = "ROCKNUM"; - - this->num_tables = record.getItem().get(0); - this->hyst_mode = hysteresis(record.getItem().getTrimmedString(0)); - this->m_water_compaction = DeckItem::to_bool(record.getItem().getTrimmedString(0)); - - this->m_active = true; - if (this->hyst_mode == Hysteresis::NONE && !this->m_water_compaction) - this->m_active = false; - } + return this->m_comp; } +const std::string& RockConfig::rocknum_property() const +{ + return this->num_property; +} -bool RockConfig::operator==(const RockConfig& other) const { - return this->num_property == other.num_property && - this->m_comp == other.m_comp && - this->num_tables == other.num_tables && - this->m_active == other.m_active && - this->m_water_compaction == other.m_water_compaction && - this->hyst_mode == other.hyst_mode; +std::size_t RockConfig::num_rock_tables() const +{ + return this->num_tables; +} + +RockConfig::Hysteresis RockConfig::hysteresis_mode() const +{ + return this->hyst_mode; +} + +bool RockConfig::water_compaction() const +{ + return this->m_water_compaction; +} + +bool RockConfig::operator==(const RockConfig& other) const +{ + return (this->num_property == other.num_property) + && (this->m_comp == other.m_comp) + && (this->num_tables == other.num_tables) + && (this->m_active == other.m_active) + && (this->m_water_compaction == other.m_water_compaction) + && (this->hyst_mode == other.hyst_mode); } } //namespace Opm diff --git a/src/opm/input/eclipse/EclipseState/Tables/Tables.cpp b/src/opm/input/eclipse/EclipseState/Tables/Tables.cpp index 0868503bd..9b58faf0f 100644 --- a/src/opm/input/eclipse/EclipseState/Tables/Tables.cpp +++ b/src/opm/input/eclipse/EclipseState/Tables/Tables.cpp @@ -27,6 +27,7 @@ #include #include #include +#include #include #include #include @@ -110,7 +111,6 @@ namespace Opm PvtgTable::PvtgTable(const DeckKeyword& keyword, size_t tableIdx) : PvtxTable("P") { - m_underSaturatedSchema.addColumn(ColumnSchema("RV", Table::STRICTLY_DECREASING, Table::DEFAULT_NONE)); m_underSaturatedSchema.addColumn(ColumnSchema("BG", Table::RANDOM, Table::DEFAULT_LINEAR)); m_underSaturatedSchema.addColumn(ColumnSchema("MUG", Table::RANDOM, Table::DEFAULT_LINEAR)); @@ -141,7 +141,6 @@ PvtgTable::operator==(const PvtgTable& data) const PvtgwTable::PvtgwTable(const DeckKeyword& keyword, size_t tableIdx) : PvtxTable("P") { - m_underSaturatedSchema.addColumn(ColumnSchema("RW", Table::STRICTLY_DECREASING, Table::DEFAULT_NONE)); m_underSaturatedSchema.addColumn(ColumnSchema("BG", Table::RANDOM, Table::DEFAULT_LINEAR)); m_underSaturatedSchema.addColumn(ColumnSchema("MUG", Table::RANDOM, Table::DEFAULT_LINEAR)); @@ -1968,7 +1967,6 @@ TracerVdTable::getTracerConcentration() const namespace { - /* * Create a compile-time sequence of integers [0,N). In C++14 this can be * replaced by std::index_sequence. @@ -1977,46 +1975,47 @@ namespace struct seq { using type = seq; }; + template struct mkseq : mkseq { }; + template struct mkseq<0u, Is...> : seq { using type = seq; }; - /* - * Convenince function for creating a 'flat table', e.g. PVTW and DENSITY. - * Assumes the following: - * - * 1. The table has vector semantics with no other to enforce - * 2. That the following struct is implemented: - * struct record { - * static constexpr std::size_t size = [number-of-members] - * double members ...; - * }; - * 3. The table is declared as - * struct table : public FlatTable< table > { - * using FlatTable< table >::FlatTable; - * } - * - * If some field can *not* be defaulted, e.g. 0, specialise the flat_props - * struct (in this namespace) as such: - * template<> struct< record, 0 > { - * static constexpr bool can_default() { return false; } - * static constexpr const char* errmsg() { "error message"; } - * }; - * and the parser will throw std::invalid_argument if the field is defaulted in - * the input. - * - */ + // Convenince function for creating a 'flat table', e.g. PVTW and + // DENSITY. Assumes the following: + // + // 1. The table has vector semantics with no other to enforce + // 2. That the following struct is implemented: + // struct record { + // static constexpr std::size_t size = [number-of-members] + // double members ...; + // }; + // 3. The table is declared as + // struct table : public FlatTable< table > { + // using FlatTable< table >::FlatTable; + // } + // + // If some field can *not* be defaulted, e.g. 0, specialise the + // flat_props struct (in this namespace) as such: + // template<> struct { + // static constexpr bool can_default() { return false; } + // static constexpr const char* errmsg() { "error message"; } + // }; + // and the parser will throw std::invalid_argument if the field is + // defaulted in the input. template - struct flat_props { + struct flat_props + { static constexpr bool can_default() { return true; } + static constexpr const char* errmsg() { return ""; @@ -2029,7 +2028,9 @@ namespace const auto& item = rec.getItem(N); if (item.defaultApplied(0) && !flat_props::can_default()) { - throw std::invalid_argument {flat_props::errmsg()}; + throw std::invalid_argument { + flat_props::errmsg() + }; } return item.getSIDouble(0); @@ -2038,7 +2039,7 @@ namespace template T flat_get(const DeckRecord& record, seq) { - return {flat_get(record)...}; + return { flat_get(record)... }; } template @@ -2055,11 +2056,13 @@ namespace } template <> - struct flat_props { + struct flat_props + { static constexpr bool can_default() { return false; } + static constexpr const char* errmsg() { return "PVTW reference pressure cannot be defaulted"; @@ -2075,11 +2078,13 @@ namespace }; template - struct flat_props { + struct flat_props + { static constexpr bool can_default() { return false; } + static constexpr const char* errmsg() { return pvcdo_err[N]; @@ -2099,13 +2104,15 @@ namespace // ------------------------------------------------------------------------ template -FlatTableWithCopy::FlatTableWithCopy(const DeckKeyword& kw, std::string_view expect) +FlatTableWithCopy::FlatTableWithCopy(const DeckKeyword& kw, + std::string_view expect) { if (!expect.empty() && (kw.name() != expect)) { - throw std::invalid_argument {fmt::format("Keyword {} cannot be used to " - "initialise {} table structures", - kw.name(), - expect)}; + throw std::invalid_argument { + fmt::format("Keyword {} cannot be used to " + "initialise {} table structures", + kw.name(), expect) + }; } this->table_.reserve(kw.size()); @@ -2116,53 +2123,55 @@ FlatTableWithCopy::FlatTableWithCopy(const DeckKeyword& kw, std::str // table in region R-1. Table must not be defaulted in region 1 // (i.e., when PVTNUM=1). if (this->table_.empty()) { - throw OpmInputError {"First record cannot be defaulted", kw.location()}; + throw OpmInputError { + "First record cannot be defaulted", + kw.location() + }; } this->table_.push_back(this->table_.back()); - } else { - this->table_.push_back(flat_get(record, mkseq {})); + } + else { + this->table_.push_back(flat_get(record, mkseq{})); } } } template -FlatTableWithCopy::FlatTableWithCopy(std::initializer_list records) - : table_ {records} -{ -} +FlatTableWithCopy:: +FlatTableWithCopy(std::initializer_list records) + : table_{ records } +{} // ------------------------------------------------------------------------ GravityTable::GravityTable(const DeckKeyword& kw) : FlatTableWithCopy(kw, ParserKeywords::GRAVITY::keywordName) -{ -} +{} GravityTable::GravityTable(std::initializer_list records) : FlatTableWithCopy(records) -{ -} +{} // ------------------------------------------------------------------------ DensityTable::DensityTable(const DeckKeyword& kw) : FlatTableWithCopy(kw, ParserKeywords::DENSITY::keywordName) -{ -} +{} DensityTable::DensityTable(std::initializer_list records) : FlatTableWithCopy(records) -{ -} +{} DensityTable::DensityTable(const GravityTable& gravity) { this->table_.reserve(gravity.size()); - constexpr auto default_air_density = 1.22 * unit::kilogram / unit::cubic(unit::meter); + constexpr auto default_air_density = + 1.22 * unit::kilogram / unit::cubic(unit::meter); - constexpr auto default_water_density = 1000.0 * unit::kilogram / unit::cubic(unit::meter); + constexpr auto default_water_density = + 1000.0 * unit::kilogram / unit::cubic(unit::meter); // Degrees API defined as // @@ -2170,10 +2179,15 @@ DensityTable::DensityTable(const GravityTable& gravity) // // with SG being the specific gravity of oil relative to pure water. - std::transform(gravity.begin(), gravity.end(), std::back_inserter(this->table_), [](const GRAVITYRecord& record) { - return DENSITYRecord {(141.5 / (record.oil_api + 131.5)) * default_water_density, - record.water_sg * default_water_density, - record.gas_sg * default_air_density}; + std::transform(gravity.begin(), gravity.end(), + std::back_inserter(this->table_), + [](const GRAVITYRecord& record) + { + return DENSITYRecord { + (141.5 / (record.oil_api + 131.5)) * default_water_density, + record.water_sg * default_water_density, + record.gas_sg * default_air_density + }; }); } @@ -2181,34 +2195,40 @@ DensityTable::DensityTable(const GravityTable& gravity) PvtwTable::PvtwTable(const DeckKeyword& kw) : FlatTableWithCopy(kw, ParserKeywords::PVTW::keywordName) -{ -} +{} PvtwTable::PvtwTable(std::initializer_list records) : FlatTableWithCopy(records) -{ -} +{} + +// ------------------------------------------------------------------------ + +RockTable::RockTable(const DeckKeyword& kw) + : FlatTableWithCopy(kw, ParserKeywords::ROCK::keywordName) +{} + +RockTable::RockTable(std::initializer_list records) + : FlatTableWithCopy(records) +{} // ------------------------------------------------------------------------ template FlatTable::FlatTable(const DeckKeyword& kw) - : std::vector(flat_records(kw, mkseq {})) -{ -} + : std::vector(flat_records(kw, mkseq{})) +{} template FlatTable::FlatTable(const DeckKeyword&); template FlatTable::FlatTable(const DeckKeyword&); template FlatTable::FlatTable(const DeckKeyword&); template FlatTable::FlatTable(const DeckKeyword&); -template FlatTable::FlatTable(const DeckKeyword&); -template FlatTable::FlatTable(const DeckKeyword&); -template FlatTable::FlatTable(const DeckKeyword&); template FlatTable::FlatTable(const DeckKeyword&); +template FlatTable::FlatTable(const DeckKeyword&); +template FlatTable::FlatTable(const DeckKeyword&); template FlatTable::FlatTable(const DeckKeyword&); template FlatTable::FlatTable(const DeckKeyword&); template FlatTable::FlatTable(const DeckKeyword&); +template FlatTable::FlatTable(const DeckKeyword&); template FlatTable::FlatTable(const DeckKeyword&); -template FlatTable::FlatTable(const DeckKeyword&); } // namespace Opm diff --git a/tests/parser/SimulationConfigTest.cpp b/tests/parser/SimulationConfigTest.cpp index 9a24ccc0b..30dc78e7c 100644 --- a/tests/parser/SimulationConfigTest.cpp +++ b/tests/parser/SimulationConfigTest.cpp @@ -17,19 +17,18 @@ along with OPM. If not, see . */ - - #define BOOST_TEST_MODULE SimulationConfigTests #include +#include + #include #include #include #include #include -#include #include #include @@ -38,6 +37,9 @@ #include #include +#include +#include + using namespace Opm; const std::string& inputStr = "RUNSPEC\n" @@ -280,7 +282,8 @@ BOOST_AUTO_TEST_CASE(SimulationConfig_TEMP_THERMAL) } -BOOST_AUTO_TEST_CASE(TESTRockConfig) { +BOOST_AUTO_TEST_CASE(TESTRockConfig_Standard) +{ const std::string deck_string = R"( RUNSPEC @@ -310,3 +313,45 @@ ROCKOPTS const auto& comp = rc.comp(); BOOST_CHECK_EQUAL(comp.size(), 3U); } + +BOOST_AUTO_TEST_CASE(TESTRockConfig_Default) +{ + const auto deck = Parser{}.parseString(R"( +RUNSPEC + +TABDIMS + 3 / -- NTSFUN = 3 + +PROPS + +ROCKOPTS + 1* 1* SATNUM / + +ROCK +123.4 0.40E-05 / +/ +271.8 1.61e-05 / + +)"); + + const auto grid = EclipseGrid { 10, 10, 10 }; + const auto fp = FieldPropsManager { + deck, Phases{true, true, true}, grid, TableManager() + }; + + const auto rc = RockConfig { deck, fp }; + + BOOST_CHECK_EQUAL(rc.rocknum_property(), "SATNUM"); + + const auto& comp = rc.comp(); + BOOST_REQUIRE_EQUAL(comp.size(), std::size_t{3}); + + BOOST_CHECK_CLOSE(comp[0].pref, 123.4*1.0e5, 1.0e-8); + BOOST_CHECK_CLOSE(comp[0].compressibility, 0.4e-5/1.0e5, 1.0e-8); + + BOOST_CHECK_CLOSE(comp[1].pref, 123.4*1.0e5, 1.0e-8); + BOOST_CHECK_CLOSE(comp[1].compressibility, 0.4e-5/1.0e5, 1.0e-8); + + BOOST_CHECK_CLOSE(comp[2].pref, 271.8*1.0e5, 1.0e-8); + BOOST_CHECK_CLOSE(comp[2].compressibility, 1.61e-5/1.0e5, 1.0e-8); +} diff --git a/tests/parser/TableManagerTests.cpp b/tests/parser/TableManagerTests.cpp index a70df923e..2971a6160 100644 --- a/tests/parser/TableManagerTests.cpp +++ b/tests/parser/TableManagerTests.cpp @@ -2398,6 +2398,28 @@ BOOST_AUTO_TEST_CASE( TestParseROCK ) { BOOST_CHECK_EQUAL( 8U , tables.numFIPRegions( )); } +BOOST_AUTO_TEST_CASE( TestParseROCK_WithDefault ) +{ + const auto deck = Opm::Parser{}.parseString(R"(RUNSPEC +TABDIMS + 1* 2 / +PROPS +ROCK + 1.1 1.2 / +/ -- Copy from region 1 +)"); + + const auto tables = Opm::TableManager { deck }; + const auto& rock = tables.getRockTable(); + BOOST_CHECK_EQUAL(rock.size(), std::size_t{2}); + + BOOST_CHECK_CLOSE(1.1e5, rock[0].reference_pressure, 1.0e-8); + BOOST_CHECK_CLOSE(1.2e-5, rock[0].compressibility, 1.0e-8); + + BOOST_CHECK_CLOSE(1.1e5, rock[1].reference_pressure, 1.0e-8); + BOOST_CHECK_CLOSE(1.2e-5, rock[1].compressibility, 1.0e-8); +} + BOOST_AUTO_TEST_CASE( TestParsePVCDO ) { const std::string data = R"( TABDIMS