]> https://gitweb.dealii.org/ - dealii.git/commitdiff
SolverMPGMRES: add another test
authorMatthias Maier <tamiko@43-1.org>
Mon, 21 Apr 2025 19:23:27 +0000 (14:23 -0500)
committerMatthias Maier <tamiko@43-1.org>
Mon, 5 May 2025 21:10:57 +0000 (16:10 -0500)
include/deal.II/lac/solver_gmres.h
tests/lac/solver_mpgmres_02.cc [new file with mode: 0644]
tests/lac/solver_mpgmres_02.output [new file with mode: 0644]

index 632ffe27d168b6f2c38ddbe01d2648d719904ac9..0eabb9e01884eed04588b18d0a68c349b21fb137 100644 (file)
@@ -2240,6 +2240,13 @@ void SolverMPGMRES<VectorType>::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 (file)
index 0000000..6e98ca3
--- /dev/null
@@ -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 <deal.II/lac/dynamic_sparsity_pattern.h>
+#include <deal.II/lac/solver_gmres.h>
+#include <deal.II/lac/sparse_matrix.h>
+
+#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<double> 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<double> rhs(M);
+  Vector<double> sol(M);
+  rhs = 1.;
+
+  SolverControl control(M, 1.e-2);
+
+  deallog.get_file_stream() << "FGMRES:\n";
+  SolverFGMRES<Vector<double>> 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<Vector<double>> 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 (file)
index 0000000..7a7440e
--- /dev/null
@@ -0,0 +1,6 @@
+
+FGMRES:
+DEAL::Solver stopped within 1 - 10 iterations
+
+MPGMRES:
+DEAL::Solver stopped within 1 - 10 iterations

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.