From: Christoph Goering Date: Fri, 15 Dec 2017 17:42:00 +0000 (+0100) Subject: added DynamicSparsityPattern::compute_mmult_pattern(..) X-Git-Tag: v9.0.0-rc1~633^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=2b61e4e9a7c0e25a2444db2d919f7e385e8e8e83;p=dealii.git added DynamicSparsityPattern::compute_mmult_pattern(..) --- diff --git a/include/deal.II/lac/dynamic_sparsity_pattern.h b/include/deal.II/lac/dynamic_sparsity_pattern.h index 4d96e1d44a..b2ce716150 100644 --- a/include/deal.II/lac/dynamic_sparsity_pattern.h +++ b/include/deal.II/lac/dynamic_sparsity_pattern.h @@ -420,6 +420,14 @@ public: */ void symmetrize (); + /** + * Construct and store in this object the sparsity pattern corresponding to + * the product of @p left and @p right sparsity pattern. + */ + template + void compute_mmult_pattern(const SparsityPatternTypeLeft &left, + const SparsityPatternTypeRight &right); + /** * Print the sparsity pattern. The output consists of one line per row of * the format [i,j1,j2,j3,...]. i is the row number and diff --git a/source/lac/dynamic_sparsity_pattern.cc b/source/lac/dynamic_sparsity_pattern.cc index 18fce9577e..e41f1c05de 100644 --- a/source/lac/dynamic_sparsity_pattern.cc +++ b/source/lac/dynamic_sparsity_pattern.cc @@ -14,6 +14,7 @@ // --------------------------------------------------------------------- #include +#include #include #include @@ -407,6 +408,37 @@ DynamicSparsityPattern::symmetrize () +template +void +DynamicSparsityPattern::compute_mmult_pattern ( + const SparsityPatternTypeLeft &left, + const SparsityPatternTypeRight &right) +{ + Assert (left.n_cols() == right.n_rows(), + ExcDimensionMismatch(left.n_cols(), right.n_rows())); + + this->reinit(left.n_rows(),right.n_cols()); + + typename SparsityPatternTypeLeft::iterator + it_left = left.begin(), + end_left = left.end(); + for (; it_left!=end_left; ++it_left) + { + const unsigned int j = it_left->column(); + + // We are sitting on entry (i,j) of the left sparsity pattern. We then need + // to add all entries (i,k) to the final sparsity pattern where (j,k) exists + // in the right sparsity pattern -- i.e., we need to iterate over row j. + typename SparsityPatternTypeRight::iterator + it_right = right.begin(j), + end_right = right.end(j); + for (; it_right!=end_right; ++it_right) + this->add(it_left->row(),it_right->column()); + } +} + + + void DynamicSparsityPattern::print (std::ostream &out) const { @@ -519,4 +551,17 @@ add_entries(std::vector::iterator, const bool); #endif +template void DynamicSparsityPattern::compute_mmult_pattern( + const DynamicSparsityPattern &, + const DynamicSparsityPattern &); +template void DynamicSparsityPattern::compute_mmult_pattern( + const DynamicSparsityPattern &, + const SparsityPattern &); +template void DynamicSparsityPattern::compute_mmult_pattern( + const SparsityPattern &, + const DynamicSparsityPattern &); +template void DynamicSparsityPattern::compute_mmult_pattern( + const SparsityPattern &, + const SparsityPattern &); + DEAL_II_NAMESPACE_CLOSE