]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Augment a test to also check sparse vanka. 10396/head
authorDavid Wells <drwells@email.unc.edu>
Wed, 27 May 2020 20:47:15 +0000 (16:47 -0400)
committerDavid Wells <drwells@email.unc.edu>
Wed, 27 May 2020 20:47:15 +0000 (16:47 -0400)
tests/lac/sparse_matrices.cc
tests/lac/sparse_matrices.output

index 73b91514ad52daeacd816c9ff6b20d9dad072dd8..f05f2c8aa6d51f93f476c1b2d987f7621ab8d829 100644 (file)
@@ -25,6 +25,7 @@
 #include <deal.II/lac/sparse_matrix.h>
 #include <deal.II/lac/sparse_matrix_ez.h>
 #include <deal.II/lac/sparse_matrix_ez.templates.h>
+#include <deal.II/lac/sparse_vanka.h>
 #include <deal.II/lac/vector.h>
 
 #include "../tests.h"
@@ -122,6 +123,45 @@ check_vmult_quadratic(std::vector<double> &residuals,
   PREC_CHECK(prich, Tsolve, block_ssor);
   PREC_CHECK(prich, Tsolve, block_sor);
   PREC_CHECK(prich, Tsolve, block_psor);
+
+  // vanka needs to match the type
+  using value_type = typename MatrixType::value_type;
+  if (std::is_same<MatrixType, SparseMatrix<value_type>>::value)
+    {
+      deallog << "Vanka" << std::endl;
+      const SparseMatrix<value_type> &B =
+        reinterpret_cast<const SparseMatrix<value_type> &>(A);
+
+      std::vector<bool> selected_dofs(B.m(), false);
+      // tag every third dof for vanka
+      for (unsigned int i = 0; i < B.m(); i += 3)
+        selected_dofs[i] = true;
+      SparseVanka<value_type> vanka;
+      vanka.initialize(
+        B, typename SparseVanka<value_type>::AdditionalData(selected_dofs));
+
+      SparseBlockVanka<value_type> block_vanka_1(
+        B,
+        selected_dofs,
+        n_blocks,
+        SparseBlockVanka<value_type>::index_intervals);
+      SparseBlockVanka<value_type> block_vanka_2(
+        B, selected_dofs, n_blocks, SparseBlockVanka<value_type>::adaptive);
+
+      PREC_CHECK(prich, solve, vanka);
+      PREC_CHECK(prich, solve, block_vanka_1);
+      PREC_CHECK(prich, solve, block_vanka_2);
+
+      // since SparseMatrixEZ doesn't support this format remove it from the
+      // residuals:
+      deallog << "vanka residual: " << residuals.back() << std::endl;
+      residuals.pop_back();
+      deallog << "vanka residual: " << residuals.back() << std::endl;
+      residuals.pop_back();
+      deallog << "vanka residual: " << residuals.back() << std::endl;
+      residuals.pop_back();
+    }
+
   deallog.pop();
 }
 
index 338d0d7d3a79f7c02de6e80471d7c42e2e2ed828..05f28facd33319c1b80dd47bedd28dbda66ccf68 100644 (file)
@@ -97,6 +97,16 @@ DEAL:5-SparseMatrix<double>:block_sor:RichardsonT::Starting value 4.00
 DEAL:5-SparseMatrix<double>:block_sor:RichardsonT::Failure step 10 value 1.07e-06
 DEAL:5-SparseMatrix<double>:block_psor:RichardsonT::Starting value 4.00
 DEAL:5-SparseMatrix<double>:block_psor:RichardsonT::Failure step 10 value 8.62e-07
+DEAL:5-SparseMatrix<double>::Vanka
+DEAL:5-SparseMatrix<double>:vanka:Richardson::Starting value 4.00
+DEAL:5-SparseMatrix<double>:vanka:Richardson::Failure step 10 value 0.000523
+DEAL:5-SparseMatrix<double>:block_vanka_1:Richardson::Starting value 4.00
+DEAL:5-SparseMatrix<double>:block_vanka_1:Richardson::Failure step 10 value 6.00
+DEAL:5-SparseMatrix<double>:block_vanka_2:Richardson::Starting value 4.00
+DEAL:5-SparseMatrix<double>:block_vanka_2:Richardson::Failure step 10 value 0.0179
+DEAL:5-SparseMatrix<double>::vanka residual: 0.0179
+DEAL:5-SparseMatrix<double>::vanka residual: 6.00
+DEAL:5-SparseMatrix<double>::vanka residual: 0.000523
 DEAL:: 0       0       5.00
 DEAL:: 0       1       -2.00
 DEAL:: 0       4       -1.00
@@ -395,6 +405,16 @@ DEAL:9-SparseMatrix<double>:block_sor:RichardsonT::Starting value 4.00
 DEAL:9-SparseMatrix<double>:block_sor:RichardsonT::Failure step 10 value 2.94e-06
 DEAL:9-SparseMatrix<double>:block_psor:RichardsonT::Starting value 4.00
 DEAL:9-SparseMatrix<double>:block_psor:RichardsonT::Failure step 10 value 1.86e-06
+DEAL:9-SparseMatrix<double>::Vanka
+DEAL:9-SparseMatrix<double>:vanka:Richardson::Starting value 4.00
+DEAL:9-SparseMatrix<double>:vanka:Richardson::Failure step 10 value 1.75e-10
+DEAL:9-SparseMatrix<double>:block_vanka_1:Richardson::Starting value 4.00
+DEAL:9-SparseMatrix<double>:block_vanka_1:Richardson::Failure step 10 value 5.61
+DEAL:9-SparseMatrix<double>:block_vanka_2:Richardson::Starting value 4.00
+DEAL:9-SparseMatrix<double>:block_vanka_2:Richardson::Failure step 10 value 0.00189
+DEAL:9-SparseMatrix<double>::vanka residual: 0.00189
+DEAL:9-SparseMatrix<double>::vanka residual: 5.61
+DEAL:9-SparseMatrix<double>::vanka residual: 1.75e-10
 DEAL::Assemble
 DEAL:9-SparseMatrixEZ<double>:identity:Richardson::Starting value 4.00
 DEAL:9-SparseMatrixEZ<double>:identity:Richardson::Failure step 10 value 2.41

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.