]> https://gitweb.dealii.org/ - dealii.git/commitdiff
added test for MPGMRES
authorWyatt Smith <wsmith97@tamu.edu>
Tue, 1 Apr 2025 19:30:17 +0000 (14:30 -0500)
committerMatthias Maier <tamiko@43-1.org>
Mon, 5 May 2025 21:10:57 +0000 (16:10 -0500)
tests/lac/solver_mpgmres_01.cc [new file with mode: 0644]
tests/lac/solver_mpgmres_01.output [new file with mode: 0644]

diff --git a/tests/lac/solver_mpgmres_01.cc b/tests/lac/solver_mpgmres_01.cc
new file mode 100644 (file)
index 0000000..28cafb3
--- /dev/null
@@ -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 <deal.II/lac/dynamic_sparsity_pattern.h>
+#include <deal.II/lac/precondition.h>
+#include <deal.II/lac/solver_gmres.h>
+#include <deal.II/lac/sparse_matrix.h>
+
+#include <iostream>
+
+#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 <typename Preconditioner>
+class MyPreconditioner
+{
+public:
+  MyPreconditioner(const std::string    &name,
+                   const Preconditioner &preconditioner)
+    : name(name)
+    , preconditioner(preconditioner)
+    , print_counter(0)
+  {}
+
+  void
+  vmult(Vector<double> &dst, const Vector<double> &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<double> 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<SparseMatrix<double>> prec_a;
+  prec_a.initialize(matrix);
+  MyPreconditioner wrapper_a("A", prec_a);
+
+  PreconditionJacobi<SparseMatrix<double>> prec_b;
+  prec_b.initialize(matrix);
+  MyPreconditioner wrapper_b("B", prec_b);
+
+  PreconditionSOR<SparseMatrix<double>> prec_c;
+  prec_c.initialize(matrix);
+  MyPreconditioner wrapper_c("C", prec_c);
+
+  Vector<double> rhs(M);
+  Vector<double> sol(M);
+  rhs = 1.;
+
+  SolverControl control(M, 1e-8);
+
+  deallog.get_file_stream() << "FGMRES:\n";
+  SolverFGMRES<Vector<double>> 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<Vector<double>>::AdditionalData data;
+  data.use_truncated_mpgmres_strategy = true;
+  SolverMPGMRES<Vector<double>> 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<Vector<double>> 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 (file)
index 0000000..c50f13e
--- /dev/null
@@ -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

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.