Fix bug for time-stepping when loading restart deck

CMakeLists_files.cmake

@@ -495,6 +495,7 @@ if(ENABLE_ECL_OUTPUT)
           tests/include_swof.txt
           tests/include_grid_3x5x4.grdecl
           tests/SPE1CASE2.DATA
+          tests/SPE1CASE2_RESTART.DATA
           tests/SPE1CASE2_RESTART_SKIPREST.DATA
           tests/SPE1CASE2.X0060
           tests/PYACTION.DATA

src/opm/parser/eclipse/EclipseState/Schedule/ScheduleDeck.cpp

@@ -202,15 +202,20 @@ ScheduleDeck::ScheduleDeck(const Deck& deck, const ScheduleRestartInfo& rst_info
     this->skiprest = rst_info.skiprest;
     if (this->m_restart_offset > 0) {
         for (std::size_t it = 0; it < this->m_restart_offset; it++) {
+            if (it == 0)
+                this->m_blocks.emplace_back(KeywordLocation{}, ScheduleTimeType::START, start_time);
+            else
                 this->m_blocks.emplace_back(KeywordLocation{}, ScheduleTimeType::RESTART, start_time);
-            if (it < this->m_restart_offset - 1)
             this->m_blocks.back().end_time(start_time);
         }
+        if (!this->skiprest) {
+            this->m_blocks.back().end_time(this->m_restart_time);
+            this->m_blocks.emplace_back(KeywordLocation{}, ScheduleTimeType::RESTART, this->m_restart_time);
+        }
     } else
         this->m_blocks.emplace_back(KeywordLocation{}, ScheduleTimeType::START, start_time);
 
-
+    ScheduleDeckContext context(this->skiprest, this->m_blocks.back().start_time());
-    ScheduleDeckContext context(this->m_restart_offset > 0, start_time);
     for( const auto& keyword : SCHEDULESection(deck)) {
         if (keyword.name() == "DATES") {
             for (size_t recordIndex = 0; recordIndex < keyword.size(); recordIndex++) {

tests/SPE1CASE2_RESTART.DATA

@@ -0,0 +1,370 @@
+-- This reservoir simulation deck is made available under the Open Database
+-- License: http://opendatacommons.org/licenses/odbl/1.0/. Any rights in
+-- individual contents of the database are licensed under the Database Contents
+-- License: http://opendatacommons.org/licenses/dbcl/1.0/
+
+-- Copyright (C) 2015 Statoil
+
+-- This simulation is based on the data given in 
+-- 'Comparison of Solutions to a Three-Dimensional
+-- Black-Oil Reservoir Simulation Problem' by Aziz S. Odeh,
+-- Journal of Petroleum Technology, January 1981
+
+---------------------------------------------------------------------------
+------------------------ SPE1 - CASE 2 ------------------------------------
+---------------------------------------------------------------------------
+
+RUNSPEC
+-- -------------------------------------------------------------------------
+
+TITLE
+   SPE1 - CASE 2
+
+DIMENS
+   10 10 3 /
+
+-- The number of equilibration regions is inferred from the EQLDIMS
+-- keyword.
+EQLDIMS
+/
+
+-- The number of PVTW tables is inferred from the TABDIMS keyword;
+-- when no data is included in the keyword the default values are used.
+TABDIMS
+/
+
+
+OIL
+GAS
+WATER
+DISGAS
+-- As seen from figure 4 in Odeh, GOR is increasing with time,
+-- which means that dissolved gas is present
+
+
+FIELD
+
+START
+   1 'JAN' 2015 /
+
+WELLDIMS
+-- Item 1: maximum number of wells in the model
+-- 	   - there are two wells in the problem; injector and producer
+-- Item 2: maximum number of grid blocks connected to any one well
+-- 	   - must be one as the wells are located at specific grid blocks
+-- Item 3: maximum number of groups in the model
+-- 	   - we are dealing with only one 'group'
+-- Item 4: maximum number of wells in any one group
+-- 	   - there must be two wells in a group as there are two wells in total
+   2 1 1 2 /
+
+UNIFOUT
+
+GRID
+
+-- The INIT keyword is used to request an .INIT file. The .INIT file
+-- is written before the simulation actually starts, and contains grid
+-- properties and saturation tables as inferred from the input
+-- deck. There are no other keywords which can be used to configure
+-- exactly what is written to the .INIT file.
+INIT
+
+-- -------------------------------------------------------------------------
+NOECHO
+
+DX 
+-- There are in total 300 cells with length 1000ft in x-direction	
+   	300*1000 /
+DY
+-- There are in total 300 cells with length 1000ft in y-direction	
+	300*1000 /
+DZ
+-- The layers are 20, 30 and 50 ft thick, in each layer there are 100 cells
+	100*20 100*30 100*50 /
+
+TOPS
+-- The depth of the top of each grid block
+	100*8325 /
+
+PORO
+-- Constant porosity of 0.3 throughout all 300 grid cells
+   	300*0.3 /
+
+PERMX
+-- The layers have perm. 500mD, 50mD and 200mD, respectively.
+	100*500 100*50 100*200 /
+
+PERMY
+-- Equal to PERMX
+	100*500 100*50 100*200 /
+
+PERMZ
+-- Cannot find perm. in z-direction in Odeh's paper
+-- For the time being, we will assume PERMZ equal to PERMX and PERMY:
+	100*500 100*50 100*200 /
+ECHO
+
+PROPS
+-- -------------------------------------------------------------------------
+
+PVTW
+-- Item 1: pressure reference (psia)
+-- Item 2: water FVF (rb per bbl or rb per stb)
+-- Item 3: water compressibility (psi^{-1})
+-- Item 4: water viscosity (cp)
+-- Item 5: water 'viscosibility' (psi^{-1})
+
+-- Using values from Norne:
+-- In METRIC units:
+-- 	277.0 1.038 4.67E-5 0.318 0.0 /
+-- In FIELD units:
+    	4017.55 1.038 3.22E-6 0.318 0.0 /
+
+ROCK
+-- Item 1: reference pressure (psia)
+-- Item 2: rock compressibility (psi^{-1})
+
+-- Using values from table 1 in Odeh:
+	14.7 3E-6 /
+
+SWOF
+-- Column 1: water saturation
+--   	     - this has been set to (almost) equally spaced values from 0.12 to 1
+-- Column 2: water relative permeability
+--   	     - generated from the Corey-type approx. formula
+--	       the coeffisient is set to 10e-5, S_{orw}=0 and S_{wi}=0.12
+-- Column 3: oil relative permeability when only oil and water are present
+--	     - we will use the same values as in column 3 in SGOF.
+-- 	       This is not really correct, but since only the first 
+--	       two values are of importance, this does not really matter
+-- Column 4: corresponding water-oil capillary pressure (psi) 
+
+0.12	0    		 	1	0
+0.18	4.64876033057851E-008	1	0
+0.24	0.000000186		0.997	0
+0.3	4.18388429752066E-007	0.98	0
+0.36	7.43801652892562E-007	0.7	0
+0.42	1.16219008264463E-006	0.35	0
+0.48	1.67355371900826E-006	0.2	0
+0.54	2.27789256198347E-006	0.09	0
+0.6	2.97520661157025E-006	0.021	0
+0.66	3.7654958677686E-006	0.01	0
+0.72	4.64876033057851E-006	0.001	0
+0.78	0.000005625		0.0001	0
+0.84	6.69421487603306E-006	0	0
+0.91	8.05914256198347E-006	0	0
+1	0.00001			0	0 /
+
+
+SGOF
+-- Column 1: gas saturation
+-- Column 2: gas relative permeability
+-- Column 3: oil relative permeability when oil, gas and connate water are present
+-- Column 4: oil-gas capillary pressure (psi)
+-- 	     - stated to be zero in Odeh's paper
+
+-- Values in column 1-3 are taken from table 3 in Odeh's paper:
+0	0	1	0
+0.001	0	1	0
+0.02	0	0.997	0
+0.05	0.005	0.980	0
+0.12	0.025	0.700	0
+0.2	0.075	0.350	0
+0.25	0.125	0.200	0
+0.3	0.190	0.090	0
+0.4	0.410	0.021	0
+0.45	0.60	0.010	0
+0.5	0.72	0.001	0
+0.6	0.87	0.0001	0
+0.7	0.94	0.000	0
+0.85	0.98	0.000	0 
+0.88	0.984	0.000	0 /
+--1.00	1.0	0.000	0 /
+-- Warning from Eclipse: first sat. value in SWOF + last sat. value in SGOF
+-- 	   		 must not be greater than 1, but Eclipse still runs
+-- Flow needs the sum to be excactly 1 so I added a row with gas sat. =  0.88
+-- The corresponding krg value was estimated by assuming linear rel. between
+-- gas sat. and krw. between gas sat. 0.85 and 1.00 (the last two values given)
+
+DENSITY
+-- Density (lb per ft³) at surface cond. of 
+-- oil, water and gas, respectively (in that order)
+
+-- Using values from Norne:
+-- In METRIC units:
+--      859.5 1033.0 0.854 /
+-- In FIELD units:
+      	53.66 64.49 0.0533 /
+
+PVDG
+-- Column 1: gas phase pressure (psia)
+-- Column 2: gas formation volume factor (rb per Mscf)
+-- 	     - in Odeh's paper the units are said to be given in rb per bbl, 
+-- 	       but this is assumed to be a mistake: FVF-values in Odeh's paper 
+--	       are given in rb per scf, not rb per bbl. This will be in 
+--	       agreement with conventions
+-- Column 3: gas viscosity (cP)
+
+-- Using values from lower right table in Odeh's table 2:
+14.700	166.666	0.008000
+264.70	12.0930	0.009600
+514.70	6.27400	0.011200
+1014.7	3.19700	0.014000
+2014.7	1.61400	0.018900
+2514.7	1.29400	0.020800
+3014.7	1.08000	0.022800
+4014.7	0.81100	0.026800
+5014.7	0.64900	0.030900
+9014.7	0.38600	0.047000 /
+
+PVTO
+-- Column 1: dissolved gas-oil ratio (Mscf per stb)
+-- Column 2: bubble point pressure (psia)
+-- Column 3: oil FVF for saturated oil (rb per stb)
+-- Column 4: oil viscosity for saturated oil (cP)
+
+-- Use values from top left table in Odeh's table 2:
+0.0010	14.7	1.0620	1.0400 /
+0.0905	264.7	1.1500	0.9750 /
+0.1800	514.7	1.2070	0.9100 /
+0.3710	1014.7	1.2950	0.8300 /
+0.6360	2014.7	1.4350	0.6950 /
+0.7750	2514.7	1.5000	0.6410 /
+0.9300	3014.7	1.5650	0.5940 /
+1.2700	4014.7	1.6950	0.5100 
+	9014.7	1.5790	0.7400 /
+1.6180	5014.7	1.8270	0.4490 
+	9014.7	1.7370	0.6310 /	
+-- It is required to enter data for undersaturated oil for the highest GOR
+-- (i.e. the last row) in the PVTO table.
+-- In order to fulfill this requirement, values for oil FVF and viscosity
+-- at 9014.7psia and GOR=1.618 for undersaturated oil have been approximated:
+-- It has been assumed that there is a linear relation between the GOR
+-- and the FVF when keeping the pressure constant at 9014.7psia.
+-- From Odeh we know that (at 9014.7psia) the FVF is 2.357 at GOR=2.984
+-- for saturated oil and that the FVF is 1.579 at GOR=1.27 for undersaturated oil,
+-- so it is possible to use the assumption described above. 
+-- An equivalent approximation for the viscosity has been used.
+/
+
+SOLUTION
+-- -------------------------------------------------------------------------
+
+RESTART
+   'SPE1CASE2' 60 /
+
+RSVD
+-- Dissolved GOR is initially constant with depth through the reservoir.
+-- The reason is that the initial reservoir pressure given is higher 
+---than the bubble point presssure of 4014.7psia, meaning that there is no 
+-- free gas initially present.
+8300 1.270
+8450 1.270 /
+
+SUMMARY
+-- -------------------------------------------------------------------------	 
+
+-- 1a) Oil rate vs time
+FOPR
+-- Field Oil Production Rate
+
+-- 1b) GOR vs time
+WGOR
+-- Well Gas-Oil Ratio
+   'PROD'
+/
+-- Using FGOR instead of WGOR:PROD results in the same graph
+FGOR
+
+-- 2a) Pressures of the cell where the injector and producer are located
+BPR
+1  1  1 /
+10 10 3 /
+/
+
+-- 2b) Gas saturation at grid points given in Odeh's paper
+BGSAT
+1  1  1 /
+1  1  2 /
+1  1  3 /
+10 1  1 /
+10 1  2 /
+10 1  3 /
+10 10 1 /
+10 10 2 /
+10 10 3 /
+/
+
+-- In order to compare Eclipse with Flow:
+WBHP
+  'INJ'
+  'PROD'
+/
+WGIR
+  'INJ'
+  'PROD'
+/
+WGIT
+  'INJ'
+  'PROD'
+/
+WGPR
+  'INJ'
+  'PROD'
+/
+WGPT
+  'INJ'
+  'PROD'
+/
+WOIR
+  'INJ'
+  'PROD'
+/
+WOIT
+  'INJ'
+  'PROD'
+/
+WOPR
+  'INJ'
+  'PROD'
+/
+WOPT
+  'INJ'
+  'PROD'
+/
+WWIR
+  'INJ'
+  'PROD'
+/
+WWIT
+  'INJ'
+  'PROD'
+/
+WWPR
+  'INJ'
+  'PROD'
+/
+WWPT
+  'INJ'
+  'PROD'
+/
+SCHEDULE
+-- -------------------------------------------------------------------------
+RPTSCHED
+	'PRES' 'SGAS' 'RS' 'WELLS' /
+
+RPTRST
+	'BASIC=1' /
+
+
+TSTEP 
+31 28 31 30 31 30 31 31 30 31 30 31 
+31 28 31 30 31 30 31 31 30 31 30 31 
+31 28 31 30 31 30 31 31 30 31 30 31 
+31 28 31 30 31 30 31 31 30 31 30 31 
+31 28 31 30 31 30 31 31 30 31 30 31 /
+
+--Advance the simulator once a year for TEN years:
+--10*365 /
+
+END

tests/parser/ScheduleRestartTests.cpp

@@ -93,19 +93,42 @@ BOOST_AUTO_TEST_CASE(LoadRST) {
     }
 }
 
-BOOST_AUTO_TEST_CASE(LoadRestartSim) {
+void compare_sched(const std::string& base_deck,
-    auto python = std::make_shared<Python>();
+                   const std::string& rst_deck,
+                   const std::string& rst_fname,
+                   std::size_t restart_step)
+{
     Parser parser;
-    auto deck = parser.parseFile("SPE1CASE2.DATA");
+    auto python = std::make_shared<Python>();
+    auto deck = parser.parseFile(base_deck);
     EclipseState ecl_state(deck);
     Schedule sched(deck, ecl_state, python);
 
-    auto restart_deck = parser.parseFile("SPE1CASE2_RESTART_SKIPREST.DATA");
+    auto restart_deck = parser.parseFile(rst_deck);
-    auto rst_file = std::make_shared<EclIO::ERst>("SPE1CASE2.X0060");
+    auto rst_file = std::make_shared<EclIO::ERst>(rst_fname);
-    auto rst_view = std::make_shared<EclIO::RestartFileView>(std::move(rst_file), 60);
+    auto rst_view = std::make_shared<EclIO::RestartFileView>(std::move(rst_file), restart_step);
     auto rst_state = RestartIO::RstState::load(std::move(rst_view));
     EclipseState ecl_state_restart(restart_deck);
     Schedule restart_sched(restart_deck, ecl_state_restart, python, {}, &rst_state);
 
-    // Verify that sched and restart_sched are identical from report_step 60 and onwords.
+    BOOST_CHECK_EQUAL(restart_sched.size(), sched.size());
+    for (std::size_t report_step=restart_step; report_step < sched.size(); report_step++) {
+        const auto& base = sched[report_step];
+        auto rst = restart_sched[report_step];
+
+        BOOST_CHECK(base.start_time() == rst.start_time());
+        if (report_step < sched.size() - 1)
+            BOOST_CHECK(base.end_time() == rst.end_time());
+
+        // Should ideally do a base == rst check here, but for now the members
+        // wells, rft_config, m_first_in_year and m_first_in_month fail.
+        // BOOST_CHECK(base == rst);
+    }
+}
+
+
+
+BOOST_AUTO_TEST_CASE(LoadRestartSim) {
+    compare_sched("SPE1CASE2.DATA", "SPE1CASE2_RESTART_SKIPREST.DATA", "SPE1CASE2.X0060", 60);
+    compare_sched("SPE1CASE2.DATA", "SPE1CASE2_RESTART.DATA", "SPE1CASE2.X0060", 60);
 }