From 372deb52aefb65775ee76a85f35e18bf7b5448bf Mon Sep 17 00:00:00 2001 From: Wyatt Smith Date: Tue, 1 Apr 2025 14:30:17 -0500 Subject: [PATCH] added test for MPGMRES --- tests/lac/solver_mpgmres_01.cc | 141 +++++++++++++++++++++++++++++ tests/lac/solver_mpgmres_01.output | 34 +++++++ 2 files changed, 175 insertions(+) create mode 100644 tests/lac/solver_mpgmres_01.cc create mode 100644 tests/lac/solver_mpgmres_01.output diff --git a/tests/lac/solver_mpgmres_01.cc b/tests/lac/solver_mpgmres_01.cc new file mode 100644 index 0000000000..28cafb3cc3 --- /dev/null +++ b/tests/lac/solver_mpgmres_01.cc @@ -0,0 +1,141 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2022 - 2024 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. +// +// ------------------------------------------------------------------------ + +// Check that the preconditioners are being applied cyclically + +#include +#include +#include +#include + +#include + +#include "../tests.h" + + +// We create a wrapper around the preconditioner that will print out a name with +// each vmult so we can check that the preconditioners are being cycled +// correctly +template +class MyPreconditioner +{ +public: + MyPreconditioner(const std::string &name, + const Preconditioner &preconditioner) + : name(name) + , preconditioner(preconditioner) + , print_counter(0) + {} + + void + vmult(Vector &dst, const Vector &src) const + { + if (print_counter < 3) + deallog.get_file_stream() + << "Applying preconditioner " << name << std::endl; + preconditioner.vmult(dst, src); + ++print_counter; + } + + void + reset_counter() + { + print_counter = 0; + } + +private: + const std::string name; + const Preconditioner &preconditioner; + mutable unsigned int print_counter; +}; + + +int +main() +{ + initlog(); + + const unsigned int N = 500; + const unsigned int M = 2 * N; + + SparsityPattern sparsity_pattern(M, M, 2); + for (unsigned int i = 0; i < M - 1; ++i) + sparsity_pattern.add(i, i + 1); + sparsity_pattern.compress(); + SparseMatrix matrix(sparsity_pattern); + + for (unsigned int i = 0; i < N; ++i) + { + matrix.diag_element(i) = 1.0 + 1 * i; + matrix.set(i, i + 1, 1.); + } + for (unsigned int i = N; i < M; ++i) + matrix.diag_element(i) = 1.; + + PreconditionChebyshev> prec_a; + prec_a.initialize(matrix); + MyPreconditioner wrapper_a("A", prec_a); + + PreconditionJacobi> prec_b; + prec_b.initialize(matrix); + MyPreconditioner wrapper_b("B", prec_b); + + PreconditionSOR> prec_c; + prec_c.initialize(matrix); + MyPreconditioner wrapper_c("C", prec_c); + + Vector rhs(M); + Vector sol(M); + rhs = 1.; + + SolverControl control(M, 1e-8); + + deallog.get_file_stream() << "FGMRES:\n"; + SolverFGMRES> solver_fgmres(control); + check_solver_within_range( + solver_fgmres.solve(matrix, sol, rhs, wrapper_a, wrapper_b, wrapper_c), + control.last_step(), + 7 - 6, + 7 + 6); + + deallog.get_file_stream() << "\nMPGMRES: (truncated)\n"; + sol = 0.; + wrapper_a.reset_counter(); + wrapper_b.reset_counter(); + wrapper_c.reset_counter(); + + SolverMPGMRES>::AdditionalData data; + data.use_truncated_mpgmres_strategy = true; + SolverMPGMRES> solver_trmpgmres(control, data); + check_solver_within_range( + solver_trmpgmres.solve(matrix, sol, rhs, wrapper_a, wrapper_b, wrapper_c), + control.last_step(), + 19 - 6, + 19 + 6); + + + deallog.get_file_stream() << "\nMPGMRES: (full)\n"; + sol = 0.; + wrapper_a.reset_counter(); + wrapper_b.reset_counter(); + wrapper_c.reset_counter(); + data.use_truncated_mpgmres_strategy = false; + SolverMPGMRES> solver_mpgmres(control, data); + + check_solver_within_range( + solver_mpgmres.solve(matrix, sol, rhs, wrapper_a, wrapper_b, wrapper_c), + control.last_step(), + 24 - 6, + 24 + 6); +} diff --git a/tests/lac/solver_mpgmres_01.output b/tests/lac/solver_mpgmres_01.output new file mode 100644 index 0000000000..c50f13ed6c --- /dev/null +++ b/tests/lac/solver_mpgmres_01.output @@ -0,0 +1,34 @@ + +FGMRES: +Applying preconditioner A +Applying preconditioner B +Applying preconditioner C +Applying preconditioner A +Applying preconditioner B +Applying preconditioner C +Applying preconditioner A +DEAL::Solver stopped within 1 - 13 iterations + +MPGMRES: (truncated) +Applying preconditioner A +Applying preconditioner B +Applying preconditioner C +Applying preconditioner A +Applying preconditioner B +Applying preconditioner C +Applying preconditioner A +Applying preconditioner B +Applying preconditioner C +DEAL::Solver stopped within 13 - 25 iterations + +MPGMRES: (full) +Applying preconditioner A +Applying preconditioner B +Applying preconditioner C +Applying preconditioner A +Applying preconditioner B +Applying preconditioner C +Applying preconditioner A +Applying preconditioner B +Applying preconditioner C +DEAL::Solver stopped within 18 - 30 iterations -- 2.39.5