From 4bb5213f3765802e27e6eb2b2a40aabc05d1b3fe Mon Sep 17 00:00:00 2001
From: Tor Harald Sandve <tor.harald.sandve@iris.no>
Date: Wed, 8 Nov 2017 12:17:31 +0100
Subject: [PATCH] BUGFIX Update the solution variable in ebos

Updating the solution variable in updateState() was conceptually wrong
and lead to incorrect results if the linear solver failed before the
updateState() method was called. i.e. in the first iteration.
---
 opm/autodiff/BlackoilModelEbos.hpp | 31 +++++++++++++++---------------
 1 file changed, 15 insertions(+), 16 deletions(-)

diff --git a/opm/autodiff/BlackoilModelEbos.hpp b/opm/autodiff/BlackoilModelEbos.hpp
index 237b624b8..379e4402e 100644
--- a/opm/autodiff/BlackoilModelEbos.hpp
+++ b/opm/autodiff/BlackoilModelEbos.hpp
@@ -203,7 +203,7 @@ namespace Opm {
         /// \param[in] timer                  simulation timer
         /// \param[in, out] reservoir_state   reservoir state variables
         /// \param[in, out] well_state        well state variables
-        void prepareStep(const SimulatorTimerInterface& /*timer*/,
+        void prepareStep(const SimulatorTimerInterface& timer,
                          const ReservoirState& /*reservoir_state*/,
                          const WellState& /* well_state */)
         {
@@ -211,6 +211,18 @@ namespace Opm {
                 updateRateConverter();
             }
 
+            // update the solution variables in ebos
+
+            // if the last time step failed we need to update the curent solution
+            // and recalculate the Intesive Quantities.
+            if ( timer.lastStepFailed() ) {
+                ebosSimulator_.model().solution( 0 /* timeIdx */ ) = ebosSimulator_.model().solution( 1 /* timeIdx */ );
+                ebosSimulator_.model().invalidateIntensiveQuantitiesCache(/*timeIdx=*/0);
+            } else {
+                // set the initial solution.
+                ebosSimulator_.model().solution( 1 /* timeIdx */ ) = ebosSimulator_.model().solution( 0 /* timeIdx */ );
+            }
+
             unsigned numDof = ebosSimulator_.model().numGridDof();
             wasSwitched_.resize(numDof);
             std::fill(wasSwitched_.begin(), wasSwitched_.end(), false);
@@ -329,7 +341,7 @@ namespace Opm {
 
                 // Apply the update, with considering model-dependent limitations and
                 // chopping of the update.
-                updateState(x,iteration);
+                updateState(x);
 
                 report.update_time += perfTimer.stop();
             }
@@ -605,8 +617,7 @@ namespace Opm {
         /// \param[in]      dx                updates to apply to primary variables
         /// \param[in, out] reservoir_state   reservoir state variables
         /// \param[in, out] well_state        well state variables
-        void updateState(const BVector& dx,
-                         const int iterationIdx)
+        void updateState(const BVector& dx)
         {
             using namespace Opm::AutoDiffGrid;
 
@@ -618,12 +629,6 @@ namespace Opm {
             const auto& elemEndIt = gridView.template end</*codim=*/0>();
             SolutionVector& solution = ebosSimulator_.model().solution( 0 /* timeIdx */ );
 
-            // Store the initial solution.
-            if( iterationIdx == 0 )
-            {
-                ebosSimulator_.model().solution( 1 /* timeIdx */ ) = solution;
-            }
-
             for (auto elemIt = gridView.template begin</*codim=*/0>();
                  elemIt != elemEndIt;
                  ++elemIt)
@@ -1611,12 +1616,6 @@ namespace Opm {
             {
                 ebosSimulator_.problem().beginTimeStep();
             }
-            // if the last time step failed we need to update the solution varables in ebos
-            // and recalculate the Intesive Quantities.
-            if ( timer.lastStepFailed() && iterationIdx == 0  ) {
-                ebosSimulator_.model().solution( 0 /* timeIdx */ ) = ebosSimulator_.model().solution( 1 /* timeIdx */ );
-                ebosSimulator_.model().invalidateIntensiveQuantitiesCache(/*timeIdx=*/0);
-            }
 
             ebosSimulator_.problem().beginIteration();
             ebosSimulator_.model().linearizer().linearize();