Merge pull request #1134 from andlaus/report-failure

flow_ebos: print statistics about failed time steps

opm/autodiff/BlackoilModelBase.hpp

@@ -295,6 +295,13 @@ namespace Opm {
         void applyVREPGroupControl(const ReservoirState& reservoir_state,
                                    WellState& well_state);
 
+        /// return the statistics if the nonlinearIteration() method failed.
+        ///
+        /// NOTE: for the flow_legacy simulator family this method is a stub, i.e. the
+        /// failure report object will *not* contain any meaningful data.
+        const SimulatorReport& failureReport() const
+        { return failureReport_; }
+
     protected:
 
         // ---------  Types and enums  ---------
@@ -307,6 +314,7 @@ namespace Opm {
 
         // ---------  Data members  ---------
 
+        SimulatorReport failureReport_;
         const Grid&         grid_;
         const BlackoilPropsAdFromDeck& fluid_;
         const DerivedGeology&           geo_;

opm/autodiff/BlackoilModelEbos.hpp

@@ -54,7 +54,6 @@
 #include <opm/common/OpmLog/OpmLog.hpp>
 #include <opm/parser/eclipse/Units/Units.hpp>
 #include <opm/core/well_controls.h>
-#include <opm/core/simulator/SimulatorReport.hpp>
 #include <opm/simulators/timestepping/SimulatorTimer.hpp>
 #include <opm/core/utility/parameters/ParameterGroup.hpp>
 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
@@ -236,6 +235,7 @@ namespace Opm {
                                            WellState& well_state)
         {
             SimulatorReport report;
+            failureReport_ = SimulatorReport();
             Dune::Timer perfTimer;
 
             perfTimer.start();
@@ -255,6 +255,7 @@ namespace Opm {
             }
             catch (...) {
                 report.assemble_time += perfTimer.stop();
+                failureReport_ += report;
                 // todo (?): make the report an attribute of the class
                 throw; // continue throwing the stick
             }
@@ -294,7 +295,8 @@ namespace Opm {
                 catch (...) {
                     report.linear_solve_time += perfTimer.stop();
                     report.total_linear_iterations += linearIterationsLastSolve();
-                    // todo (?): make the report an attribute of the class
+
+                    failureReport_ += report;
                     throw; // re-throw up
                 }
 
@@ -1429,6 +1431,10 @@ namespace Opm {
         const Simulator& ebosSimulator() const
         { return ebosSimulator_; }
 
+        /// return the statistics if the nonlinearIteration() method failed
+        const SimulatorReport& failureReport() const
+        { return failureReport_; }
+
     protected:
         const ISTLSolverType& istlSolver() const
         {
@@ -1452,6 +1458,7 @@ namespace Opm {
         const bool has_vapoil_;
 
         ModelParameters                 param_;
+        SimulatorReport failureReport_;
 
         // Well Model
         StandardWellsDense<FluidSystem, BlackoilIndices, ElementContext, MaterialLaw> well_model_;

opm/autodiff/BlackoilSequentialModel.hpp

@@ -282,11 +282,19 @@ namespace Opm {
             return transport_solver_.model().getFIPData();
         }
 
+        /// return the statistics if the nonlinearIteration() method failed.
+        ///
+        /// NOTE: for the flow_legacy simulator family this method is a stub, i.e. the
+        /// failure report object will *not* contain any meaningful data.
+        const SimulatorReport& failureReport() const
+        { return failureReport_; }
 
     protected:
         typedef NonlinearSolver<PressureModel> PressureSolver;
         typedef NonlinearSolver<TransportModel> TransportSolver;
 
+        SimulatorReport failureReport_;
+
         std::unique_ptr<PressureModel> pressure_model_;
         std::unique_ptr<TransportModel> transport_model_;
         PressureSolver pressure_solver_;

opm/autodiff/FlowMainEbos.hpp

@@ -643,12 +643,13 @@ namespace Opm
                     OpmLog::info(msg);
                 }
 
-                SimulatorReport fullReport = simulator_->run(simtimer, *state_);
+                SimulatorReport successReport = simulator_->run(simtimer, *state_);
+                SimulatorReport failureReport = simulator_->failureReport();
 
                 if (output_cout_) {
                     std::ostringstream ss;
                     ss << "\n\n================    End of simulation     ===============\n\n";
-                    fullReport.reportFullyImplicit(ss);
+                    successReport.reportFullyImplicit(ss, &failureReport);
                     OpmLog::info(ss.str());
                     if (param_.anyUnused()) {
                         // This allows a user to catch typos and misunderstandings in the
@@ -662,7 +663,7 @@ namespace Opm
                 if (output_to_files_) {
                     std::string filename = output_dir_ + "/walltime.txt";
                     std::fstream tot_os(filename.c_str(), std::fstream::trunc | std::fstream::out);
-                    fullReport.reportParam(tot_os);
+                    successReport.reportParam(tot_os);
                 }
             } else {
                 if (output_cout_) {

opm/autodiff/NonlinearSolver.hpp

@@ -99,6 +99,9 @@ namespace Opm {
              ReservoirState& reservoir_state,
              WellState& well_state);
 
+        /// return the statistics if the step() method failed
+        const SimulatorReport& failureReport() const
+        { return failureReport_; }
 
         /// Number of linearizations used in all calls to step().
         int linearizations() const;
@@ -173,6 +176,7 @@ namespace Opm {
 
     private:
         // ---------  Data members  ---------
+        SimulatorReport failureReport_;
         SolverParameters param_;
         std::unique_ptr<PhysicalModel> model_;
         int linearizations_;

opm/autodiff/NonlinearSolver_impl.hpp

@@ -124,6 +124,7 @@ namespace Opm
     {
         SimulatorReport iterReport;
         SimulatorReport report;
+        failureReport_ = SimulatorReport();
 
         // Do model-specific once-per-step calculations.
         model_->prepareStep(timer, initial_reservoir_state, initial_well_state);
@@ -137,6 +138,7 @@ namespace Opm
 
         // ----------  Main nonlinear solver loop  ----------
         do {
+            try {
                 // Do the nonlinear step. If we are in a converged state, the
                 // model will usually do an early return without an expensive
                 // solve, unless the minIter() count has not been reached yet.
@@ -147,9 +149,19 @@ namespace Opm
 
                 converged = report.converged;
                 iteration += 1;
+            }
+            catch (...) {
+                // if an iteration fails during a time step, all previous iterations
+                // count as a failure as well
+                failureReport_ += report;
+                failureReport_ += model_->failureReport();
+                throw;
+            }
         } while ( (!converged && (iteration <= maxIter())) || (iteration <= minIter()));
 
         if (!converged) {
+            failureReport_ += report;
+
             std::string msg = "Failed to complete a time step within " + std::to_string(maxIter()) + " iterations.";
             if (model_->terminalOutputEnabled()) {
                 OpmLog::problem(msg);

opm/autodiff/SimulatorFullyImplicitBlackoilEbos.hpp

@@ -147,6 +147,9 @@ public:
         WellState prev_well_state;
 
         ExtraData extra;
+
+        failureReport_ = SimulatorReport();
+
         if (output_writer_.isRestart()) {
             // This is a restart, populate WellState and ReservoirState state objects from restart file
             output_writer_.initFromRestartFile(props_.phaseUsage(), grid(), state, prev_well_state, extra);
@@ -319,11 +322,13 @@ public:
                 stepReport = adaptiveTimeStepping->step( timer, *solver, state, well_state, event, output_writer_,
                                                          output_writer_.requireFIPNUM() ? &fipnum : nullptr );
                 report += stepReport;
+                failureReport_ += adaptiveTimeStepping->failureReport();
             }
             else {
                 // solve for complete report step
                 stepReport = solver->step(timer, state, well_state);
                 report += stepReport;
+                failureReport_ += solver->failureReport();
 
                 if( terminal_output_ )
                 {
@@ -407,6 +412,10 @@ public:
         return report;
     }
 
+    /** \brief Returns the simulator report for the failed substeps of the simulation.
+     */
+    const SimulatorReport& failureReport() const { return failureReport_; };
+
     const Grid& grid() const
     { return ebosSimulator_.gridManager().grid(); }
 
@@ -830,6 +839,8 @@ protected:
     typedef RateConverter::SurfaceToReservoirVoidage<FluidSystem, std::vector<int> > RateConverterType;
     typedef typename Solver::SolverParameters SolverParameters;
 
+    SimulatorReport failureReport_;
+
     const parameter::ParameterGroup param_;
     ModelParameters model_param_;
     SolverParameters solver_param_;

opm/polymer/SimulatorCompressiblePolymer.hpp

@@ -20,6 +20,7 @@
 #ifndef OPM_SIMULATORCOMPRESSIBLEPOLYMER_HEADER_INCLUDED
 #define OPM_SIMULATORCOMPRESSIBLEPOLYMER_HEADER_INCLUDED
 
+#include <opm/core/simulator/SimulatorReport.hpp>
 #include <boost/shared_ptr.hpp>
 #include <vector>
 
@@ -91,8 +92,18 @@ namespace Opm
                             PolymerBlackoilState& state,
                             WellState& well_state);
 
+        /// return the statistics if the nonlinearIteration() method failed.
+        ///
+        /// NOTE: for the flow_legacy simulator family this method is a stub, i.e. the
+        /// failure report object will *not* contain any meaningful data.
+        const SimulatorReport& failureReport() const
+        { return failureReport_; }
+
     private:
         class Impl;
+
+        SimulatorReport failureReport_;
+
         // Using shared_ptr instead of scoped_ptr since scoped_ptr requires complete type for Impl.
         boost::shared_ptr<Impl> pimpl_;
     };

opm/polymer/SimulatorPolymer.hpp

@@ -20,6 +20,7 @@
 #ifndef OPM_SIMULATORPOLYMER_HEADER_INCLUDED
 #define OPM_SIMULATORPOLYMER_HEADER_INCLUDED
 
+#include <opm/core/simulator/SimulatorReport.hpp>
 #include <boost/shared_ptr.hpp>
 #include <vector>
 
@@ -95,8 +96,18 @@ namespace Opm
                             PolymerState& state,
                             WellState& well_state);
 
+        /// return the statistics if the nonlinearIteration() method failed.
+        ///
+        /// NOTE: for the flow_legacy simulator family this method is a stub, i.e. the
+        /// failure report object will *not* contain any meaningful data.
+        const SimulatorReport& failureReport() const
+        { return failureReport_; }
+
     private:
         class Impl;
+
+        SimulatorReport failureReport_;
+
         // Using shared_ptr instead of scoped_ptr since scoped_ptr requires complete type for Impl.
         boost::shared_ptr<Impl> pimpl_;
     };

opm/polymer/fullyimplicit/FullyImplicitCompressiblePolymerSolver.hpp

@@ -30,6 +30,7 @@
 #include <opm/polymer/PolymerProperties.hpp>
 #include <opm/polymer/fullyimplicit/WellStateFullyImplicitBlackoilPolymer.hpp>
 #include <opm/polymer/fullyimplicit/PolymerPropsAd.hpp>
+#include <opm/core/simulator/SimulatorReport.hpp>
 #include <opm/core/utility/parameters/ParameterGroup.hpp>
 #include <opm/simulators/timestepping/SimulatorTimerInterface.hpp>
 #include <opm/common/data/SimulationDataContainer.hpp>
@@ -175,6 +176,13 @@ namespace Opm {
         computeFluidInPlace(const PolymerBlackoilState& x,
                             const std::vector<int>& fipnum);
 
+        /// return the statistics if the nonlinearIteration() method failed.
+        ///
+        /// NOTE: for the flow_legacy simulator family this method is a stub, i.e. the
+        /// failure report object will *not* contain any meaningful data.
+        const SimulatorReport& failureReport() const
+        { return failureReport_; }
+
     private:
 
         struct SolutionState {
@@ -197,6 +205,7 @@ namespace Opm {
                Oil   = Opm::Oil  };
 
         // Member data
+        SimulatorReport failureReport_;
         const UnstructuredGrid&         grid_;
         const BlackoilPropsAdFromDeck& fluid_;
         const DerivedGeology&           geo_;

opm/polymer/fullyimplicit/SimulatorFullyImplicitCompressiblePolymer.hpp

@@ -21,6 +21,7 @@
 #ifndef OPM_SIMULATORFULLYIMPLICITCOMPRESSIBLEPOLYMER_HEADER_INCLUDED
 #define OPM_SIMULATORFULLYIMPLICITCOMPRESSIBLEPOLYMER_HEADER_INCLUDED
 
+#include <opm/core/simulator/SimulatorReport.hpp>
 #include <opm/core/utility/parameters/ParameterGroup.hpp>
 #include <opm/common/ErrorMacros.hpp>
 
@@ -125,7 +126,17 @@ namespace Opm
                                    const Wells* wells,
                                    const WellState& well_state,
                                    DynamicListEconLimited& list_econ_limited) const;
+
+        /// return the statistics if the nonlinearIteration() method failed.
+        ///
+        /// NOTE: for the flow_legacy simulator family this method is a stub, i.e. the
+        /// failure report object will *not* contain any meaningful data.
+        const SimulatorReport& failureReport() const
+        { return failureReport_; }
+
 private:
+        SimulatorReport failureReport_;
+
         const Deck& deck_;
         const PolymerPropsAd& polymer_props_;
 

opm/simulators/timestepping/AdaptiveTimeStepping.hpp

@@ -88,6 +88,11 @@ namespace Opm {
                               Output& outputWriter,
                               const std::vector<int>* fipnum = nullptr);
 
+        /** \brief Returns the simulator report for the failed substeps of the last
+         *         report step.
+         */
+        const SimulatorReport& failureReport() const { return failureReport_; };
+
         double suggestedNextStep() const { return suggested_next_timestep_; }
 
         void setSuggestedNextStep(const double x) { suggested_next_timestep_ = x; }
@@ -104,6 +109,7 @@ namespace Opm {
 
         typedef std::unique_ptr< TimeStepControlInterface > TimeStepControlType;
 
+        SimulatorReport failureReport_;       //!< statistics for the failed substeps of the last timestep
         TimeStepControlType timeStepControl_; //!< time step control object
         const double restart_factor_;         //!< factor to multiply time step with when solver fails to converge
         const double growth_factor_;          //!< factor to multiply time step when solver recovered from failed convergence

opm/simulators/timestepping/AdaptiveTimeStepping_impl.hpp

@@ -211,6 +211,9 @@ namespace Opm {
         State  last_state( state );
         WState last_well_state( well_state );
 
+        // reset the statistics for the failed substeps
+        failureReport_ = SimulatorReport();
+
         // counter for solver restarts
         int restarts = 0;
 
@@ -238,18 +241,26 @@ namespace Opm {
                 }
             }
             catch (const Opm::NumericalProblem& e) {
+                substepReport += solver.failureReport();
+
                 detail::logException(e, solver_verbose_);
                 // since linearIterations is < 0 this will restart the solver
             }
             catch (const std::runtime_error& e) {
+                substepReport += solver.failureReport();
+
                 detail::logException(e, solver_verbose_);
                 // also catch linear solver not converged
             }
             catch (const Dune::ISTLError& e) {
+                substepReport += solver.failureReport();
+
                 detail::logException(e, solver_verbose_);
                 // also catch errors in ISTL AMG that occur when time step is too large
             }
             catch (const Dune::MatrixBlockError& e) {
+                substepReport += solver.failureReport();
+
                 detail::logException(e, solver_verbose_);
                 // this can be thrown by ISTL's ILU0 in block mode, yet is not an ISTLError
             }
@@ -323,8 +334,10 @@ namespace Opm {
             }
             else // in case of no convergence (linearIterations < 0)
             {
-                report.converged = false;
                 substepTimer.setLastStepFailed(true);
+
+                failureReport_ += substepReport;
+
                 // increase restart counter
                 if( restarts >= solver_restart_max_ ) {
                     const auto msg = std::string("Solver failed to converge after cutting timestep ")