From 02f208a0b2480956ddbaaa6fce3bed8455f8093e Mon Sep 17 00:00:00 2001 From: Vladimir Ivannikov Date: Thu, 28 Dec 2023 09:41:23 +0100 Subject: [PATCH] Add tests --- tests/trilinos/nox_solver_05.cc | 144 ++++++++++++++++++ ...lver_05.with_trilinos_with_nox=true.output | 28 ++++ tests/trilinos/nox_solver_06.cc | 140 +++++++++++++++++ ...lver_06.with_trilinos_with_nox=true.output | 21 +++ 4 files changed, 333 insertions(+) create mode 100644 tests/trilinos/nox_solver_05.cc create mode 100644 tests/trilinos/nox_solver_05.with_trilinos_with_nox=true.output create mode 100644 tests/trilinos/nox_solver_06.cc create mode 100644 tests/trilinos/nox_solver_06.with_trilinos_with_nox=true.output diff --git a/tests/trilinos/nox_solver_05.cc b/tests/trilinos/nox_solver_05.cc new file mode 100644 index 0000000000..fd41079f17 --- /dev/null +++ b/tests/trilinos/nox_solver_05.cc @@ -0,0 +1,144 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2022 - 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. +// +// --------------------------------------------------------------------- + + + +// Check TrilinosWrappers::NOXSolver. This test solves the same +// problem as nox_solver_03, but it throws an exception in solve_with_jacobian. +// This callback is special, since with default settings NOX tries to perform a +// recovery step once an exception has been raised there. We need to treat such +// situations accordingly. This test verifies that the recovery functionality of +// NOX works. + + +#include + +#include +#include + +#include + +#include "../tests.h" + + +int +main(int argc, char **argv) +{ + initlog(); + + Utilities::MPI::MPI_InitFinalize mpi_init(argc, argv, 1); + + using Number = double; + using VectorType = Vector; + + // set up solver control + const unsigned int n_max_iterations = 10; + const double abs_tolerance = 1e-9; + const double rel_tolerance = 1e-5; + + TrilinosWrappers::NOXSolver::AdditionalData additional_data( + n_max_iterations, abs_tolerance, rel_tolerance); + + // set up parameters + Teuchos::RCP non_linear_parameters = + Teuchos::rcp(new Teuchos::ParameterList); + + non_linear_parameters->set("Nonlinear Solver", "Line Search Based"); + non_linear_parameters->sublist("Printing").set("Output Information", 15); + non_linear_parameters->sublist("Direction").set("Method", "Newton"); + non_linear_parameters->sublist("Line Search").set("Method", "Polynomial"); + + + // set up solver + TrilinosWrappers::NOXSolver solver(additional_data, + non_linear_parameters); + + // ... helper functions + solver.residual = [](const VectorType &u, VectorType &F) { + deallog << "Evaluating the solution at u=(" << u[0] << ',' << u[1] << ')' + << std::endl; + + F(0) = std::cos(u[0] + u[1]) - 1 + 2 * u[0]; + F(1) = std::sin(u[0] - u[1]) + 2 * u[1]; + }; + + FullMatrix J(2, 2); + FullMatrix J_inverse(2, 2); + + solver.setup_jacobian = [&J, &J_inverse](const VectorType &u) { + // We don't do any kind of set-up in this program, but we can at least + // say that we're here + deallog << "Setting up Jacobian system at u=(" << u[0] << ',' << u[1] << ')' + << std::endl; + + J(0, 0) = -std::sin(u[0] + u[1]) + 2; + J(0, 1) = -std::sin(u[0] + u[1]); + J(1, 0) = std::cos(u[0] - u[1]); + J(1, 1) = -std::cos(u[0] - u[1]) + 2; + + J_inverse.invert(J); + }; + + solver.apply_jacobian = [&](const VectorType &src, VectorType &dst) { + J.vmult(dst, src); + }; + + unsigned int counter = 0; + + solver.solve_with_jacobian = [&J_inverse, + &counter](const VectorType &rhs, + VectorType &dst, + const double /*tolerance*/) { + deallog << "Solving Jacobian system with rhs=(" << rhs[0] << ',' << rhs[1] + << ')' << std::endl; + + J_inverse.vmult(dst, rhs); + + counter++; + + // Fail linear solver step but NOX should perform a recovery step + if (counter == 2) + throw RecoverableUserCallbackError( + "Recoverable failure in solve_with_jacobian."); + + // Throw an irrecoverable exception then + if (counter == 5) + throw ExcMessage("Irrecoverable failure in solve_with_jacobian."); + }; + + // initial guess + const unsigned int N = 2; + VectorType solution(N); + solution[0] = 0.5; + solution[1] = 1.234; + + // solve with the given initial guess + try + { + auto niter = solver.solve(solution); + solution.print(deallog.get_file_stream()); + deallog << "Converged in " << niter << " iterations." << std::endl; + } + catch (const std::exception &exc) + { + deallog << "Caught an irrecoverable exception in a callback:" << std::endl + << exc.what() << std::endl; + } + catch (const char *s) + { + deallog << "Trilinos threw its own exception of type char*:" << std::endl + << s << std::endl; + } +} diff --git a/tests/trilinos/nox_solver_05.with_trilinos_with_nox=true.output b/tests/trilinos/nox_solver_05.with_trilinos_with_nox=true.output new file mode 100644 index 0000000000..bd6fab75d1 --- /dev/null +++ b/tests/trilinos/nox_solver_05.with_trilinos_with_nox=true.output @@ -0,0 +1,28 @@ + +DEAL::Evaluating the solution at u=(0.500000,1.23400) +DEAL::Setting up Jacobian system at u=(0.500000,1.23400) +DEAL::Solving Jacobian system with rhs=(-0.162480,1.79816) +DEAL::Evaluating the solution at u=(-0.282294,0.265967) +DEAL::Setting up Jacobian system at u=(-0.282294,0.265967) +DEAL::Solving Jacobian system with rhs=(-0.564722,0.0107297) +DEAL::Evaluating the solution at u=(-0.564143,0.485115) +DEAL::Evaluating the solution at u=(-0.564143,0.485115) +DEAL::Setting up Jacobian system at u=(-0.564143,0.485115) +DEAL::Solving Jacobian system with rhs=(-1.13141,0.103177) +DEAL::Evaluating the solution at u=(-0.0103358,0.232689) +DEAL::Setting up Jacobian system at u=(-0.0103358,0.232689) +DEAL::Solving Jacobian system with rhs=(-0.0452903,0.224738) +DEAL::Evaluating the solution at u=(-0.0117725,0.0157208) +DEAL::Setting up Jacobian system at u=(-0.0117725,0.0157208) +DEAL::Solving Jacobian system with rhs=(-0.0235528,0.00395177) +DEAL::Caught an irrecoverable exception in a callback: +DEAL:: +-------------------------------------------------------- +An error occurred in line <0> of file <> in function + +The violated condition was: + +Additional information: + Irrecoverable failure in solve_with_jacobian. +-------------------------------------------------------- + diff --git a/tests/trilinos/nox_solver_06.cc b/tests/trilinos/nox_solver_06.cc new file mode 100644 index 0000000000..eb5373c79f --- /dev/null +++ b/tests/trilinos/nox_solver_06.cc @@ -0,0 +1,140 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2022 - 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. +// +// --------------------------------------------------------------------- + + + +// Check TrilinosWrappers::NOXSolver. This test solves the same +// problem as nox_solver_03 and checks the same special case of +// solve_with_jacobian as nox_solver_05, but this time we disable the recovery +// functionality. So the NOX wrapper is expected to behave like in +// nox_solver_04. + + +#include + +#include +#include + +#include + +#include "../tests.h" + + +int +main(int argc, char **argv) +{ + initlog(); + + Utilities::MPI::MPI_InitFinalize mpi_init(argc, argv, 1); + + using Number = double; + using VectorType = Vector; + + // set up solver control + const unsigned int n_max_iterations = 10; + const double abs_tolerance = 1e-9; + const double rel_tolerance = 1e-5; + + TrilinosWrappers::NOXSolver::AdditionalData additional_data( + n_max_iterations, abs_tolerance, rel_tolerance); + + // set up parameters + Teuchos::RCP non_linear_parameters = + Teuchos::rcp(new Teuchos::ParameterList); + + non_linear_parameters->set("Nonlinear Solver", "Line Search Based"); + non_linear_parameters->sublist("Printing").set("Output Information", 15); + non_linear_parameters->sublist("Direction").set("Method", "Newton"); + non_linear_parameters->sublist("Direction") + .sublist("Newton") + .set("Rescue Bad Newton Solve", false); + non_linear_parameters->sublist("Line Search").set("Method", "Polynomial"); + + + // set up solver + TrilinosWrappers::NOXSolver solver(additional_data, + non_linear_parameters); + + // ... helper functions + solver.residual = [](const VectorType &u, VectorType &F) { + deallog << "Evaluating the solution at u=(" << u[0] << ',' << u[1] << ')' + << std::endl; + + F(0) = std::cos(u[0] + u[1]) - 1 + 2 * u[0]; + F(1) = std::sin(u[0] - u[1]) + 2 * u[1]; + }; + + FullMatrix J(2, 2); + FullMatrix J_inverse(2, 2); + + solver.setup_jacobian = [&J, &J_inverse](const VectorType &u) { + // We don't do any kind of set-up in this program, but we can at least + // say that we're here + deallog << "Setting up Jacobian system at u=(" << u[0] << ',' << u[1] << ')' + << std::endl; + + J(0, 0) = -std::sin(u[0] + u[1]) + 2; + J(0, 1) = -std::sin(u[0] + u[1]); + J(1, 0) = std::cos(u[0] - u[1]); + J(1, 1) = -std::cos(u[0] - u[1]) + 2; + + J_inverse.invert(J); + }; + + solver.apply_jacobian = [&](const VectorType &src, VectorType &dst) { + J.vmult(dst, src); + }; + + unsigned int counter = 0; + + solver.solve_with_jacobian = [&J_inverse, + &counter](const VectorType &rhs, + VectorType &dst, + const double /*tolerance*/) { + deallog << "Solving Jacobian system with rhs=(" << rhs[0] << ',' << rhs[1] + << ')' << std::endl; + + J_inverse.vmult(dst, rhs); + + counter++; + + if (counter == 3) + throw ExcMessage("Irrecoverable failure in solve_with_jacobian."); + }; + + // initial guess + const unsigned int N = 2; + VectorType solution(N); + solution[0] = 0.5; + solution[1] = 1.234; + + // solve with the given initial guess + try + { + auto niter = solver.solve(solution); + solution.print(deallog.get_file_stream()); + deallog << "Converged in " << niter << " iterations." << std::endl; + } + catch (const std::exception &exc) + { + deallog << "Caught an irrecoverable exception in a callback:" << std::endl + << exc.what() << std::endl; + } + catch (const char *s) + { + deallog << "Trilinos threw its own exception of type char*:" << std::endl + << s << std::endl; + } +} diff --git a/tests/trilinos/nox_solver_06.with_trilinos_with_nox=true.output b/tests/trilinos/nox_solver_06.with_trilinos_with_nox=true.output new file mode 100644 index 0000000000..eb8e195fac --- /dev/null +++ b/tests/trilinos/nox_solver_06.with_trilinos_with_nox=true.output @@ -0,0 +1,21 @@ + +DEAL::Evaluating the solution at u=(0.500000,1.23400) +DEAL::Setting up Jacobian system at u=(0.500000,1.23400) +DEAL::Solving Jacobian system with rhs=(-0.162480,1.79816) +DEAL::Evaluating the solution at u=(-0.282294,0.265967) +DEAL::Setting up Jacobian system at u=(-0.282294,0.265967) +DEAL::Solving Jacobian system with rhs=(-0.564722,0.0107297) +DEAL::Evaluating the solution at u=(-0.000445202,0.0468183) +DEAL::Setting up Jacobian system at u=(-0.000445202,0.0468183) +DEAL::Solving Jacobian system with rhs=(-0.00196544,0.0463907) +DEAL::Caught an irrecoverable exception in a callback: +DEAL:: +-------------------------------------------------------- +An error occurred in line <0> of file <> in function + +The violated condition was: + +Additional information: + Irrecoverable failure in solve_with_jacobian. +-------------------------------------------------------- + -- 2.39.5