From: Martin Kronbichler <kronbichler@lnm.mw.tum.de> Date: Thu, 24 May 2018 11:10:10 +0000 (+0200) Subject: Do call plain update_ghost_values()/compress() for block vector with many blocks. X-Git-Tag: v9.1.0-rc1~1101^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=ab12185ea15a5d238dc9941bd8ac13fda91d0e1f;p=dealii.git Do call plain update_ghost_values()/compress() for block vector with many blocks. --- diff --git a/include/deal.II/matrix_free/matrix_free.h b/include/deal.II/matrix_free/matrix_free.h index dc8ce0b827..2a4652aab4 100644 --- a/include/deal.II/matrix_free/matrix_free.h +++ b/include/deal.II/matrix_free/matrix_free.h @@ -2919,130 +2919,155 @@ namespace internal + // if the input vector did not have ghosts imported, clear them here again + // in order to avoid subsequent operations e.g. in linear solvers to work + // with ghosts all the time template <int dim, typename VectorStruct, typename Number> inline - void update_ghost_values_start (const VectorStruct &vec, - VectorDataExchange<dim,Number> &exchanger, - const unsigned int channel = 0) + void reset_ghost_values (const VectorStruct &vec, + VectorDataExchange<dim,Number> &exchanger) { - update_ghost_values_start_block(vec, channel, - std::integral_constant<bool, - IsBlockVector<VectorStruct>::value>(), - exchanger); + reset_ghost_values_block(vec, + std::integral_constant<bool, + IsBlockVector<VectorStruct>::value>(), + exchanger); } template <int dim, typename Number, typename Number2> inline - void update_ghost_values_start (const LinearAlgebra::distributed::Vector<Number> &vec, - VectorDataExchange<dim,Number2> &exchanger, - const unsigned int channel = 0) + void reset_ghost_values (const LinearAlgebra::distributed::Vector<Number> &vec, + VectorDataExchange<dim,Number2> &exchanger) { - exchanger.update_ghost_values_start(channel, vec); + exchanger.reset_ghost_values(vec); } template <int dim, typename VectorStruct, typename Number> inline - void update_ghost_values_start (const std::vector<VectorStruct> &vec, - VectorDataExchange<dim,Number> &exchanger) + void reset_ghost_values (const std::vector<VectorStruct> &vec, + VectorDataExchange<dim,Number> &exchanger) { - unsigned int component_index = 0; for (unsigned int comp=0; comp<vec.size(); comp++) - { - update_ghost_values_start(vec[comp], exchanger, component_index); - component_index += n_components(vec[comp]); - } + reset_ghost_values(vec[comp], exchanger); } template <int dim, typename VectorStruct, typename Number> inline - void update_ghost_values_start (const std::vector<VectorStruct *> &vec, - VectorDataExchange<dim,Number> &exchanger) + void reset_ghost_values (const std::vector<VectorStruct *> &vec, + VectorDataExchange<dim,Number> &exchanger) { - unsigned int component_index = 0; for (unsigned int comp=0; comp<vec.size(); comp++) - { - update_ghost_values_start(*vec[comp], exchanger, component_index); - component_index += n_components(*vec[comp]); - } + reset_ghost_values(*vec[comp], exchanger); } template <int dim, typename VectorStruct, typename Number> inline - void update_ghost_values_start_block (const VectorStruct &vec, - const unsigned int channel, - std::integral_constant<bool,true>, - VectorDataExchange<dim,Number> &exchanger) + void reset_ghost_values_block (const VectorStruct &vec, + std::integral_constant<bool,true>, + VectorDataExchange<dim,Number> &exchanger) { for (unsigned int i=0; i<vec.n_blocks(); ++i) - update_ghost_values_start(vec.block(i), exchanger, channel+i); + reset_ghost_values(vec.block(i), exchanger); + } + + + + // A helper function to identify block vectors with many components where we + // should not try to overlap computations and communcation because there + // would be too many outstanding communcation requests. This is the base case + // for generic vectors + template <typename VectorStruct> + constexpr unsigned int get_communication_block_size(const VectorStruct &) + { + return numbers::invalid_unsigned_int; + } + + + + // Specialized case for the block vector which as the additional member + // variable + template <typename Number> + constexpr unsigned int get_communication_block_size(const LinearAlgebra::distributed::BlockVector<Number> &) + { + return LinearAlgebra::distributed::BlockVector<Number>::communication_block_size; } - // if the input vector did not have ghosts imported, clear them here again - // in order to avoid subsequent operations e.g. in linear solvers to work - // with ghosts all the time template <int dim, typename VectorStruct, typename Number> inline - void reset_ghost_values (const VectorStruct &vec, - VectorDataExchange<dim,Number> &exchanger) + void update_ghost_values_start (const VectorStruct &vec, + VectorDataExchange<dim,Number> &exchanger, + const unsigned int channel = 0) { - reset_ghost_values_block(vec, - std::integral_constant<bool, - IsBlockVector<VectorStruct>::value>(), - exchanger); + update_ghost_values_start_block(vec, channel, + std::integral_constant<bool, + IsBlockVector<VectorStruct>::value>(), + exchanger); } template <int dim, typename Number, typename Number2> inline - void reset_ghost_values (const LinearAlgebra::distributed::Vector<Number> &vec, - VectorDataExchange<dim,Number2> &exchanger) + void update_ghost_values_start (const LinearAlgebra::distributed::Vector<Number> &vec, + VectorDataExchange<dim,Number2> &exchanger, + const unsigned int channel = 0) { - exchanger.reset_ghost_values(vec); + exchanger.update_ghost_values_start(channel, vec); } template <int dim, typename VectorStruct, typename Number> inline - void reset_ghost_values (const std::vector<VectorStruct> &vec, - VectorDataExchange<dim,Number> &exchanger) + void update_ghost_values_start (const std::vector<VectorStruct> &vec, + VectorDataExchange<dim,Number> &exchanger) { + unsigned int component_index = 0; for (unsigned int comp=0; comp<vec.size(); comp++) - reset_ghost_values(vec[comp], exchanger); + { + update_ghost_values_start(vec[comp], exchanger, component_index); + component_index += n_components(vec[comp]); + } } template <int dim, typename VectorStruct, typename Number> inline - void reset_ghost_values (const std::vector<VectorStruct *> &vec, - VectorDataExchange<dim,Number> &exchanger) + void update_ghost_values_start (const std::vector<VectorStruct *> &vec, + VectorDataExchange<dim,Number> &exchanger) { + unsigned int component_index = 0; for (unsigned int comp=0; comp<vec.size(); comp++) - reset_ghost_values(*vec[comp], exchanger); + { + update_ghost_values_start(*vec[comp], exchanger, component_index); + component_index += n_components(*vec[comp]); + } } template <int dim, typename VectorStruct, typename Number> inline - void reset_ghost_values_block (const VectorStruct &vec, - std::integral_constant<bool,true>, - VectorDataExchange<dim,Number> &exchanger) + void update_ghost_values_start_block (const VectorStruct &vec, + const unsigned int channel, + std::integral_constant<bool,true>, + VectorDataExchange<dim,Number> &exchanger) { - for (unsigned int i=0; i<vec.n_blocks(); ++i) - reset_ghost_values(vec.block(i), exchanger); + if (get_communication_block_size(vec) < vec.n_blocks()) + vec.update_ghost_values(); + else + for (unsigned int i=0; i<vec.n_blocks(); ++i) + update_ghost_values_start(vec.block(i), exchanger, channel+i); } @@ -3109,8 +3134,14 @@ namespace internal std::integral_constant<bool,true>, VectorDataExchange<dim,Number> &exchanger) { - for (unsigned int i=0; i<vec.n_blocks(); ++i) - update_ghost_values_finish(vec.block(i), exchanger, channel+i); + if (get_communication_block_size(vec) < vec.n_blocks()) + { + // do nothing, everything has already been completed in the _start() + // call + } + else + for (unsigned int i=0; i<vec.n_blocks(); ++i) + update_ghost_values_finish(vec.block(i), exchanger, channel+i); } @@ -3177,8 +3208,11 @@ namespace internal std::integral_constant<bool,true>, VectorDataExchange<dim, Number> &exchanger) { - for (unsigned int i=0; i<vec.n_blocks(); ++i) - compress_start(vec.block(i), exchanger, channel+i); + if (get_communication_block_size(vec) < vec.n_blocks()) + vec.compress(dealii::VectorOperation::add); + else + for (unsigned int i=0; i<vec.n_blocks(); ++i) + compress_start(vec.block(i), exchanger, channel+i); } @@ -3245,8 +3279,14 @@ namespace internal std::integral_constant<bool,true>, VectorDataExchange<dim, Number> &exchanger) { - for (unsigned int i=0; i<vec.n_blocks(); ++i) - compress_finish(vec.block(i), exchanger, channel+i); + if (get_communication_block_size(vec) < vec.n_blocks()) + { + // do nothing, everything has already been completed in the _start() + // call + } + else + for (unsigned int i=0; i<vec.n_blocks(); ++i) + compress_finish(vec.block(i), exchanger, channel+i); }