From: Matthias Maier Date: Mon, 21 Apr 2025 19:23:27 +0000 (-0500) Subject: SolverMPGMRES: add another test X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=839b8ff446c180d7694cdfc01c6de89606225dc2;p=dealii.git SolverMPGMRES: add another test --- diff --git a/include/deal.II/lac/solver_gmres.h b/include/deal.II/lac/solver_gmres.h index 632ffe27d1..0eabb9e018 100644 --- a/include/deal.II/lac/solver_gmres.h +++ b/include/deal.II/lac/solver_gmres.h @@ -2240,6 +2240,13 @@ void SolverMPGMRES::solve_internal( const auto previous_vector_index = [this, n_preconditioners, indexing_strategy]( unsigned int i) -> unsigned int { + // In the special case of no preconditioners we simply fall back to the + // FGMRES indexing strategy. + if (n_preconditioners == 0) + { + return i; + } + switch (indexing_strategy) { case IndexingStrategy::fgmres: diff --git a/tests/lac/solver_mpgmres_02.cc b/tests/lac/solver_mpgmres_02.cc new file mode 100644 index 0000000000..6e98ca333d --- /dev/null +++ b/tests/lac/solver_mpgmres_02.cc @@ -0,0 +1,64 @@ +// ------------------------------------------------------------------------ +// +// 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 we can call FGMRES::solve() and MPGMRES::solve() without a +// preconditioner + +#include +#include +#include + +#include "../tests.h" + +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; + for (unsigned int i = N; i < M; ++i) + matrix.diag_element(i) = 1.; + + Vector rhs(M); + Vector sol(M); + rhs = 1.; + + SolverControl control(M, 1.e-2); + + deallog.get_file_stream() << "FGMRES:\n"; + SolverFGMRES> solver_fgmres(control); + check_solver_within_range(solver_fgmres.solve(matrix, sol, rhs), + control.last_step(), + 1, + 10); + + deallog.get_file_stream() << "\nMPGMRES:\n"; + sol = 0.; + + SolverMPGMRES> solver_mpgmres(control); + check_solver_within_range(solver_mpgmres.solve(matrix, sol, rhs), + control.last_step(), + 1, + 10); +} diff --git a/tests/lac/solver_mpgmres_02.output b/tests/lac/solver_mpgmres_02.output new file mode 100644 index 0000000000..7a7440e9fa --- /dev/null +++ b/tests/lac/solver_mpgmres_02.output @@ -0,0 +1,6 @@ + +FGMRES: +DEAL::Solver stopped within 1 - 10 iterations + +MPGMRES: +DEAL::Solver stopped within 1 - 10 iterations