]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add a testcase for SNES. 18443/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Fri, 9 May 2025 23:54:04 +0000 (17:54 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Mon, 19 May 2025 16:19:32 +0000 (10:19 -0600)
tests/petsc/petsc_snes_08.cc [new file with mode: 0644]
tests/petsc/petsc_snes_08.output [new file with mode: 0644]

diff --git a/tests/petsc/petsc_snes_08.cc b/tests/petsc/petsc_snes_08.cc
new file mode 100644 (file)
index 0000000..e598a7c
--- /dev/null
@@ -0,0 +1,107 @@
+// ------------------------------------------------------------------------
+//
+// SPDX-License-Identifier: LGPL-2.1-or-later
+// Copyright (C) 2017 - 2023 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// Part of the source code is dual licensed under Apache-2.0 WITH
+// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
+// governing the source code and code contributions can be found in
+// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
+//
+// ------------------------------------------------------------------------
+
+#include <deal.II/base/utilities.h>
+
+#include <deal.II/lac/petsc_vector.h>
+
+#include <deal.II/numerics/nonlinear.h>
+
+#include "../tests.h"
+
+
+// Like the sundials/kinsol_06_v3 test, but for SNES. After a few
+// tries with the given residual/Jacobian, SNES gives up and
+// terminates because it can't find a solution. We need to make sure
+// that we throw a catchable exception.
+//
+// This testcase is a variation of the kinsol_06 test, modified by
+// Simon Wiesheier, and posted on the mailing list. Then further
+// adapted by Wolfgang Bangerth.
+
+
+int
+main(int argc, char *argv[])
+{
+  initlog();
+
+  using VectorType = PETScWrappers::MPI::Vector;
+
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+  // Size of the problem
+  const unsigned int N = 1;
+
+  NonlinearSolverSelector<VectorType>::AdditionalData additional_data;
+  additional_data.solver_type =
+    NonlinearSolverSelector<VectorType>::AdditionalData::SolverType::petsc_snes;
+  additional_data.maximum_non_linear_iterations = 1;
+
+  NonlinearSolverSelector<VectorType> nonlinear_solver(additional_data);
+
+  nonlinear_solver.reinit_vector = [N](VectorType &v) {
+    v.reinit(MPI_COMM_WORLD, N, N);
+  };
+
+  nonlinear_solver.residual = [&](const VectorType &u, VectorType &F) {
+    // Count how many times this function has been called.
+    static int count = 0;
+    deallog << "Computing residual for the " << count + 1
+            << "th time, at u=" << u[0] << std::endl;
+    if ((u[0] < -10) || (u[0] > 20) /* || (count > 0) */)
+      {
+        deallog << "Reporting recoverable failure." << std::endl;
+        throw RecoverableUserCallbackError();
+      }
+    count++;
+
+
+    F.reinit(u);
+    F[0] = std::atan(u[0]) - 0.5;
+    F.compress(VectorOperation::values::insert);
+  };
+
+  double J_inverse;
+
+  nonlinear_solver.setup_jacobian = [&J_inverse](const VectorType &u) {
+    deallog << "Setting up Jacobian system at u=" << u[0] << std::endl;
+
+    const double J = 1. / (1 + u[0] * u[0]);
+    J_inverse      = 1. / J;
+  };
+
+
+  nonlinear_solver.solve_with_jacobian =
+    [&](const VectorType &rhs, VectorType &dst, double) {
+      dst[0] = J_inverse * rhs[0];
+      dst.compress(VectorOperation::values::insert);
+    };
+
+
+  VectorType v(MPI_COMM_WORLD, N, N);
+  v[0] = 10;
+  v.compress(VectorOperation::values::insert);
+
+  try
+    {
+      nonlinear_solver.solve(v);
+    }
+  catch (const ExceptionBase &e)
+    {
+      deallog
+        << "Nonlinear Solver threw an exception with the following message:"
+        << std::endl;
+      e.print_info(deallog.get_file_stream());
+    }
+}
diff --git a/tests/petsc/petsc_snes_08.output b/tests/petsc/petsc_snes_08.output
new file mode 100644 (file)
index 0000000..13eb5f6
--- /dev/null
@@ -0,0 +1,14 @@
+
+DEAL::Computing residual for the 1th time, at u=10.0000
+DEAL::Setting up Jacobian system at u=10.0000
+DEAL::Computing residual for the 2th time, at u=10.0000
+DEAL::Computing residual for the 3th time, at u=-88.0839
+DEAL::Reporting recoverable failure.
+DEAL::Computing residual for the 3th time, at u=-39.0419
+DEAL::Reporting recoverable failure.
+DEAL::Computing residual for the 3th time, at u=-14.5210
+DEAL::Reporting recoverable failure.
+DEAL::Computing residual for the 3th time, at u=-2.26049
+DEAL::Computing residual for the 4th time, at u=8.77395
+DEAL::Nonlinear Solver threw an exception with the following message:
+    SNES solver did not converge after 1 iterations with reason DIVERGED_MAX_IT

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.