From: Peter Munch Date: Tue, 30 Aug 2022 07:33:09 +0000 (+0200) Subject: Enable SolverIDR for block vectors X-Git-Tag: v9.5.0-rc1~996^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F14232%2Fhead;p=dealii.git Enable SolverIDR for block vectors --- diff --git a/include/deal.II/lac/solver_idr.h b/include/deal.II/lac/solver_idr.h index 64bc3c6793..3c6c66daba 100644 --- a/include/deal.II/lac/solver_idr.h +++ b/include/deal.II/lac/solver_idr.h @@ -24,6 +24,7 @@ #include #include +#include #include #include #include @@ -226,6 +227,53 @@ namespace internal } return *data[i]; } + + + + template ::value, VectorType> + * = nullptr> + unsigned int + n_blocks(const VectorType &) + { + return 1; + } + + + + template ::value, VectorType> * = + nullptr> + unsigned int + n_blocks(const VectorType &vector) + { + return vector.n_blocks(); + } + + + + template ::value, VectorType> + * = nullptr> + VectorType & + block(VectorType &vector, const unsigned int b) + { + AssertDimension(b, 0); + (void)b; + return vector; + } + + + + template ::value, VectorType> * = + nullptr> + typename VectorType::BlockType & + block(VectorType &vector, const unsigned int b) + { + return vector.block(b); + } + } // namespace SolverIDRImplementation } // namespace internal @@ -328,8 +376,13 @@ SolverIDR::solve(const MatrixType & A, VectorType &tmp_q = Q(i, x); if (i != 0) { - for (auto indx : tmp_q.locally_owned_elements()) - tmp_q(indx) = normal_distribution(rng); + for (unsigned int b = 0; + b < internal::SolverIDRImplementation::n_blocks(tmp_q); + ++b) + for (auto indx : internal::SolverIDRImplementation::block(tmp_q, b) + .locally_owned_elements()) + internal::SolverIDRImplementation::block(tmp_q, b)(indx) = + normal_distribution(rng); tmp_q.compress(VectorOperation::insert); } else