# include <deal.II/base/table.h>
# include <deal.II/lac/block_matrix_base.h>
+# include <deal.II/lac/block_sparsity_pattern.h>
# include <deal.II/lac/petsc_parallel_sparse_matrix.h>
# include <deal.II/lac/petsc_parallel_block_vector.h>
# include <deal.II/lac/exceptions.h>
void reinit (const unsigned int n_block_rows,
const unsigned int n_block_columns);
+
+ /**
+ * Efficiently reinit the block matrix for a parallel computation.
+ * Only the BlockSparsityPattern of the Simple type can efficiently
+ * store large sparsity patterns in parallel, so this is the only
+ * supported argument.
+ * The IndexSets describe the locally owned range of DoFs for each block.
+ * Note that each IndexSet needs to contiguous. For a symmetric structure
+ * hand in the same vector for the first two arguments.
+ */
+ void reinit(const std::vector<IndexSet> &rows,
+ const std::vector<IndexSet> &cols,
+ const BlockCompressedSimpleSparsityPattern &bcsp,
+ const MPI_Comm &com);
+
+
+ /**
+ * Same as above but only for a symmetric structure only.
+ */
+ void reinit(const std::vector<IndexSet> &sizes,
+ const BlockCompressedSimpleSparsityPattern &bcsp,
+ const MPI_Comm &com);
+
+
+
/**
* Matrix-vector multiplication:
* let $dst = M*src$ with $M$
}
}
+ void
+ BlockSparseMatrix::
+ reinit(const std::vector<IndexSet> &rows,
+ const std::vector<IndexSet> &cols,
+ const BlockCompressedSimpleSparsityPattern &bcsp,
+ const MPI_Comm &com)
+ {
+ Assert(rows.size() == bcsp.n_block_rows(), ExcMessage("invalid size"));
+ Assert(cols.size() == bcsp.n_block_cols(), ExcMessage("invalid size"));
+
+
+ clear();
+ this->sub_objects.reinit (bcsp.n_block_rows(),
+ bcsp.n_block_cols());
+
+ std::vector<unsigned int> row_sizes;
+ for (unsigned int r=0; r<this->n_block_rows(); ++r)
+ row_sizes.push_back( bcsp.block(r,0).n_rows() );
+ this->row_block_indices.reinit (row_sizes);
+
+ std::vector<unsigned int> col_sizes;
+ for (unsigned int c=0; c<this->n_block_cols(); ++c)
+ col_sizes.push_back( bcsp.block(0,c).n_cols() );
+ this->column_block_indices.reinit (col_sizes);
+
+ for (unsigned int r=0; r<this->n_block_rows(); ++r)
+ for (unsigned int c=0; c<this->n_block_cols(); ++c)
+ {
+ Assert(rows[r].size() == bcsp.block(r,c).n_rows(), ExcMessage("invalid size"));
+ Assert(cols[r].size() == bcsp.block(r,c).n_cols(), ExcMessage("invalid size"));
+
+ BlockType *p = new BlockType();
+ p->reinit(rows[r],
+ cols[c],
+ bcsp.block(r,c),
+ com);
+ this->sub_objects[r][c] = p;
+ }
+
+ collect_sizes();
+ }
+
+ void
+ BlockSparseMatrix::
+ reinit(const std::vector<IndexSet> &sizes,
+ const BlockCompressedSimpleSparsityPattern &bcsp,
+ const MPI_Comm &com)
+ {
+ reinit(sizes, sizes, bcsp, com);
+ }
+
void