From 4deb00fd765325d7a22b7e3d1f1ec4c6a8f0a228 Mon Sep 17 00:00:00 2001
From: Dmitri Naumov <github@naumov.de>
Date: Mon, 14 Jan 2019 01:07:22 +0100
Subject: [PATCH] [NL] Solvers now return a NonlinearSolverStatus.

It contains a boolean for met error norms and a
number of iterations of the solver made.

The number of iterations is needed for
a timestepping scheme.
---
 NumLib/ODESolver/NonlinearSolver.cpp      |  8 +--
 NumLib/ODESolver/NonlinearSolver.h        | 10 ++-
 NumLib/ODESolver/NonlinearSolverStatus.h  | 20 ++++++
 ProcessLib/ProcessData.h                  |  6 +-
 ProcessLib/UncoupledProcessesTimeLoop.cpp | 77 +++++++++++------------
 ProcessLib/UncoupledProcessesTimeLoop.h   |  8 +--
 6 files changed, 76 insertions(+), 53 deletions(-)
 create mode 100644 NumLib/ODESolver/NonlinearSolverStatus.h

diff --git a/NumLib/ODESolver/NonlinearSolver.cpp b/NumLib/ODESolver/NonlinearSolver.cpp
index 34b68924b23..6b1642cae70 100644
--- a/NumLib/ODESolver/NonlinearSolver.cpp
+++ b/NumLib/ODESolver/NonlinearSolver.cpp
@@ -26,7 +26,7 @@ void NonlinearSolver<NonlinearSolverTag::Picard>::assemble(
     _equation_system->assemble(x);
 }
 
-bool NonlinearSolver<NonlinearSolverTag::Picard>::solve(
+NonlinearSolverStatus NonlinearSolver<NonlinearSolverTag::Picard>::solve(
     GlobalVector& x,
     std::function<void(int, GlobalVector const&)> const&
         postIterationCallback)
@@ -161,7 +161,7 @@ bool NonlinearSolver<NonlinearSolverTag::Picard>::solve(
     NumLib::GlobalVectorProvider::provider.releaseVector(rhs);
     NumLib::GlobalVectorProvider::provider.releaseVector(x_new);
 
-    return error_norms_met;
+    return {error_norms_met, iteration};
 }
 
 void NonlinearSolver<NonlinearSolverTag::Newton>::assemble(
@@ -173,7 +173,7 @@ void NonlinearSolver<NonlinearSolverTag::Newton>::assemble(
     //      equation every time and could not forget it.
 }
 
-bool NonlinearSolver<NonlinearSolverTag::Newton>::solve(
+NonlinearSolverStatus NonlinearSolver<NonlinearSolverTag::Newton>::solve(
     GlobalVector& x,
     std::function<void(int, GlobalVector const&)> const&
         postIterationCallback)
@@ -312,7 +312,7 @@ bool NonlinearSolver<NonlinearSolverTag::Newton>::solve(
     NumLib::GlobalVectorProvider::provider.releaseVector(
         minus_delta_x);
 
-    return error_norms_met;
+    return {error_norms_met, iteration};
 }
 
 std::pair<std::unique_ptr<NonlinearSolverBase>, NonlinearSolverTag>
diff --git a/NumLib/ODESolver/NonlinearSolver.h b/NumLib/ODESolver/NonlinearSolver.h
index d2eb306f08c..a77f7025318 100644
--- a/NumLib/ODESolver/NonlinearSolver.h
+++ b/NumLib/ODESolver/NonlinearSolver.h
@@ -14,6 +14,7 @@
 #include <logog/include/logog.hpp>
 
 #include "ConvergenceCriterion.h"
+#include "NonlinearSolverStatus.h"
 #include "NonlinearSystem.h"
 #include "Types.h"
 
@@ -53,7 +54,8 @@ public:
      * \retval true if the equation system could be solved
      * \retval false otherwise
      */
-    virtual bool solve(GlobalVector& x,
+    virtual NonlinearSolverStatus solve(
+        GlobalVector& x,
         std::function<void(int, GlobalVector const&)> const&
             postIterationCallback) = 0;
 
@@ -105,7 +107,8 @@ public:
 
     void assemble(GlobalVector const& x) const override;
 
-    bool solve(GlobalVector& x,
+    NonlinearSolverStatus solve(
+        GlobalVector& x,
         std::function<void(int, GlobalVector const&)> const&
             postIterationCallback) override;
 
@@ -164,7 +167,8 @@ public:
 
     void assemble(GlobalVector const& x) const override;
 
-    bool solve(GlobalVector& x,
+    NonlinearSolverStatus solve(
+        GlobalVector& x,
         std::function<void(int, GlobalVector const&)> const&
             postIterationCallback) override;
 
diff --git a/NumLib/ODESolver/NonlinearSolverStatus.h b/NumLib/ODESolver/NonlinearSolverStatus.h
new file mode 100644
index 00000000000..daddb29209d
--- /dev/null
+++ b/NumLib/ODESolver/NonlinearSolverStatus.h
@@ -0,0 +1,20 @@
+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2019, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ */
+
+#pragma once
+
+namespace NumLib
+{
+/// Status of the non-linear solver.
+struct NonlinearSolverStatus
+{
+    bool error_norms_met;
+    int number_iterations;
+};
+}  // namespace NumLib
diff --git a/ProcessLib/ProcessData.h b/ProcessLib/ProcessData.h
index 540e15964cf..e1d9321c616 100644
--- a/ProcessLib/ProcessData.h
+++ b/ProcessLib/ProcessData.h
@@ -34,7 +34,7 @@ struct ProcessData
         : timestepper(std::move(timestepper_)),
           nonlinear_solver_tag(NLTag),
           nonlinear_solver(nonlinear_solver),
-          nonlinear_solver_converged(true),
+          nonlinear_solver_status{true, 0},
           conv_crit(std::move(conv_crit_)),
           time_disc(std::move(time_disc_)),
           process(process_)
@@ -45,7 +45,7 @@ struct ProcessData
         : timestepper(std::move(pd.timestepper)),
           nonlinear_solver_tag(pd.nonlinear_solver_tag),
           nonlinear_solver(pd.nonlinear_solver),
-          nonlinear_solver_converged(pd.nonlinear_solver_converged),
+          nonlinear_solver_status(pd.nonlinear_solver_status),
           conv_crit(std::move(pd.conv_crit)),
           time_disc(std::move(pd.time_disc)),
           tdisc_ode_sys(std::move(pd.tdisc_ode_sys)),
@@ -66,7 +66,7 @@ struct ProcessData
     //! other members of this struct to their concrety types.
     NumLib::NonlinearSolverTag const nonlinear_solver_tag;
     NumLib::NonlinearSolverBase& nonlinear_solver;
-    bool nonlinear_solver_converged;
+    NumLib::NonlinearSolverStatus nonlinear_solver_status;
     std::unique_ptr<NumLib::ConvergenceCriterion> conv_crit;
 
     std::unique_ptr<NumLib::TimeDiscretization> time_disc;
diff --git a/ProcessLib/UncoupledProcessesTimeLoop.cpp b/ProcessLib/UncoupledProcessesTimeLoop.cpp
index 31ce7793784..7cab9ddaece 100644
--- a/ProcessLib/UncoupledProcessesTimeLoop.cpp
+++ b/ProcessLib/UncoupledProcessesTimeLoop.cpp
@@ -235,11 +235,10 @@ std::vector<GlobalVector*> setInitialConditions(
     return process_solutions;
 }
 
-bool solveOneTimeStepOneProcess(int const process_id, GlobalVector& x,
-                                std::size_t const timestep, double const t,
-                                double const delta_t,
-                                ProcessData const& process_data,
-                                Output& output_control)
+NumLib::NonlinearSolverStatus solveOneTimeStepOneProcess(
+    int const process_id, GlobalVector& x, std::size_t const timestep,
+    double const t, double const delta_t, ProcessData const& process_data,
+    Output& output_control)
 {
     auto& process = process_data.process;
     auto& time_disc = *process_data.time_disc;
@@ -263,15 +262,15 @@ bool solveOneTimeStepOneProcess(int const process_id, GlobalVector& x,
                                                   timestep, t, x, iteration);
     };
 
-    bool nonlinear_solver_succeeded =
+    auto const nonlinear_solver_status =
         nonlinear_solver.solve(x, post_iteration_callback);
 
-    if (nonlinear_solver_succeeded)
+    if (nonlinear_solver_status.error_norms_met)
     {
         process.postNonLinearSolver(x, t, process_id);
     }
 
-    return nonlinear_solver_succeeded;
+    return nonlinear_solver_status;
 }
 
 UncoupledProcessesTimeLoop::UncoupledProcessesTimeLoop(
@@ -352,7 +351,7 @@ double UncoupledProcessesTimeLoop::computeTimeStepping(
                              x, norm_type))
                 : 0.;
 
-        if (!ppd.nonlinear_solver_converged)
+        if (!ppd.nonlinear_solver_status.error_norms_met)
         {
             timestepper->setAcceptedOrNot(false);
         }
@@ -366,7 +365,7 @@ double UncoupledProcessesTimeLoop::computeTimeStepping(
             all_process_steps_accepted = false;
         }
 
-        if (!ppd.nonlinear_solver_converged)
+        if (!ppd.nonlinear_solver_status.error_norms_met)
         {
             WARN("Time step will be rejected due to nonlinear solver diverged");
             all_process_steps_accepted = false;
@@ -525,7 +524,7 @@ bool UncoupledProcessesTimeLoop::loop()
     double t = _start_time;
     std::size_t accepted_steps = 0;
     std::size_t rejected_steps = 0;
-    bool nonlinear_solver_succeeded = true;
+    NumLib::NonlinearSolverStatus nonlinear_solver_status{true, 0};
 
     double dt = computeTimeStepping(0.0, t, accepted_steps, rejected_steps);
 
@@ -552,12 +551,12 @@ bool UncoupledProcessesTimeLoop::loop()
 
         if (is_staggered_coupling)
         {
-            nonlinear_solver_succeeded =
+            nonlinear_solver_status =
                 solveCoupledEquationSystemsByStaggeredScheme(t, dt, timesteps);
         }
         else
         {
-            nonlinear_solver_succeeded =
+            nonlinear_solver_status =
                 solveUncoupledEquationSystems(t, dt, timesteps);
         }
 
@@ -616,7 +615,7 @@ bool UncoupledProcessesTimeLoop::loop()
         accepted_steps + rejected_steps, accepted_steps, rejected_steps);
 
     // output last time step
-    if (nonlinear_solver_succeeded)
+    if (nonlinear_solver_status.error_norms_met)
     {
         const bool output_initial_condition = false;
         outputSolutions(output_initial_condition, is_staggered_coupling,
@@ -624,17 +623,18 @@ bool UncoupledProcessesTimeLoop::loop()
                         &Output::doOutputLastTimestep);
     }
 
-    return nonlinear_solver_succeeded;
+    return nonlinear_solver_status.error_norms_met;
 }
 
 static std::string const nonlinear_fixed_dt_fails_info =
-    "Nonlinear solver fails. Because the time stepper"
-    " FixedTimeStepping is used, the program has to be"
-    " terminated ";
+    "Nonlinear solver fails. Because the time stepper FixedTimeStepping is "
+    "used, the program has to be terminated.";
 
-bool UncoupledProcessesTimeLoop::solveUncoupledEquationSystems(
+NumLib::NonlinearSolverStatus
+UncoupledProcessesTimeLoop::solveUncoupledEquationSystems(
     const double t, const double dt, const std::size_t timestep_id)
 {
+    NumLib::NonlinearSolverStatus nonlinear_solver_status;
     // TODO(wenqing): use process name
     unsigned process_id = 0;
     for (auto& process_data : _per_process_data)
@@ -653,19 +653,19 @@ bool UncoupledProcessesTimeLoop::solveUncoupledEquationSystems(
         auto& pcs = process_data->process;
         pcs.preTimestep(x, t, dt, process_id);
 
-        const auto nonlinear_solver_succeeded = solveOneTimeStepOneProcess(
+        nonlinear_solver_status = solveOneTimeStepOneProcess(
             process_id, x, timestep_id, t, dt, *process_data, *_output);
-        process_data->nonlinear_solver_converged = nonlinear_solver_succeeded;
+        process_data->nonlinear_solver_status = nonlinear_solver_status;
         pcs.postTimestep(x, t, dt, process_id);
         pcs.computeSecondaryVariable(t, x, process_id);
 
         INFO("[time] Solving process #%u took %g s in time step #%u ",
              process_id, time_timestep_process.elapsed(), timestep_id);
 
-        if (!nonlinear_solver_succeeded)
+        if (!nonlinear_solver_status.error_norms_met)
         {
-            ERR("The nonlinear solver failed in time step #%u at t = %g "
-                "s for process #%u.",
+            ERR("The nonlinear solver failed in time step #%u at t = %g s for "
+                "process #%u.",
                 timestep_id, t, process_id);
 
             if (!process_data->timestepper->isSolutionErrorComputationNeeded())
@@ -675,16 +675,17 @@ bool UncoupledProcessesTimeLoop::solveUncoupledEquationSystems(
                 OGS_FATAL(nonlinear_fixed_dt_fails_info.data());
             }
 
-            return false;
+            return nonlinear_solver_status;
         }
 
         ++process_id;
     }  // end of for (auto& process_data : _per_process_data)
 
-    return true;
+    return nonlinear_solver_status;
 }
 
-bool UncoupledProcessesTimeLoop::solveCoupledEquationSystemsByStaggeredScheme(
+NumLib::NonlinearSolverStatus
+UncoupledProcessesTimeLoop::solveCoupledEquationSystemsByStaggeredScheme(
     const double t, const double dt, const std::size_t timestep_id)
 {
     // Coupling iteration
@@ -703,13 +704,13 @@ bool UncoupledProcessesTimeLoop::solveCoupledEquationSystemsByStaggeredScheme(
         }
     };
 
+    NumLib::NonlinearSolverStatus nonlinear_solver_status{true, 0};
     bool coupling_iteration_converged = true;
     for (int global_coupling_iteration = 0;
          global_coupling_iteration < _global_coupling_max_iterations;
          global_coupling_iteration++, resetCouplingConvergenceCriteria())
     {
         // TODO(wenqing): use process name
-        bool nonlinear_solver_succeeded = true;
         coupling_iteration_converged = true;
         int process_id = 0;
         int const last_process_id = _per_process_data.size() - 1;
@@ -741,10 +742,9 @@ bool UncoupledProcessesTimeLoop::solveCoupledEquationSystemsByStaggeredScheme(
             process_data->process.setCoupledSolutionsForStaggeredScheme(
                 &coupled_solutions);
 
-            const auto nonlinear_solver_succeeded = solveOneTimeStepOneProcess(
+            nonlinear_solver_status = solveOneTimeStepOneProcess(
                 process_id, x, timestep_id, t, dt, *process_data, *_output);
-            process_data->nonlinear_solver_converged =
-                nonlinear_solver_succeeded;
+            process_data->nonlinear_solver_status = nonlinear_solver_status;
 
             INFO(
                 "[time] Solving process #%u took %g s in time step #%u "
@@ -752,11 +752,10 @@ bool UncoupledProcessesTimeLoop::solveCoupledEquationSystemsByStaggeredScheme(
                 process_id, time_timestep_process.elapsed(), timestep_id,
                 global_coupling_iteration);
 
-            if (!nonlinear_solver_succeeded)
+            if (!nonlinear_solver_status.error_norms_met)
             {
-                ERR("The nonlinear solver failed in time step #%u at t = %g "
-                    "s"
-                    " for process #%u.",
+                ERR("The nonlinear solver failed in time step #%u at t = %g s "
+                    "for process #%u.",
                     timestep_id, t, process_id);
 
                 if (!process_data->timestepper
@@ -797,9 +796,9 @@ bool UncoupledProcessesTimeLoop::solveCoupledEquationSystemsByStaggeredScheme(
             break;
         }
 
-        if (!nonlinear_solver_succeeded)
+        if (!nonlinear_solver_status.error_norms_met)
         {
-            return false;
+            return nonlinear_solver_status;
         }
     }
 
@@ -833,7 +832,7 @@ bool UncoupledProcessesTimeLoop::solveCoupledEquationSystemsByStaggeredScheme(
         ++process_id;
     }
 
-    return true;
+    return nonlinear_solver_status;
 }
 
 template <typename OutputClass, typename OutputClassMember>
@@ -848,7 +847,7 @@ void UncoupledProcessesTimeLoop::outputSolutions(
         auto& pcs = process_data->process;
         // If nonlinear solver diverged, the solution has already been
         // saved.
-        if ((!process_data->nonlinear_solver_converged) ||
+        if ((!process_data->nonlinear_solver_status.error_norms_met) ||
             process_data->skip_time_stepping)
         {
             ++process_id;
diff --git a/ProcessLib/UncoupledProcessesTimeLoop.h b/ProcessLib/UncoupledProcessesTimeLoop.h
index 9620a315ec9..92a1351cb35 100644
--- a/ProcessLib/UncoupledProcessesTimeLoop.h
+++ b/ProcessLib/UncoupledProcessesTimeLoop.h
@@ -35,7 +35,7 @@ struct ProcessData;
 class UncoupledProcessesTimeLoop
 {
 public:
-    explicit UncoupledProcessesTimeLoop(
+    UncoupledProcessesTimeLoop(
         std::unique_ptr<Output>&& output,
         std::vector<std::unique_ptr<ProcessData>>&& per_process_data,
         const int global_coupling_max_iterations,
@@ -98,8 +98,8 @@ private:
      * @return            true:  if all nonlinear solvers convergence.
      *                    false: if any of nonlinear solvers divergences.
      */
-    bool solveUncoupledEquationSystems(const double t, const double dt,
-                                       const std::size_t timestep_id);
+    NumLib::NonlinearSolverStatus solveUncoupledEquationSystems(
+        const double t, const double dt, const std::size_t timestep_id);
 
     /**
      * \brief Member to solver coupled systems of equations by the staggered
@@ -111,7 +111,7 @@ private:
      * @return            true:   if all nonlinear solvers convergence.
      *                    false:  if any of nonlinear solvers divergences.
      */
-    bool solveCoupledEquationSystemsByStaggeredScheme(
+    NumLib::NonlinearSolverStatus solveCoupledEquationSystemsByStaggeredScheme(
         const double t, const double dt, const std::size_t timestep_id);
 
     /**
-- 
GitLab