]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add SNES, address comments, various improvements 15160/head
authorTimo Heister <timo.heister@gmail.com>
Mon, 15 May 2023 16:18:51 +0000 (12:18 -0400)
committerTimo Heister <timo.heister@gmail.com>
Sun, 11 Jun 2023 19:41:05 +0000 (15:41 -0400)
include/deal.II/numerics/nonlinear.h
tests/numerics/nonlinear_solver_selector_01.cc
tests/numerics/nonlinear_solver_selector_01_nox.cc [moved from tests/numerics/nonlinear_solver_selector_02.cc with 100% similarity]
tests/numerics/nonlinear_solver_selector_01_nox.with_trilinos=on.output [moved from tests/numerics/nonlinear_solver_selector_02.with_trilinos=on.output with 100% similarity]
tests/numerics/nonlinear_solver_selector_03.cc
tests/numerics/nonlinear_solver_selector_03_nox.cc [moved from tests/numerics/nonlinear_solver_selector_04.cc with 93% similarity]
tests/numerics/nonlinear_solver_selector_03_nox.with_trilinos=on.mpirun=4.output [moved from tests/numerics/nonlinear_solver_selector_04.with_trilinos=on.mpirun=4.output with 100% similarity]
tests/numerics/nonlinear_solver_selector_03_petsc.cc [new file with mode: 0644]
tests/numerics/nonlinear_solver_selector_03_petsc.with_petsc=on.mpirun=4.output [new file with mode: 0644]

index 11bb956a8c3a10d3734b592e16f3e598eb42cbb2..77d007599dc67ab4a8fdd5c37aee194707f140bc 100644 (file)
@@ -17,6 +17,8 @@
 
 #include <deal.II/base/exceptions.h>
 
+#include <deal.II/lac/petsc_snes.h>
+
 #include <deal.II/sundials/kinsol.h>
 
 #include <deal.II/trilinos/nox.h>
 DEAL_II_NAMESPACE_OPEN
 
 /**
- * Selects a nonlinear solver by choosing between KINSOL or NOX.
- * KINSOL and NOX are nonlinear solvers included in the SUNDIALS
- * package and the Trilinos package, respectively. If no solver type is
- * specified it will automaticlaly choose a solver based on what
- * deal.II was configured with, KINSOL having priority.
+ * This class offers a unified interface to several nonlinear solver
+ * implementations like KINSOL, NOX, and SNES.
+ *
+ * KINSOL nonlinear solvers are part of the SUNDIALS package, while
+ * NOX is part of Trilinos and SNES is part of PETSc, respectively.
+ * If no solver is manually specified, this class will automaticlaly
+ * choose one of the available solvers depending on the enabled
+ * dependencies.
  *
- * By calling the @p solve function of this @p NonlinearSolverSelector, it selects the
- * @p solve function of that @p Solver that was specified in the constructor
- * of this class, similar to the SolverSelector class.
+ * By calling the @p solve function of this @p
+ * NonlinearSolverSelector, it selects the @p solve function of that
+ * @p Solver that was specified in the constructor of this class,
+ * similar to the SolverSelector class.
  *
  * <h3>Usage</h3>
  * An example of code one would run with this class is the
@@ -48,20 +54,19 @@ DEAL_II_NAMESPACE_OPEN
  *
  * // Functions that are required for solving a nonlinear problem,
  * // which are utilized in both @p KINSOL and @p NOX.
- * nonlinear_solver.reinit_vector = [&](Vector<double> &x) {...}
+ * nonlinear_solver.reinit_vector = [&](Vector<double> &x) {...};
  *
  * nonlinear_solver.residual =
  *           [&](const Vector<double> &evaluation_point,
- *               Vector<double> &      residual) {...}
+ *               Vector<double> &      residual) {...};
  *
  * nonlinear_solver.setup_jacobian =
- *           [&](const Vector<double> &current_u,
- *               const Vector<double> &current_f) {...}
+ *           [&](const Vector<double> &current_u) {...};
  *
  * nonlinear_solver.solve_with_jacobian =
  *           [&](const Vector<double> &rhs,
  *               Vector<double> &      dst,
- *               const double tolerance) {...}
+ *               const double tolerance) {...};
  *
  * // Calling the @p solve function with an initial guess.
  * nonlinear_solver.solve(current_solution);
@@ -92,7 +97,7 @@ public:
        */
       newton,
       /**
-       * Newton iteration with linesearch.
+       * Newton iteration with line search.
        */
       linesearch,
       /**
@@ -116,20 +121,25 @@ public:
        * NOX nonlinear solver, part of the TRILINOS package.
        */
       nox,
+      /*
+       * Use the PETSc SNES solver.
+       */
+      petsc_snes
     };
 
     /**
      * Initialization parameters for NonlinearSolverSelector.
      *
-     * @param solver_type Nonlinear solver type, can be 'auto', 'kinsol', or 'nox'.
-     * @param strategy Method of solving nonlinear problem, can be 'newton',
-     * 'linesearch', or 'picard'.
+     * @param solver_type Nonlinear solver type.
+     * @param strategy Method of solving the nonlinear problem.
      * @param maximum_non_linear_iterations Maximum number of nonlinear
-     * iterations. This parameter is shared between KINSOL and NOX.
-     * @param function_tolerance Function norm stopping tolerance.
-     * @param relative_tolerance Relative function norm stopping tolerance.
-     * @param step_tolerance Step tolerance for solution update.
-     * @param anderson_subspace_size Anderson acceleration subspace size
+     *   iterations.
+     * @param function_tolerance Absolute stopping tolerance for the norm
+     *   of the residual $F(u)$.
+     * @param relative_tolerance Relative stopping tolerance.
+     * @param step_tolerance Tolerance for minimum scaled step length
+     * @param anderson_subspace_size Size of the Anderson acceleration
+     *   subspace, use 0 to disable.
      */
     AdditionalData(const SolverType &      solver_type = automatic,
                    const SolutionStrategy &strategy    = linesearch,
@@ -140,9 +150,9 @@ public:
                    const unsigned int      anderson_subspace_size        = 0);
 
     /**
-     * The type of nonlinear solver to use. The default value is set to 'auto',
-     * which will use either KINSOL or NOX depending on which package is
-     * installed, KINSOL having priority.
+     * The type of nonlinear solver to use. If the default 'automatic' is used,
+     * it will choose the first available package in the order KINSOL, NOX, or
+     * SNES.
      */
     SolverType solver_type;
 
@@ -163,29 +173,29 @@ public:
      * A scalar used as a stopping tolerance on the scaled
      * maximum norm of the system function $F(u)$ or $G(u)$.
      *
-     * If set to zero, default values provided by KINSOL will be used.
+     * If set to zero, default values will be used.
      */
     double function_tolerance;
 
     /**
-     * A scalar used as a stopping tolerance on the minimum
-     * scaled step length.
+     * Relative $l_2$ tolerance of the residual to be reached.
      *
-     * If set to zero, default values provided by KINSOL will be used.
+     * @note Solver terminates successfully if either the function
+     * tolerance or the relative tolerance has been reached.
      */
-    double step_tolerance = 0.0;
+    const double relative_tolerance;
 
     /**
-     * Relative $l_2$ tolerance of the residual to be reached.
+     * A scalar used as a stopping tolerance on the minimum
+     * scaled step length.
      *
-     * @note Solver terminates successfully if either the absolute or
-     * the relative tolerance has been reached.
+     * If set to zero, default values will be used.
      */
-    const double relative_tolerance;
+    double step_tolerance;
 
     /**
      * The size of the subspace used with Anderson acceleration
-     * in conjunction with Picard or fixed-point iteration.
+     * in conjunction with Picard iteration.
      *
      * If you set this to 0, no acceleration is used.
      */
@@ -204,11 +214,10 @@ public:
   NonlinearSolverSelector(const AdditionalData &additional_data);
 
   /**
-   * Constructor
+   * Constructor.
    *
    * @param additional_data NonlinearSolverSelector configuration data
-   * @param mpi_communicator MPI communicator over which logging operations are
-   * computer.
+   * @param mpi_communicator MPI communicator used by the nonlinear solver.
    */
   NonlinearSolverSelector(const AdditionalData &additional_data,
                           const MPI_Comm &      mpi_communicator);
@@ -220,8 +229,15 @@ public:
   void
   select(const typename AdditionalData::SolverType &type);
 
+  /**
+   * Set the generic additional data for all nonlinear solvers.
+   */
+  void
+  set_data(const AdditionalData &additional_data);
+
 /**
- * Set the additional data. For more information see the @p Solver class.
+ * Set the additional data for NOX. See TrilinosWrappers::NOXSolver
+ * for more information.
  */
 #ifdef DEAL_II_WITH_TRILINOS
   void
@@ -233,7 +249,8 @@ public:
 #endif
 
 /**
- * Set the additional data. For more information see the @p Solver class.
+ * Set the additional data for KINSOL. See SUNDIALS::KINSOL for more
+ * information.
  */
 #ifdef DEAL_II_WITH_SUNDIALS
   void
@@ -241,13 +258,20 @@ public:
              &additional_data);
 #endif
 
+/**
+ * Set the additional data for SNES. See PETScWrappers::NonlinearSolverData
+ * for more information.
+ */
+#ifdef DEAL_II_WITH_PETSC
+  void
+  set_data(const typename PETScWrappers::NonlinearSolverData &additional_data);
+#endif
+
   /**
-   * Solve the nonlinear system. KINSOL uses the content of
-   * `initial_guess_and_solution` as an initial guess, and
-   * stores the final solution in the same vector.
+   * Solve the nonlinear system.
    *
-   * The functions herein are nearly identical in setup to what can be found
-   * in the KINSOL and NOX documentation.
+   * The content of `initial_guess_and_solution` is used as an initial guess
+   * and the final solution is stored in the same vector.
    */
   void
   solve(VectorType &initial_guess_and_solution);
@@ -283,14 +307,6 @@ public:
    * F/\partial u$. If strategy = SolutionStrategy::picard, $A$ is the
    * approximate Jacobian matrix $L$.
    *
-   * The setup_jacobian() function may call a user-supplied function, or a
-   * function within the linear solver module, to compute Jacobian-related
-   * data that is required by the linear solver. It may also preprocess that
-   * data as needed for solve_with_jacobian(), which may involve calling a
-   * generic function (such as for LU factorization) or, more generally,
-   * build preconditioners from the assembled Jacobian. In any case, the
-   * data so generated may then be used whenever a linear system is solved.
-   *
    * @param current_u Current value of $u$
    */
   std::function<int(const VectorType &current_u)> setup_jacobian;
@@ -299,12 +315,11 @@ public:
    * A function object that users may supply and that is intended to solve
    * a linear system with the Jacobian matrix.
    *
-   * Specific details on this function can be found in the KINSOL
-   *
-   * Arguments to the function are:
+   * @note PETSc SNES does not provide us with a linear tolerance to solve
+   *       linear system with. We will provide a value of 1e-6 in that case.
    *
    * @param[in] rhs The system right hand side to solve for.
-   * @param[out] dst The solution of $J^{-1} * src$.
+   * @param[out] dst The solution of $J^{-1} * \texttt{src}$.
    * @param[in] tolerance The tolerance with which to solve the linear system
    *   of equations.
    */
@@ -312,7 +327,7 @@ public:
     int(const VectorType &rhs, VectorType &dst, const double tolerance)>
     solve_with_jacobian;
 
-protected:
+private:
   /**
    * NonlinearSolverSelector configuration data.
    */
@@ -341,16 +356,28 @@ private:
     Teuchos::rcp(new Teuchos::ParameterList);
 #endif
 
+/**
+ * PETSc SNES configuration data
+ */
+#ifdef DEAL_II_WITH_PETSC
+  typename PETScWrappers::NonlinearSolverData additional_data_petsc_snes;
+#endif
+
   /**
-   * Data transfer function
+   * Solve with PETSc SNES. Internal functions specialized for PETSc
+   * Vectors. This is necessary to only instantiate the SNES solvers
+   * with PETSc Vectors, as this is the only vector type currently
+   * supported.
    */
   void
-  data_transfer(const AdditionalData &additional_data);
+  solve_with_petsc(VectorType &initial_guess_and_solution);
 };
 
+
+
 template <typename VectorType>
 void
-NonlinearSolverSelector<VectorType>::data_transfer(
+NonlinearSolverSelector<VectorType>::set_data(
   const AdditionalData &additional_data)
 {
 #ifdef DEAL_II_WITH_SUNDIALS
@@ -381,7 +408,6 @@ NonlinearSolverSelector<VectorType>::data_transfer(
 
 // Do the same thing we did above but with NOX
 #ifdef DEAL_II_WITH_TRILINOS
-  // Some default settings for parameters.
   parameters_nox->set("Nonlinear Solver", "Line Search Based");
   Teuchos::ParameterList &Line_Search = parameters_nox->sublist("Line Search");
   Line_Search.set("Method", "Full Step");
@@ -390,19 +416,62 @@ NonlinearSolverSelector<VectorType>::data_transfer(
   additional_data_nox.abs_tol  = additional_data.function_tolerance;
   additional_data_nox.rel_tol  = additional_data.relative_tolerance;
 #endif
+
+#ifdef DEAL_II_WITH_PETSC
+  additional_data_petsc_snes.options_prefix = "";
+
+  if (additional_data.anderson_subspace_size > 0 &&
+      additional_data.strategy ==
+        NonlinearSolverSelector<VectorType>::AdditionalData::picard)
+    {
+      additional_data_petsc_snes.snes_type = "anderson";
+      // TODO additional_data.anderson_subspace_size;
+    }
+  else if (additional_data.strategy ==
+           NonlinearSolverSelector<VectorType>::AdditionalData::linesearch)
+    {
+      additional_data_petsc_snes.snes_type            = "newtonls";
+      additional_data_petsc_snes.snes_linesearch_type = "bt";
+    }
+  else if (additional_data.strategy ==
+           NonlinearSolverSelector<VectorType>::AdditionalData::newton)
+    {
+      additional_data_petsc_snes.snes_linesearch_type = "newtonls";
+      additional_data_petsc_snes.snes_linesearch_type = "basic";
+    }
+  else if (additional_data.strategy ==
+           NonlinearSolverSelector<VectorType>::AdditionalData::picard)
+    additional_data_petsc_snes.snes_type = "nrichardson";
+
+  additional_data_petsc_snes.absolute_tolerance =
+    additional_data.function_tolerance;
+  additional_data_petsc_snes.relative_tolerance =
+    additional_data.relative_tolerance;
+  additional_data_petsc_snes.step_tolerance = additional_data.step_tolerance;
+  additional_data_petsc_snes.maximum_non_linear_iterations =
+    additional_data.maximum_non_linear_iterations;
+  additional_data_petsc_snes.max_n_function_evaluations = -1;
+
+#endif
 }
 
+
+
 template <typename VectorType>
 NonlinearSolverSelector<VectorType>::NonlinearSolverSelector() = default;
 
+
+
 template <typename VectorType>
 NonlinearSolverSelector<VectorType>::NonlinearSolverSelector(
   const AdditionalData &additional_data)
   : additional_data(additional_data)
 {
-  data_transfer(additional_data);
+  set_data(additional_data);
 }
 
+
+
 template <typename VectorType>
 NonlinearSolverSelector<VectorType>::NonlinearSolverSelector(
   const AdditionalData &additional_data,
@@ -410,10 +479,11 @@ NonlinearSolverSelector<VectorType>::NonlinearSolverSelector(
   : additional_data(additional_data)
   , mpi_communicator(mpi_communicator)
 {
-  data_transfer(additional_data);
+  set_data(additional_data);
 }
 
 
+
 template <typename VectorType>
 void
 NonlinearSolverSelector<VectorType>::select(
@@ -422,6 +492,8 @@ NonlinearSolverSelector<VectorType>::select(
   additional_data.solver_type = type;
 }
 
+
+
 template <typename VectorType>
 NonlinearSolverSelector<VectorType>::AdditionalData::AdditionalData(
   const SolverType &      solver_type,
@@ -440,6 +512,8 @@ NonlinearSolverSelector<VectorType>::AdditionalData::AdditionalData(
   , anderson_subspace_size(anderson_subspace_size)
 {}
 
+
+
 #ifdef DEAL_II_WITH_TRILINOS
 template <typename VectorType>
 void
@@ -453,6 +527,8 @@ NonlinearSolverSelector<VectorType>::set_data(
 }
 #endif
 
+
+
 #ifdef DEAL_II_WITH_SUNDIALS
 template <typename VectorType>
 void
@@ -463,6 +539,47 @@ NonlinearSolverSelector<VectorType>::set_data(
 }
 #endif
 
+
+
+template <typename VectorType>
+void
+NonlinearSolverSelector<VectorType>::solve_with_petsc(
+  VectorType & /*initial_guess_and_solution*/)
+{
+  AssertThrow(false,
+              ExcMessage("PETSc SNES requires you to use PETSc vectors."));
+}
+
+
+
+#ifdef DEAL_II_WITH_PETSC
+template <>
+void
+NonlinearSolverSelector<PETScWrappers::MPI::Vector>::solve_with_petsc(
+  PETScWrappers::MPI::Vector &initial_guess_and_solution)
+{
+  PETScWrappers::NonlinearSolver<PETScWrappers::MPI::Vector> nonlinear_solver(
+    additional_data_petsc_snes, mpi_communicator);
+
+  nonlinear_solver.residual = residual;
+
+  nonlinear_solver.setup_jacobian = setup_jacobian;
+
+  nonlinear_solver.solve_with_jacobian =
+    [&](const PETScWrappers::MPI::Vector &src,
+        PETScWrappers::MPI::Vector &      dst) -> int {
+    // PETSc does not gives a tolerance, so we have to choose something
+    // reasonable to provide to the user:
+    const double tolerance = 1e-6;
+    return solve_with_jacobian(src, dst, tolerance);
+  };
+
+  nonlinear_solver.solve(initial_guess_and_solution);
+}
+#endif
+
+
+
 template <typename VectorType>
 void
 NonlinearSolverSelector<VectorType>::solve(
@@ -472,6 +589,9 @@ NonlinearSolverSelector<VectorType>::solve(
   // available then we will use NOX.
   if (additional_data.solver_type == AdditionalData::SolverType::automatic)
     {
+#ifdef DEAL_II_WITH_PETSC
+      additional_data.solver_type = AdditionalData::SolverType::petsc_snes;
+#endif
 #ifdef DEAL_II_WITH_TRILINOS
       additional_data.solver_type = AdditionalData::SolverType::nox;
 #endif
@@ -490,16 +610,13 @@ NonlinearSolverSelector<VectorType>::solve(
       SUNDIALS::KINSOL<VectorType> nonlinear_solver(additional_data_kinsol,
                                                     mpi_communicator);
 
-      // We set the KINSOL reinit vector equal to the same function
-      // defined for NonlinearSolverSelector.
       nonlinear_solver.reinit_vector = reinit_vector;
-
-      nonlinear_solver.residual = residual;
+      nonlinear_solver.residual      = residual;
 
       // We cannot simply set these two functions equal to each other
       // because they have a different number of inputs.
       nonlinear_solver.setup_jacobian = [&](const VectorType &current_u,
-                                            const VectorType /*&current_f*/) {
+                                            const VectorType & /*current_f*/) {
         return NonlinearSolverSelector<VectorType>::setup_jacobian(current_u);
       };
 
@@ -514,6 +631,7 @@ NonlinearSolverSelector<VectorType>::solve(
   else if (additional_data.solver_type == AdditionalData::SolverType::nox)
     {
 #ifdef DEAL_II_WITH_TRILINOS
+
       TrilinosWrappers::NOXSolver<VectorType> nonlinear_solver(
         additional_data_nox, parameters_nox);
 
@@ -532,27 +650,33 @@ NonlinearSolverSelector<VectorType>::solve(
         false, ExcMessage("You do not have Trilinos configured with deal.II"));
 #endif
     }
+  else if (additional_data.solver_type ==
+           AdditionalData::SolverType::petsc_snes)
+    {
+      // Forward to internal function specializations, which throws for
+      // non-supported vector types:
+      solve_with_petsc(initial_guess_and_solution);
+    }
   else
     {
-      std::string solver1;
-      std::string solver2;
-
+      const std::string solvers =
 #ifdef DEAL_II_WITH_SUNDIALS
-      solver1 = "kinsol \n";
+        "kinsol\n"
 #endif
 #ifdef DEAL_II_WITH_TRILINOS
-      solver2 = "nox \n";
+        "NOX\n"
 #endif
+#ifdef DEAL_II_WITH_PETSC
+        "SNES\n"
+#endif
+        ;
 
-      DeclException2(InvalidNonlinearSolver,
-                     std::string,
-                     std::string,
-                     "Invalid nonlinear solver specified, you may use:\n"
-                       << arg1 << arg2);
-
-      AssertThrow(false, InvalidNonlinearSolver(solver1, solver2));
+      AssertThrow(false,
+                  ExcMessage(
+                    "Invalid nonlinear solver specified. "
+                    "The solvers available in your installation are:\n" +
+                    solvers));
     }
 }
 
-
 DEAL_II_NAMESPACE_CLOSE
index 266e6fe7fb3627db26caaf2f55bb239f8ce85c5c..e3126ba98b3af5f7ac265f51243bfab18f5ee928 100644 (file)
@@ -258,6 +258,7 @@ namespace nonlinear_solver_selector_test
     Vector<double> &      residual)
   {
     deallog << "  Computing residual vector..." << std::flush;
+    residual = 0.;
 
     const QGauss<dim> quadrature_formula(fe.degree + 1);
     FEValues<dim>     fe_values(fe,
index 036aa753ab4c86f303637d8f566fbe6ec1afc844..9026f572315f330aff127999ba99808749a6c841 100644 (file)
@@ -15,8 +15,8 @@
 
 
 
-// Tests the NonlinearSolverSelector class using an examplebased on the
-// step-77 tutorial program. This test checks the compatability of the
+// Tests the NonlinearSolverSelector class using an example based on the
+// step-77 tutorial program. This test checks the compatibility of the
 // class with MPI.
 
 #include <deal.II/base/function.h>
 // Included from step-40
 #include <deal.II/lac/generic_linear_algebra.h>
 
-#define FORCE_USE_OF_TRILINOS
-
 namespace LA
 {
-#if defined(DEAL_II_WITH_PETSC) && !defined(DEAL_II_PETSC_WITH_COMPLEX) && \
-  !(defined(DEAL_II_WITH_TRILINOS) && defined(FORCE_USE_OF_TRILINOS))
+#if defined(DEAL_II_WITH_PETSC) && \
+  (!defined(DEAL_II_WITH_TRILINOS) || defined(FORCE_USE_OF_PETSC))
   using namespace dealii::LinearAlgebraPETSc;
 #  define USE_PETSC_LA
 #elif defined(DEAL_II_WITH_TRILINOS)
@@ -190,6 +188,8 @@ namespace MPI_nonlinear_solver_selector_test
 
           for (const auto &boundary_value : boundary_values)
             current_solution(boundary_value.first) = boundary_value.second;
+
+          current_solution.compress(VectorOperation::insert);
         }
 
         {
@@ -316,6 +316,7 @@ namespace MPI_nonlinear_solver_selector_test
     LA::MPI::Vector &      residual)
   {
     deallog << "  Computing residual vector..." << std::flush;
+    residual = 0.;
 
     LA::MPI::Vector evaluation_point_1;
     evaluation_point_1.reinit(locally_owned_dofs,
@@ -370,8 +371,8 @@ namespace MPI_nonlinear_solver_selector_test
           }
       }
 
-    zero_constraints.set_zero(residual);
     residual.compress(VectorOperation::add);
+    zero_constraints.set_zero(residual);
 
     deallog << " norm=" << residual.l2_norm() << std::endl;
   }
similarity index 93%
rename from tests/numerics/nonlinear_solver_selector_04.cc
rename to tests/numerics/nonlinear_solver_selector_03_nox.cc
index 7173a5f45ed90a4a9b0e52c7fd6ff059f2d53b79..8517914b5f36c271330c33776d2a5d9ee3273e6a 100644 (file)
@@ -16,7 +16,7 @@
 
 
 // Tests the class NonlinearSolverSelector using an example based on
-// the test nonlinear_solver_selector_01. Here we use the nonlinear
+// the test nonlinear_solver_selector_03. Here we use the nonlinear
 // solver NOX instead of KINSOL with MPI.
 
 #define SOLVER NonlinearSolverSelector<LA::MPI::Vector>::AdditionalData::nox
diff --git a/tests/numerics/nonlinear_solver_selector_03_petsc.cc b/tests/numerics/nonlinear_solver_selector_03_petsc.cc
new file mode 100644 (file)
index 0000000..d6fd6e6
--- /dev/null
@@ -0,0 +1,27 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2007 - 2023 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// Tests the class NonlinearSolverSelector using an example based on
+// the test nonlinear_solver_selector_03. Here we use the nonlinear
+// PETSc solver SNES instead of KINSOL with MPI.
+
+#define SOLVER \
+  NonlinearSolverSelector<LA::MPI::Vector>::AdditionalData::petsc_snes
+
+#define FORCE_USE_OF_PETSC
+
+#include "nonlinear_solver_selector_03.cc"
diff --git a/tests/numerics/nonlinear_solver_selector_03_petsc.with_petsc=on.mpirun=4.output b/tests/numerics/nonlinear_solver_selector_03_petsc.with_petsc=on.mpirun=4.output
new file mode 100644 (file)
index 0000000..f46bfa6
--- /dev/null
@@ -0,0 +1,254 @@
+
+DEAL::Running with PETSc on 4 MPI rank(s)...
+DEAL::  Target_tolerance: 0.00100000
+DEAL::
+DEAL::  Computing residual vector... norm=0.170956
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 2.30748
+DEAL::Convergence step 9 value 3.07002e-13
+DEAL::   Solved in 9 iterations.
+DEAL::  Computing residual vector... norm=0.170956
+DEAL::  Computing residual vector... norm=0.129347
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 1.67096
+DEAL::Convergence step 9 value 2.69526e-13
+DEAL::   Solved in 9 iterations.
+DEAL::  Computing residual vector... norm=0.129347
+DEAL::  Computing residual vector... norm=0.100710
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 1.37283
+DEAL::Convergence step 9 value 1.64162e-13
+DEAL::   Solved in 9 iterations.
+DEAL::  Computing residual vector... norm=0.100710
+DEAL::  Computing residual vector... norm=0.0776513
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.802384
+DEAL::Convergence step 9 value 2.11842e-13
+DEAL::   Solved in 9 iterations.
+DEAL::  Computing residual vector... norm=0.0776513
+DEAL::  Computing residual vector... norm=0.0593571
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.409330
+DEAL::Convergence step 9 value 6.54958e-14
+DEAL::   Solved in 9 iterations.
+DEAL::  Computing residual vector... norm=0.0593571
+DEAL::  Computing residual vector... norm=0.0459878
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.290541
+DEAL::Convergence step 9 value 7.18389e-14
+DEAL::   Solved in 9 iterations.
+DEAL::  Computing residual vector... norm=0.0459878
+DEAL::  Computing residual vector... norm=0.0355434
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.222945
+DEAL::Convergence step 9 value 5.02084e-14
+DEAL::   Solved in 9 iterations.
+DEAL::  Computing residual vector... norm=0.0355434
+DEAL::  Computing residual vector... norm=0.0274503
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.171509
+DEAL::Convergence step 9 value 5.52278e-14
+DEAL::   Solved in 9 iterations.
+DEAL::  Computing residual vector... norm=0.0274503
+DEAL::  Computing residual vector... norm=0.0212126
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.134719
+DEAL::Convergence step 9 value 4.72412e-14
+DEAL::   Solved in 9 iterations.
+DEAL::  Computing residual vector... norm=0.0212126
+DEAL::  Computing residual vector... norm=0.0164293
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.106579
+DEAL::Convergence step 8 value 9.61276e-13
+DEAL::   Solved in 8 iterations.
+DEAL::  Computing residual vector... norm=0.0164293
+DEAL::  Computing residual vector... norm=0.0127897
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.0845457
+DEAL::Convergence step 8 value 6.79505e-13
+DEAL::   Solved in 8 iterations.
+DEAL::  Computing residual vector... norm=0.0127897
+DEAL::  Computing residual vector... norm=0.0100451
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.0671822
+DEAL::Convergence step 8 value 5.29774e-13
+DEAL::   Solved in 8 iterations.
+DEAL::  Computing residual vector... norm=0.0100451
+DEAL::  Computing residual vector... norm=0.00798892
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.0531083
+DEAL::Convergence step 8 value 6.78081e-13
+DEAL::   Solved in 8 iterations.
+DEAL::  Computing residual vector... norm=0.00798892
+DEAL::  Computing residual vector... norm=0.00645103
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.0417547
+DEAL::Convergence step 8 value 2.20078e-13
+DEAL::   Solved in 8 iterations.
+DEAL::  Computing residual vector... norm=0.00645103
+DEAL::  Computing residual vector... norm=0.00529638
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.0326761
+DEAL::Convergence step 8 value 8.07466e-14
+DEAL::   Solved in 8 iterations.
+DEAL::  Computing residual vector... norm=0.00529638
+DEAL::  Computing residual vector... norm=0.00442194
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.0254397
+DEAL::Convergence step 8 value 2.27458e-13
+DEAL::   Solved in 8 iterations.
+DEAL::  Computing residual vector... norm=0.00442194
+DEAL::  Computing residual vector... norm=0.00375154
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.0197527
+DEAL::Convergence step 8 value 1.05240e-13
+DEAL::   Solved in 8 iterations.
+DEAL::  Computing residual vector... norm=0.00375154
+DEAL::  Computing residual vector... norm=0.00323013
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.0154004
+DEAL::Convergence step 7 value 9.37496e-13
+DEAL::   Solved in 7 iterations.
+DEAL::  Computing residual vector... norm=0.00323013
+DEAL::  Computing residual vector... norm=0.00281833
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.0120033
+DEAL::Convergence step 8 value 2.93204e-14
+DEAL::   Solved in 8 iterations.
+DEAL::  Computing residual vector... norm=0.00281834
+DEAL::  Computing residual vector... norm=0.00248811
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.00940651
+DEAL::Convergence step 8 value 1.96359e-14
+DEAL::   Solved in 8 iterations.
+DEAL::  Computing residual vector... norm=0.00248811
+DEAL::  Computing residual vector... norm=0.00221939
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.00741723
+DEAL::Convergence step 7 value 6.96541e-13
+DEAL::   Solved in 7 iterations.
+DEAL::  Computing residual vector... norm=0.00221939
+DEAL::  Computing residual vector... norm=0.00199769
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.00589663
+DEAL::Convergence step 8 value 1.65160e-14
+DEAL::   Solved in 8 iterations.
+DEAL::  Computing residual vector... norm=0.00199769
+DEAL::  Computing residual vector... norm=0.00181244
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.00473229
+DEAL::Convergence step 7 value 9.22009e-13
+DEAL::   Solved in 7 iterations.
+DEAL::  Computing residual vector... norm=0.00181244
+DEAL::  Computing residual vector... norm=0.00165582
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.00383393
+DEAL::Convergence step 7 value 7.05884e-13
+DEAL::   Solved in 7 iterations.
+DEAL::  Computing residual vector... norm=0.00165583
+DEAL::  Computing residual vector... norm=0.00152199
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.00313701
+DEAL::Convergence step 7 value 5.38032e-13
+DEAL::   Solved in 7 iterations.
+DEAL::  Computing residual vector... norm=0.00152199
+DEAL::  Computing residual vector... norm=0.00140649
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.00259180
+DEAL::Convergence step 7 value 4.06858e-13
+DEAL::   Solved in 7 iterations.
+DEAL::  Computing residual vector... norm=0.00140649
+DEAL::  Computing residual vector... norm=0.00130588
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.00216116
+DEAL::Convergence step 7 value 3.07418e-13
+DEAL::   Solved in 7 iterations.
+DEAL::  Computing residual vector... norm=0.00130589
+DEAL::  Computing residual vector... norm=0.00121752
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.00181099
+DEAL::Convergence step 7 value 2.90636e-13
+DEAL::   Solved in 7 iterations.
+DEAL::  Computing residual vector... norm=0.00121752
+DEAL::  Computing residual vector... norm=0.00113928
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.00154247
+DEAL::Convergence step 7 value 2.60502e-13
+DEAL::   Solved in 7 iterations.
+DEAL::  Computing residual vector... norm=0.00113929
+DEAL::  Computing residual vector... norm=0.00106952
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.00131648
+DEAL::Convergence step 7 value 2.17370e-13
+DEAL::   Solved in 7 iterations.
+DEAL::  Computing residual vector... norm=0.00106952
+DEAL::  Computing residual vector... norm=0.00100687
+DEAL::  Computing Jacobian matrix
+DEAL::  Factorizing Jacobian matrix
+DEAL::  Solving linear system
+DEAL::Starting value 0.00113050
+DEAL::Convergence step 7 value 1.81968e-13
+DEAL::   Solved in 7 iterations.
+DEAL::  Computing residual vector... norm=0.00100688
+DEAL::  Computing residual vector... norm=0.000950282
+DEAL::

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.