Fix calling TrilinosWrappers::BlockSparseMatrix::reinit with BlockSparseMatrix

diff --git a/include/deal.II/lac/trilinos_block_sparse_matrix.h b/include/deal.II/lac/trilinos_block_sparse_matrix.h
index ef88d17c6a0f1cc7cadfb7a63fa583b238d839e6..3f448489828ad8f83d59eed6a64d99d36027a1f7 100644 (file)
--- a/include/deal.II/lac/trilinos_block_sparse_matrix.h
+++ b/include/deal.II/lac/trilinos_block_sparse_matrix.h
@@ -168,6 +168,19 @@ namespace TrilinosWrappers
     void
     reinit(const BlockSparsityPatternType &block_sparsity_pattern);
 
+    /**
+     * This function initializes the Trilinos matrix using the deal.II sparse
+     * matrix and the entries stored therein. It uses a threshold to copy only
+     * elements whose modulus is larger than the threshold (so zeros in the
+     * deal.II matrix can be filtered away).
+     */
+    void
+    reinit(
+      const std::vector<IndexSet> &              parallel_partitioning,
+      const ::dealii::BlockSparseMatrix<double> &dealii_block_sparse_matrix,
+      const MPI_Comm &                           communicator = MPI_COMM_WORLD,
+      const double                               drop_tolerance = 1e-13);
+
     /**
      * This function initializes the Trilinos matrix using the deal.II sparse
      * matrix and the entries stored therein. It uses a threshold to copy only

diff --git a/source/lac/trilinos_block_sparse_matrix.cc b/source/lac/trilinos_block_sparse_matrix.cc
index c2a02b1cc453ea9a1b4e10b80bf5d445cd7ecd00..70349aaa06495be622dec2a5ed675fc123dcd51a 100644 (file)
--- a/source/lac/trilinos_block_sparse_matrix.cc
+++ b/source/lac/trilinos_block_sparse_matrix.cc
@@ -160,6 +160,43 @@ namespace TrilinosWrappers
 
 
 
+  void
+  BlockSparseMatrix::reinit(
+    const std::vector<IndexSet> &              parallel_partitioning,
+    const ::dealii::BlockSparseMatrix<double> &dealii_block_sparse_matrix,
+    const MPI_Comm &                           communicator,
+    const double                               drop_tolerance)
+  {
+    const size_type n_block_rows = parallel_partitioning.size();
+
+    Assert(n_block_rows == dealii_block_sparse_matrix.n_block_rows(),
+           ExcDimensionMismatch(n_block_rows,
+                                dealii_block_sparse_matrix.n_block_rows()));
+    Assert(n_block_rows == dealii_block_sparse_matrix.n_block_cols(),
+           ExcDimensionMismatch(n_block_rows,
+                                dealii_block_sparse_matrix.n_block_cols()));
+
+    // Call the other basic reinit function ...
+    reinit(n_block_rows, n_block_rows);
+
+    // ... and then assign the correct
+    // data to the blocks.
+    for (size_type r = 0; r < this->n_block_rows(); ++r)
+      for (size_type c = 0; c < this->n_block_cols(); ++c)
+        {
+          this->sub_objects[r][c]->reinit(parallel_partitioning[r],
+                                          parallel_partitioning[c],
+                                          dealii_block_sparse_matrix.block(r,
+                                                                           c),
+                                          communicator,
+                                          drop_tolerance);
+        }
+
+    collect_sizes();
+  }
+
+
+
   void
   BlockSparseMatrix::reinit(
     const ::dealii::BlockSparseMatrix<double> &dealii_block_sparse_matrix,