Merge pull request #2303 from joakim-hove/serialize-rockconfig

Serialize rockconfig

ebos/eclproblem.hh

@@ -86,6 +86,7 @@
 #include <opm/material/common/Valgrind.hpp>
 #include <opm/parser/eclipse/Deck/Deck.hpp>
 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
+#include <opm/parser/eclipse/EclipseState/SimulationConfig/RockConfig.hpp>
 #include <opm/parser/eclipse/EclipseState/Tables/Eqldims.hpp>
 #include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
 #include <opm/parser/eclipse/EclipseState/Schedule/Action/ActionContext.hpp>
@@ -2162,71 +2163,34 @@ private:
     {
         const auto& simulator = this->simulator();
         const auto& vanguard = simulator.vanguard();
-        const auto& deck = vanguard.deck();
+        const auto& eclState = vanguard.eclState();
-        auto& eclState = vanguard.eclState();
+        const auto& rock_config = eclState.getSimulationConfig().rock_config();
 
         // read the rock compressibility parameters
-        if (deck.hasKeyword("ROCK")) {
+        {
-            const auto& rockKeyword = deck.getKeyword("ROCK");
+            const auto& comp = rock_config.comp();
-            rockParams_.resize(rockKeyword.size());
+            rockParams_.clear();
-            for (size_t rockRecordIdx = 0; rockRecordIdx < rockKeyword.size(); ++ rockRecordIdx) {
+            for (const auto& c : comp)
-                const auto& rockRecord = rockKeyword.getRecord(rockRecordIdx);
+                rockParams_.push_back( { c.pref, c.compressibility } );
-                rockParams_[rockRecordIdx].referencePressure =
-                        rockRecord.getItem("PREF").getSIDouble(0);
-                rockParams_[rockRecordIdx].compressibility =
-                        rockRecord.getItem("COMPRESSIBILITY").getSIDouble(0);
-            }
         }
 
         // read the parameters for water-induced rock compaction
         readRockCompactionParameters_();
 
-        // check the kind of region which is supposed to be used by checking the ROCKOPTS
+        const auto& num = eclState.fieldProps().get_global_int(rock_config.rocknum_property());
-        // keyword. note that for some funny reason, the ROCK keyword uses PVTNUM by
+        unsigned numElem = vanguard.gridView().size(0);
-        // default, *not* ROCKNUM!
+        rockTableIdx_.resize(numElem);
-        std::string propName = "PVTNUM";
+        for (size_t elemIdx = 0; elemIdx < numElem; ++ elemIdx) {
-        if (deck.hasKeyword("ROCKOPTS")) {
+            unsigned cartElemIdx = vanguard.cartesianIndex(elemIdx);
-            const auto& rockoptsKeyword = deck.getKeyword("ROCKOPTS");
-            std::string rockTableType =
-                rockoptsKeyword.getRecord(0).getItem("TABLE_TYPE").getTrimmedString(0);
-            if (rockTableType == "PVTNUM")
-                propName = "PVTNUM";
-            else if (rockTableType == "SATNUM")
-                propName = "SATNUM";
-            else if (rockTableType == "ROCKNUM")
-                propName = "ROCKNUM";
-            else {
-                throw std::runtime_error("Unknown table type '"+rockTableType
-                                         +" for the ROCKOPTS keyword given");
-            }
-        }
 
-        // If ROCKCOMP is used and ROCKNUM is specified ROCK2D ROCK2DTR ROCKTAB etc. uses ROCKNUM
+            rockTableIdx_[elemIdx] = num[cartElemIdx] - 1;
-        // to give the correct table index.
-        if (deck.hasKeyword("ROCKCOMP") && eclState.fieldProps().has_int("ROCKNUM"))
-            propName = "ROCKNUM";
-
-        if (eclState.fieldProps().has_int(propName)) {
-            const auto& tmp = eclState.fieldProps().get_global_int(propName);
-            unsigned numElem = vanguard.gridView().size(0);
-            rockTableIdx_.resize(numElem);
-            for (size_t elemIdx = 0; elemIdx < numElem; ++ elemIdx) {
-                unsigned cartElemIdx = vanguard.cartesianIndex(elemIdx);
-
-                // reminder: Eclipse uses FORTRAN-style indices
-                rockTableIdx_[elemIdx] = tmp[cartElemIdx] - 1;
-            }
         }
 
         // Store overburden pressure pr element
         const auto& overburdTables = eclState.getTableManager().getOverburdTables();
         if (!overburdTables.empty()) {
-            unsigned numElem = vanguard.gridView().size(0);
             overburdenPressure_.resize(numElem,0.0);
-
+            size_t numRocktabTables = rock_config.num_rock_tables();
-            const auto& rockcomp = deck.getKeyword("ROCKCOMP");
-            const auto& rockcompRecord = rockcomp.getRecord(0);
-            size_t numRocktabTables = rockcompRecord.getItem("NTROCC").template get< int >(0);
 
             if (overburdTables.size() != numRocktabTables)
                 throw std::runtime_error(std::to_string(numRocktabTables) +" OVERBURD tables is expected, but " + std::to_string(overburdTables.size()) +" is provided");
@@ -2251,47 +2215,36 @@ private:
     void readRockCompactionParameters_()
     {
         const auto& vanguard = this->simulator().vanguard();
-        const auto& deck = vanguard.deck();
         const auto& eclState = vanguard.eclState();
+        const auto& rock_config = eclState.getSimulationConfig().rock_config();
 
-        if (!deck.hasKeyword("ROCKCOMP"))
+        if (!rock_config.active())
             return; // deck does not enable rock compaction
 
-        const auto& rockcomp = deck.getKeyword("ROCKCOMP");
+        unsigned numElem = vanguard.gridView().size(0);
-        //for (size_t rockRecordIdx = 0; rockRecordIdx < rockcomp.size(); ++ rockRecordIdx) {
+        switch (rock_config.hysteresis_mode()) {
-        assert(rockcomp.size() == 1);
+        case RockConfig::Hysteresis::REVERS:
-        const auto& rockcompRecord = rockcomp.getRecord(0);
+            break;
-        const auto& option = rockcompRecord.getItem("HYSTERESIS").getTrimmedString(0);
+        case RockConfig::Hysteresis::IRREVERS:
-        if (option == "REVERS") {
-            // interpolate the porv volume multiplier using the pressure in the cell
-        }
-        else if (option == "IRREVERS") {
             // interpolate the porv volume multiplier using the minimum pressure in the cell
             // i.e. don't allow re-inflation.
-            unsigned numElem = vanguard.gridView().size(0);
             minOilPressure_.resize(numElem, 1e99);
+            break;
+        default:
+            throw std::runtime_error("Not support ROCKOMP hysteresis option ");
         }
-        else if (option == "NO")
-            // rock compaction turned on but disabled by ROCKCOMP option
-            return;
-        else
-            throw std::runtime_error("ROCKCOMP option " + option + " not supported for item 1");
 
-        size_t numRocktabTables = rockcompRecord.getItem("NTROCC").template get<int>(0);
+        size_t numRocktabTables = rock_config.num_rock_tables();
-        const auto& waterCompactionItem = rockcompRecord.getItem("WATER_COMPACTION").getTrimmedString(0);
+        bool waterCompaction = rock_config.water_compaction();
-        bool waterCompaction = false;
+
-        if (waterCompactionItem == "YES") {
+        if (!waterCompaction)
-            waterCompaction = true;
+            throw std::runtime_error("Only water compatction allowed");
-            unsigned numElem = vanguard.gridView().size(0);
-            maxWaterSaturation_.resize(numElem, 0.0);
-        }
-        else
-            throw std::runtime_error("ROCKCOMP option " + waterCompactionItem + " not supported for item 3. Only YES is supported");
 
         if (waterCompaction) {
             const auto& rock2dTables = eclState.getTableManager().getRock2dTables();
             const auto& rock2dtrTables = eclState.getTableManager().getRock2dtrTables();
             const auto& rockwnodTables = eclState.getTableManager().getRockwnodTables();
+            maxWaterSaturation_.resize(numElem, 0.0);
 
             if (rock2dTables.size() != numRocktabTables)
                 throw std::runtime_error("Water compation option is selected in ROCKCOMP." + std::to_string(numRocktabTables)
@@ -3022,7 +2975,7 @@ private:
             for (const auto& bcface : bcconfig) {
                 const auto& type = bcface.bctype;
                 if (type == BCType::RATE) {
-                    int compIdx = 0;
+                    int compIdx = 0; // default initialize to avoid -Wmaybe-uninitialized warning
 
                     switch (bcface.component) {
                     case BCComponent::OIL:
@@ -3047,8 +3000,7 @@ private:
                         compIdx = Indices::polymerConcentrationIdx;
                         break;
                     case BCComponent::NONE:
-                        if (type == BCType::RATE)
+                        throw std::logic_error("you need to specify the component when RATE type is set in BC");
-                            throw std::logic_error("you need to specify the component when RATE type is set in BC");
                         break;
                     }
 

opm/simulators/utils/ParallelRestart.cpp

@@ -65,6 +65,7 @@
 #include <opm/parser/eclipse/EclipseState/Schedule/Well/WList.hpp>
 #include <opm/parser/eclipse/EclipseState/Schedule/Well/WListManager.hpp>
 #include <opm/parser/eclipse/EclipseState/SimulationConfig/BCConfig.hpp>
+#include <opm/parser/eclipse/EclipseState/SimulationConfig/RockConfig.hpp>
 #include <opm/parser/eclipse/EclipseState/SimulationConfig/SimulationConfig.hpp>
 #include <opm/parser/eclipse/EclipseState/SimulationConfig/ThresholdPressure.hpp>
 #include <opm/parser/eclipse/EclipseState/SummaryConfig/SummaryConfig.hpp>
@@ -446,6 +447,7 @@ std::size_t packSize(const std::array<T,N>& data, Dune::MPIHelper::MPICommunicat
 
 HANDLE_AS_POD(Actdims)
 HANDLE_AS_POD(Aqudims)
+HANDLE_AS_POD(BCConfig::BCFace)
 HANDLE_AS_POD(data::Connection)
 HANDLE_AS_POD(data::Rates)
 HANDLE_AS_POD(data::Segment)
@@ -455,6 +457,7 @@ HANDLE_AS_POD(MLimits)
 HANDLE_AS_POD(PVTWRecord)
 HANDLE_AS_POD(PVCDORecord)
 HANDLE_AS_POD(Regdims)
+HANDLE_AS_POD(RockConfig::RockComp)
 HANDLE_AS_POD(ROCKRecord)
 HANDLE_AS_POD(SatFuncControls)
 HANDLE_AS_POD(Tabdims)
@@ -465,7 +468,6 @@ HANDLE_AS_POD(WellBrineProperties)
 HANDLE_AS_POD(Welldims)
 HANDLE_AS_POD(WellFoamProperties)
 HANDLE_AS_POD(WellSegmentDims)
-HANDLE_AS_POD(BCConfig::BCFace)
 
 std::size_t packSize(const data::Well& data, Dune::MPIHelper::MPICommunicator comm)
 {
@@ -520,6 +522,15 @@ std::size_t packSize(const BCConfig& bc, Dune::MPIHelper::MPICommunicator comm)
     return packSize(bc.faces(), comm);
 }
 
+std::size_t packSize(const RockConfig& data, Dune::MPIHelper::MPICommunicator comm)
+{
+    return packSize(data.active(), comm) +
+           packSize(data.rocknum_property(), comm) +
+           packSize(data.comp(), comm) +
+           packSize(data.num_rock_tables(), comm) +
+           packSize(data.water_compaction(), comm) +
+           packSize(data.hysteresis_mode(), comm);
+}
 
 std::size_t packSize(const NNC& data, Dune::MPIHelper::MPICommunicator comm)
 {
@@ -613,6 +624,7 @@ std::size_t packSize(const SimulationConfig& data, Dune::MPIHelper::MPICommunica
 {
     return packSize(data.getThresholdPressure(), comm) +
            packSize(data.bcconfig(), comm) +
+           packSize(data.rock_config(), comm) +
            packSize(data.useCPR(), comm) +
            packSize(data.hasDISGAS(), comm) +
            packSize(data.hasVAPOIL(), comm) +
@@ -2344,6 +2356,18 @@ void pack(const BCConfig& bc, std::vector<char>& buffer, int& position,
     pack(bc.faces(), buffer, position, comm);
 }
 
+void pack(const RockConfig& data, std::vector<char>& buffer, int& position,
+          Dune::MPIHelper::MPICommunicator comm)
+{
+    pack(data.active(), buffer, position, comm);
+    pack(data.rocknum_property(), buffer, position, comm);
+    pack(data.comp(), buffer, position, comm);
+    pack(data.num_rock_tables(), buffer, position, comm);
+    pack(data.water_compaction(), buffer, position, comm);
+    pack(data.hysteresis_mode(), buffer, position, comm);
+}
+
+
 void pack(const NNC& data, std::vector<char>& buffer, int& position,
           Dune::MPIHelper::MPICommunicator comm)
 {
@@ -2452,6 +2476,7 @@ void pack(const SimulationConfig& data, std::vector<char>& buffer, int& position
 {
     pack(data.getThresholdPressure(), buffer, position, comm);
     pack(data.bcconfig(), buffer, position, comm);
+    pack(data.rock_config(), buffer, position, comm);
     pack(data.useCPR(), buffer, position, comm);
     pack(data.hasDISGAS(), buffer, position, comm);
     pack(data.hasVAPOIL(), buffer, position, comm);
@@ -4330,6 +4355,26 @@ void unpack(RestartValue& data, std::vector<char>& buffer, int& position,
     unpack(data.extra, buffer, position, comm);
 }
 
+void unpack(RockConfig& data, std::vector<char>& buffer, int& position,
+            Dune::MPIHelper::MPICommunicator comm)
+{
+    RockConfig rock_config;
+    bool active;
+    std::vector<RockConfig::RockComp> rock_comp;
+    std::string rocknum_property;
+    std::size_t num_rock_tables;
+    bool water_compaction;
+    RockConfig::Hysteresis hyst_mode;
+
+    unpack(active, buffer, position, comm);
+    unpack(rocknum_property, buffer, position, comm);
+    unpack(rock_comp, buffer, position, comm);
+    unpack(num_rock_tables, buffer, position, comm);
+    unpack(water_compaction, buffer, position, comm);
+    unpack(hyst_mode, buffer, position, comm);
+    data = RockConfig(active, rock_comp, rocknum_property, num_rock_tables, water_compaction, hyst_mode);
+}
+
 void unpack(ThresholdPressure& data, std::vector<char>& buffer, int& position,
             Dune::MPIHelper::MPICommunicator comm)
 {
@@ -4500,14 +4545,16 @@ void unpack(SimulationConfig& data, std::vector<char>& buffer, int& position,
 {
     ThresholdPressure thresholdPressure;
     BCConfig bc;
+    RockConfig rock_config;
     bool useCPR, DISGAS, VAPOIL, isThermal;
     unpack(thresholdPressure, buffer, position, comm);
     unpack(bc, buffer, position, comm);
+    unpack(rock_config, buffer, position, comm);
     unpack(useCPR, buffer, position, comm);
     unpack(DISGAS, buffer, position, comm);
     unpack(VAPOIL, buffer, position, comm);
     unpack(isThermal, buffer, position, comm);
-    data = SimulationConfig(thresholdPressure, bc, useCPR, DISGAS, VAPOIL, isThermal);
+    data = SimulationConfig(thresholdPressure, bc, rock_config, useCPR, DISGAS, VAPOIL, isThermal);
 }
 
 void unpack(TimeMap& data, std::vector<char>& buffer, int& position,

opm/simulators/utils/ParallelRestart.hpp

@@ -123,6 +123,7 @@ class Regdims;
 class RestartConfig;
 class RestartSchedule;
 class RFTConfig;
+class RockConfig;
 class ROCKRecord;
 class RockTable;
 class Rock2dTable;
@@ -726,6 +727,8 @@ ADD_PACK_PROTOTYPES(RestartKey)
 ADD_PACK_PROTOTYPES(RestartSchedule)
 ADD_PACK_PROTOTYPES(RestartValue)
 ADD_PACK_PROTOTYPES(RFTConfig)
+ADD_PACK_PROTOTYPES(RockConfig)
+ADD_PACK_PROTOTYPES(RockConfig::RockComp)
 ADD_PACK_PROTOTYPES(ROCKRecord)
 ADD_PACK_PROTOTYPES(RockTable)
 ADD_PACK_PROTOTYPES(Rock2dTable)

tests/test_ParallelRestart.cpp

@@ -74,6 +74,7 @@
 #include <opm/parser/eclipse/EclipseState/Schedule/Well/WList.hpp>
 #include <opm/parser/eclipse/EclipseState/Schedule/Well/WListManager.hpp>
 #include <opm/parser/eclipse/EclipseState/SimulationConfig/BCConfig.hpp>
+#include <opm/parser/eclipse/EclipseState/SimulationConfig/RockConfig.hpp>
 #include <opm/parser/eclipse/EclipseState/SimulationConfig/SimulationConfig.hpp>
 #include <opm/parser/eclipse/EclipseState/SimulationConfig/ThresholdPressure.hpp>
 #include <opm/parser/eclipse/EclipseState/SummaryConfig/SummaryConfig.hpp>
@@ -183,6 +184,12 @@ Opm::ThresholdPressure getThresholdPressure()
                                   {{{1,2},{false,3.0}},{{2,3},{true,4.0}}});
 }
 
+Opm::RockConfig getRockConfig()
+{
+    return Opm::RockConfig(true, {{100, 0.25}, {200, 0.30}}, "ROCKNUM", 10, false, Opm::RockConfig::Hysteresis::HYSTER);
+}
+
+
 
 Opm::BCConfig getBCConfig()
 {
@@ -604,6 +611,16 @@ BOOST_AUTO_TEST_CASE(ThresholdPressure)
 #endif
 }
 
+BOOST_AUTO_TEST_CASE(RockConfig)
+{
+#if HAVE_MPI
+    Opm::RockConfig val1 = getRockConfig();
+    auto val2 = PackUnpack(val1);
+    BOOST_CHECK(std::get<1>(val2) == std::get<2>(val2));
+    BOOST_CHECK(val1 == std::get<0>(val2));
+#endif
+}
+
 
 BOOST_AUTO_TEST_CASE(EDITNNC)
 {
@@ -774,7 +791,7 @@ BOOST_AUTO_TEST_CASE(InitConfig)
 BOOST_AUTO_TEST_CASE(SimulationConfig)
 {
 #if HAVE_MPI
-    Opm::SimulationConfig val1(getThresholdPressure(), getBCConfig(), false, true, false, true);
+    Opm::SimulationConfig val1(getThresholdPressure(), getBCConfig(), getRockConfig(), false, true, false, true);
     auto val2 = PackUnpack(val1);
     BOOST_CHECK(std::get<1>(val2) == std::get<2>(val2));
     BOOST_CHECK(val1 == std::get<0>(val2));