From: Peter Munch <peterrmuench@gmail.com>
Date: Wed, 24 Mar 2021 10:30:16 +0000 (+0100)
Subject: MatrixFreeTools::compute_matrix() for PETScWrappers::MPI::SparseMatrix
X-Git-Tag: v9.3.0-rc1~279^2
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F11958%2Fhead;p=dealii.git

MatrixFreeTools::compute_matrix() for PETScWrappers::MPI::SparseMatrix
---

diff --git a/include/deal.II/matrix_free/tools.h b/include/deal.II/matrix_free/tools.h
index d3498fba29..7040ffffd8 100644
--- a/include/deal.II/matrix_free/tools.h
+++ b/include/deal.II/matrix_free/tools.h
@@ -826,7 +826,8 @@ namespace MatrixFreeTools
 
             for (unsigned int v = 0; v < n_filled_lanes; ++v)
               {
-                const auto cell_v = matrix_free.get_cell_iterator(cell, v);
+                const auto cell_v =
+                  matrix_free.get_cell_iterator(cell, v, dof_no);
 
                 if (matrix_free.get_mg_level() != numbers::invalid_unsigned_int)
                   cell_v->get_mg_dof_indices(dof_indices);
@@ -837,7 +838,6 @@ namespace MatrixFreeTools
                   dof_indices_mf[j] = dof_indices[lexicographic_numbering[j]];
 
                 constraints.distribute_local_to_global(matrices[v],
-                                                       dof_indices_mf,
                                                        dof_indices_mf,
                                                        dst);
               }
@@ -847,13 +847,6 @@ namespace MatrixFreeTools
       matrix);
 
     matrix.compress(VectorOperation::add);
-
-    for (auto &entry : matrix)
-      if (entry.row() == entry.column() && entry.value() == 0.0)
-        {
-          Assert(constraints.is_constrained(entry.row()), ExcNotImplemented());
-          entry.value() = 1.0;
-        }
   }
 
   template <typename CLASS,
diff --git a/tests/matrix_free/compute_diagonal_util.h b/tests/matrix_free/compute_diagonal_util.h
index 6d2a3c5ddd..6486dc4dc1 100644
--- a/tests/matrix_free/compute_diagonal_util.h
+++ b/tests/matrix_free/compute_diagonal_util.h
@@ -188,6 +188,10 @@ public:
 
         for (; a1 != A1.end(); ++a1, ++a2, ++a_ref)
           {
+            if (a1->row() == a1->column() &&
+                constraints.is_constrained(a1->row()))
+              continue;
+
             Assert(std::abs(a1->value() - a_ref->value()) < 1e-6,
                    ExcNotImplemented());
             Assert(std::abs(a2->value() - a_ref->value()) < 1e-6,