From: Wolfgang Bangerth Date: Fri, 9 May 2025 23:54:04 +0000 (-0600) Subject: Add a testcase for SNES. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=c5c7201d3ae184bf4cf534af11d06dbce178df83;p=dealii.git Add a testcase for SNES. --- diff --git a/tests/petsc/petsc_snes_08.cc b/tests/petsc/petsc_snes_08.cc new file mode 100644 index 0000000000..e598a7ca8f --- /dev/null +++ b/tests/petsc/petsc_snes_08.cc @@ -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 + +#include + +#include + +#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::AdditionalData additional_data; + additional_data.solver_type = + NonlinearSolverSelector::AdditionalData::SolverType::petsc_snes; + additional_data.maximum_non_linear_iterations = 1; + + NonlinearSolverSelector 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 index 0000000000..13eb5f6517 --- /dev/null +++ b/tests/petsc/petsc_snes_08.output @@ -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