From 80ef18ac6fed499a7bed7f4e5bd6e5ae210f8e6a Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Tue, 6 May 2025 14:38:52 -0600 Subject: [PATCH] Add a testcase. --- tests/sundials/kinsol_06_v3.cc | 102 ++++++++++++++++++++++++++ tests/sundials/kinsol_06_v3.output | 23 ++++++ tests/sundials/kinsol_06_v3_input.prm | 19 +++++ 3 files changed, 144 insertions(+) create mode 100644 tests/sundials/kinsol_06_v3.cc create mode 100644 tests/sundials/kinsol_06_v3.output create mode 100644 tests/sundials/kinsol_06_v3_input.prm diff --git a/tests/sundials/kinsol_06_v3.cc b/tests/sundials/kinsol_06_v3.cc new file mode 100644 index 0000000000..97e0afdae2 --- /dev/null +++ b/tests/sundials/kinsol_06_v3.cc @@ -0,0 +1,102 @@ +// ------------------------------------------------------------------------ +// +// 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 kinsol_06 test, but with a case where KINSOL, after a few +// tries with the given residual/Jacobian, 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() +{ + initlog(); + + using VectorType = Vector; + + // Size of the problem + const unsigned int N = 1; + + SUNDIALS::KINSOL::AdditionalData data; + ParameterHandler prm; + data.add_parameters(prm); + + std::ifstream ifile(SOURCE_DIR "/kinsol_06_v3_input.prm"); + prm.parse_input(ifile); + + SUNDIALS::KINSOL kinsol(data); + + kinsol.reinit_vector = [N](VectorType &v) { v.reinit(N); }; + + kinsol.residual = [&](const VectorType &u, VectorType &F) { + // Count how many times this function has been called. Let it call + // *every* time after the first attempt: + 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)) + { + deallog << "Reporting recoverable failure." << std::endl; + throw RecoverableUserCallbackError(); + } + count++; + + + F.reinit(u); + F[0] = std::atan(u[0]) - 0.5; + }; + + double J_inverse; + + kinsol.setup_jacobian = [&J_inverse](const VectorType &u, + const VectorType &F) { + deallog << "Setting up Jacobian system at u=" << u[0] << std::endl; + + const double J = 1. / (1 + u[0] * u[0]); + J_inverse = 1. / J; + }; + + + kinsol.solve_with_jacobian = [&](const VectorType &rhs, + VectorType &dst, + double) { dst[0] = J_inverse * rhs[0]; }; + + VectorType v(N); + v[0] = 10; + + try + { + auto niter = kinsol.solve(v); + } + catch (const ExceptionBase &e) + { + deallog << "KINSOL threw an exception with the following message:" + << std::endl; + e.print_info(deallog.get_file_stream()); + } +} diff --git a/tests/sundials/kinsol_06_v3.output b/tests/sundials/kinsol_06_v3.output new file mode 100644 index 0000000000..dda6d57d6f --- /dev/null +++ b/tests/sundials/kinsol_06_v3.output @@ -0,0 +1,23 @@ + +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::Setting up Jacobian system at u=-2.26049 +DEAL::Computing residual for the 4th time, at u=-2.26049 +DEAL::Computing residual for the 5th time, at u=7.84693 +DEAL::Computing residual for the 6th time, at u=7.84693 +DEAL::Computing residual for the 7th time, at u=2.07902 +DEAL::Setting up Jacobian system at u=2.07902 +DEAL::Computing residual for the 8th time, at u=2.07902 +DEAL::Computing residual for the 9th time, at u=-1.23396 +DEAL::Computing residual for the 10th time, at u=-1.23396 +DEAL::Computing residual for the 11th time, at u=6.16274 +DEAL::KINSOL threw an exception with the following message: + One of the SUNDIALS KINSOL internal functions returned a negative error code: -6. Please consult SUNDIALS manual. diff --git a/tests/sundials/kinsol_06_v3_input.prm b/tests/sundials/kinsol_06_v3_input.prm new file mode 100644 index 0000000000..fd9f762886 --- /dev/null +++ b/tests/sundials/kinsol_06_v3_input.prm @@ -0,0 +1,19 @@ +set Function norm stopping tolerance = 0.000000 +set Maximum number of nonlinear iterations = 5 +set Scaled step stopping tolerance = 0.000000 +set Solution strategy = newton + +subsection Fixed point and Picard parameters + set Anderson acceleration subspace size = 5 +end + +subsection Linesearch parameters + set Maximum number of beta-condition failures = 0 +end + +subsection Newton parameters + set Maximum allowable scaled length of the Newton step = 0.000000 + set Maximum iterations without matrix setup = 0 + set No initial matrix setup = false + set Relative error for different quotient computation = 0.000000 +end -- 2.39.5