]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Convert NOX to use the new exception scheme.
authorWolfgang Bangerth <bangerth@colostate.edu>
Fri, 5 May 2023 21:06:37 +0000 (15:06 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Tue, 9 May 2023 18:31:36 +0000 (12:31 -0600)
include/deal.II/trilinos/nox.h
include/deal.II/trilinos/nox.templates.h

index 7c1efc0ba4bacd579cc88ba512487aa49a614ecc..3ca2cc26126290e062fc2e7ff88301fcc0bf9b51 100644 (file)
@@ -26,6 +26,7 @@
 
 #  include <Teuchos_ParameterList.hpp>
 
+#  include <exception>
 #  include <functional>
 
 DEAL_II_NAMESPACE_OPEN
@@ -180,17 +181,29 @@ namespace TrilinosWrappers
      * A function object that users should supply and that is intended to
      * compute the residual $F(u)$.
      *
-     * @note This function should return 0 in the case of success.
+     * @note This variable represents a
+     * @ref GlossUserProvidedCallBack "user provided callback".
+     * See there for a description of how to deal with errors and other
+     * requirements and conventions. NOX can not deal
+     * with "recoverable" errors, so if a callback
+     * throws an exception of type RecoverableUserCallbackError, then this
+     * exception is treated like any other exception.
      */
-    std::function<int(const VectorType &u, VectorType &F)> residual;
+    std::function<void(const VectorType &u, VectorType &F)> residual;
 
     /**
      * A user function that sets up the Jacobian, based on the
      * current solution @p current_u.
      *
-     * @note This function should return 0 in the case of success.
+     * @note This variable represents a
+     * @ref GlossUserProvidedCallBack "user provided callback".
+     * See there for a description of how to deal with errors and other
+     * requirements and conventions. NOX can not deal
+     * with "recoverable" errors, so if a callback
+     * throws an exception of type RecoverableUserCallbackError, then this
+     * exception is treated like any other exception.
      */
-    std::function<int(const VectorType &current_u)> setup_jacobian;
+    std::function<void(const VectorType &current_u)> setup_jacobian;
 
     /**
      * A user function that sets up the preconditioner for inverting
@@ -201,9 +214,15 @@ namespace TrilinosWrappers
      * update_preconditioner_predicate and
      * AdditionalData::threshold_nonlinear_iterations).
      *
-     * @note This function should return 0 in the case of success.
+     * @note This variable represents a
+     * @ref GlossUserProvidedCallBack "user provided callback".
+     * See there for a description of how to deal with errors and other
+     * requirements and conventions. NOX can not deal
+     * with "recoverable" errors, so if a callback
+     * throws an exception of type RecoverableUserCallbackError, then this
+     * exception is treated like any other exception.
      */
-    std::function<int(const VectorType &current_u)> setup_preconditioner;
+    std::function<void(const VectorType &current_u)> setup_preconditioner;
 
     /**
      * A user function that applies the Jacobian $\nabla_u F(u)$ to
@@ -217,9 +236,15 @@ namespace TrilinosWrappers
      * chosen, whereas for the full step case (`NOX::LineSearch::FullStep`)
      * it won't be called.
      *
-     * @note This function should return 0 in the case of success.
+     * @note This variable represents a
+     * @ref GlossUserProvidedCallBack "user provided callback".
+     * See there for a description of how to deal with errors and other
+     * requirements and conventions. NOX can not deal
+     * with "recoverable" errors, so if a callback
+     * throws an exception of type RecoverableUserCallbackError, then this
+     * exception is treated like any other exception.
      */
-    std::function<int(const VectorType &x, VectorType &y)> apply_jacobian;
+    std::function<void(const VectorType &x, VectorType &y)> apply_jacobian;
 
     /**
      * A user function that applies the inverse of the Jacobian
@@ -233,10 +258,16 @@ namespace TrilinosWrappers
      * @note This function is optional and is used in the case of certain
      * configurations.
      *
-     * @note This function should return 0 in the case of success.
+     * @note This variable represents a
+     * @ref GlossUserProvidedCallBack "user provided callback".
+     * See there for a description of how to deal with errors and other
+     * requirements and conventions. NOX can not deal
+     * with "recoverable" errors, so if a callback
+     * throws an exception of type RecoverableUserCallbackError, then this
+     * exception is treated like any other exception.
      */
     std::function<
-      int(const VectorType &y, VectorType &x, const double tolerance)>
+      void(const VectorType &y, VectorType &x, const double tolerance)>
       solve_with_jacobian;
 
     /**
@@ -251,6 +282,14 @@ namespace TrilinosWrappers
      *
      * @note This function is optional and is used in the case of certain
      * configurations.
+     *
+     * @note This variable represents a
+     * @ref GlossUserProvidedCallBack "user provided callback".
+     * See there for a description of how to deal with errors and other
+     * requirements and conventions. NOX can not deal
+     * with "recoverable" errors, so if a callback
+     * throws an exception of type RecoverableUserCallbackError, then this
+     * exception is treated like any other exception.
      */
     std::function<
       int(const VectorType &y, VectorType &x, const double tolerance)>
@@ -266,6 +305,14 @@ namespace TrilinosWrappers
      * and the current residual vector @p f.
      *
      * @note This function is optional.
+     *
+     * @note This variable represents a
+     * @ref GlossUserProvidedCallBack "user provided callback".
+     * See there for a description of how to deal with errors and other
+     * requirements and conventions. NOX can not deal
+     * with "recoverable" errors, so if a callback
+     * throws an exception of type RecoverableUserCallbackError, then this
+     * exception is treated like any other exception.
      */
     std::function<SolverControl::State(const unsigned int i,
                                        const double       norm_f,
@@ -282,6 +329,14 @@ namespace TrilinosWrappers
      *
      * @note This function is optional. If no function is attached, this
      * means implicitly a return value of `false`.
+     *
+     * @note This variable represents a
+     * @ref GlossUserProvidedCallBack "user provided callback".
+     * See there for a description of how to deal with errors and other
+     * requirements and conventions. NOX can not deal
+     * with "recoverable" errors, so if a callback
+     * throws an exception of type RecoverableUserCallbackError, then this
+     * exception is treated like any other exception.
      */
     std::function<bool()> update_preconditioner_predicate;
 
@@ -317,6 +372,13 @@ namespace TrilinosWrappers
      * The number of linear iterations of the last Jacobian solve.
      */
     unsigned int n_last_linear_iterations;
+
+    /**
+     * A pointer to any exception that may have been thrown in user-defined
+     * call-backs and that we have to deal after the KINSOL function we call
+     * has returned.
+     */
+    mutable std::exception_ptr pending_exception;
   };
 } // namespace TrilinosWrappers
 
index b32f64dd12257bef42185f62be3f24e56d901b6c..69d3bdbb4e35593289f5cb9739b0fbacc3122f4d 100644 (file)
@@ -43,6 +43,46 @@ namespace TrilinosWrappers
   {
     namespace NOXWrappers
     {
+      namespace
+      {
+        /**
+         * A function that calls the function object given by its first argument
+         * with the set of arguments following at the end. If the call returns
+         * regularly, the current function returns zero to indicate success. If
+         * the call fails with an, then the current function returns with an
+         * error code of -1. In that case, the exception thrown by `f` is
+         * captured and `eptr` is set to the exception. In case of success,
+         * `eptr` is set to `nullptr`.
+         */
+        template <typename F, typename... Args>
+        int
+        call_and_possibly_capture_exception(const F &           f,
+                                            std::exception_ptr &eptr,
+                                            Args &&...args)
+        {
+          // See whether there is already something in the exception pointer
+          // variable. There is no reason why this should be so, and
+          // we should probably bail out:
+          AssertThrow(eptr == nullptr, ExcInternalError());
+
+          // Call the function and if that succeeds, return zero:
+          try
+            {
+              f(std::forward<Args>(args)...);
+              eptr = nullptr;
+              return 0;
+            }
+          // In case of an exception, capture the exception and
+          // return -1:
+          catch (...)
+            {
+              eptr = std::current_exception();
+              return -1;
+            }
+        }
+      } // namespace
+
+
       template <typename VectorType>
       class Group;
 
@@ -937,7 +977,8 @@ namespace TrilinosWrappers
         n_residual_evaluations++;
 
         // evalute residual
-        return residual(x, f);
+        return internal::NOXWrappers::call_and_possibly_capture_exception(
+          residual, pending_exception, x, f);
       },
       [&](const VectorType &x) -> int {
         Assert(
@@ -946,7 +987,8 @@ namespace TrilinosWrappers
             "No setup_jacobian function has been attached to the NOXSolver object."));
 
         // setup Jacobian
-        int flag = setup_jacobian(x);
+        int flag = internal::NOXWrappers::call_and_possibly_capture_exception(
+          setup_jacobian, pending_exception, x);
 
         if (flag != 0)
           return flag;
@@ -967,7 +1009,8 @@ namespace TrilinosWrappers
               update_preconditioner = update_preconditioner_predicate();
 
             if (update_preconditioner)
-              flag = setup_preconditioner(x);
+              flag = internal::NOXWrappers::call_and_possibly_capture_exception(
+                setup_preconditioner, pending_exception, x);
           }
 
         return flag;
@@ -981,7 +1024,8 @@ namespace TrilinosWrappers
         n_jacobian_applications++;
 
         // apply Jacobian
-        return apply_jacobian(x, v);
+        return internal::NOXWrappers::call_and_possibly_capture_exception(
+          apply_jacobian, pending_exception, x, v);
       },
       [&](const VectorType &f, VectorType &x, const double tolerance) -> int {
         n_nonlinear_iterations++;
@@ -996,7 +1040,8 @@ namespace TrilinosWrappers
                 "solve_with_jacobian_and_track_n_linear_iterations!"));
 
             // without tracking of linear iterations
-            return solve_with_jacobian(f, x, tolerance);
+            return internal::NOXWrappers::call_and_possibly_capture_exception(
+              solve_with_jacobian, pending_exception, f, x, tolerance);
           }
         else if (solve_with_jacobian_and_track_n_linear_iterations)
           {
@@ -1063,8 +1108,18 @@ namespace TrilinosWrappers
     // create non-linear solver
     const auto solver = NOX::Solver::buildSolver(group, check, parameters);
 
-    // solve
+    // Solve, then check whether an exception was thrown by one of the user
+    // callback functions. If so, exit by rethrowing the exception that
+    // we had previously saved. This also calls the destructors of all of
+    // the member variables above, so we do not have to clean things up by hand.
     const auto status = solver->solve();
+    if (pending_exception)
+      {
+        std::exception_ptr this_exception = pending_exception;
+        pending_exception                 = nullptr;
+
+        std::rethrow_exception(this_exception);
+      }
 
     AssertThrow(status == NOX::StatusTest::Converged, ExcNOXNoConvergence());
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.