From: Wolfgang Bangerth Date: Tue, 15 Feb 2022 22:45:51 +0000 (-0700) Subject: Use a template variable rather than a dummy function to infer MPI types. X-Git-Tag: v9.4.0-rc1~413^2~7 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=1c91300ad7d484ad3c1c2fd7b5752eaf74e60f92;p=dealii.git Use a template variable rather than a dummy function to infer MPI types. --- diff --git a/include/deal.II/base/mpi.h b/include/deal.II/base/mpi.h index 9762169624..4c9e5ccfec 100644 --- a/include/deal.II/base/mpi.h +++ b/include/deal.II/base/mpi.h @@ -1303,154 +1303,180 @@ namespace Utilities /* --------------------------- inline functions ------------------------- */ - /** - * Given a pointer to an object of class T, return the matching - * `MPI_Datatype` to be used for MPI communication. - * - * As an example, passing an `int*` to this function returns `MPI_INT`. - * - * @note In reality, these functions are not template functions templated - * on the parameter T, but free standing inline function overloads. This - * templated version only exists so that it shows up in the documentation. - * The `=delete` statement at the end of the declaration ensures that the - * compiler will never choose this general template and instead look - * for one of the overloads. - */ - template - inline MPI_Datatype - mpi_type_id(const T *) = delete; + namespace internal + { + namespace MPIDataTypes + { + /** + * Given a pointer to an object of class T, return the matching + * `MPI_Datatype` to be used for MPI communication. + * + * As an example, passing an `int*` to this function returns `MPI_INT`. + * + * @note In reality, these functions are not template functions templated + * on the parameter T, but free standing inline function overloads. This + * templated version only exists so that it shows up in the + * documentation. The `=delete` statement at the end of the declaration + * ensures that the compiler will never choose this general template and + * instead look for one of the overloads. + */ + template + inline MPI_Datatype + mpi_type_id(const T *) = delete; #ifndef DOXYGEN # ifdef DEAL_II_WITH_MPI - inline MPI_Datatype - mpi_type_id(const bool *) - { - return MPI_CXX_BOOL; - } + inline MPI_Datatype + mpi_type_id(const bool *) + { + return MPI_CXX_BOOL; + } - inline MPI_Datatype - mpi_type_id(const char *) - { - return MPI_CHAR; - } + inline MPI_Datatype + mpi_type_id(const char *) + { + return MPI_CHAR; + } - inline MPI_Datatype - mpi_type_id(const signed char *) - { - return MPI_SIGNED_CHAR; - } + inline MPI_Datatype + mpi_type_id(const signed char *) + { + return MPI_SIGNED_CHAR; + } - inline MPI_Datatype - mpi_type_id(const short *) - { - return MPI_SHORT; - } + inline MPI_Datatype + mpi_type_id(const short *) + { + return MPI_SHORT; + } - inline MPI_Datatype - mpi_type_id(const int *) - { - return MPI_INT; - } + inline MPI_Datatype + mpi_type_id(const int *) + { + return MPI_INT; + } - inline MPI_Datatype - mpi_type_id(const long int *) - { - return MPI_LONG; - } + inline MPI_Datatype + mpi_type_id(const long int *) + { + return MPI_LONG; + } - inline MPI_Datatype - mpi_type_id(const unsigned char *) - { - return MPI_UNSIGNED_CHAR; - } + inline MPI_Datatype + mpi_type_id(const unsigned char *) + { + return MPI_UNSIGNED_CHAR; + } - inline MPI_Datatype - mpi_type_id(const unsigned short *) - { - return MPI_UNSIGNED_SHORT; - } + inline MPI_Datatype + mpi_type_id(const unsigned short *) + { + return MPI_UNSIGNED_SHORT; + } - inline MPI_Datatype - mpi_type_id(const unsigned int *) - { - return MPI_UNSIGNED; - } + inline MPI_Datatype + mpi_type_id(const unsigned int *) + { + return MPI_UNSIGNED; + } - inline MPI_Datatype - mpi_type_id(const unsigned long int *) - { - return MPI_UNSIGNED_LONG; - } + inline MPI_Datatype + mpi_type_id(const unsigned long int *) + { + return MPI_UNSIGNED_LONG; + } - inline MPI_Datatype - mpi_type_id(const unsigned long long int *) - { - return MPI_UNSIGNED_LONG_LONG; - } + inline MPI_Datatype + mpi_type_id(const unsigned long long int *) + { + return MPI_UNSIGNED_LONG_LONG; + } - inline MPI_Datatype - mpi_type_id(const float *) - { - return MPI_FLOAT; - } + inline MPI_Datatype + mpi_type_id(const float *) + { + return MPI_FLOAT; + } - inline MPI_Datatype - mpi_type_id(const double *) - { - return MPI_DOUBLE; - } + inline MPI_Datatype + mpi_type_id(const double *) + { + return MPI_DOUBLE; + } - inline MPI_Datatype - mpi_type_id(const long double *) - { - return MPI_LONG_DOUBLE; - } + inline MPI_Datatype + mpi_type_id(const long double *) + { + return MPI_LONG_DOUBLE; + } - inline MPI_Datatype - mpi_type_id(const std::complex *) - { - return MPI_COMPLEX; - } + inline MPI_Datatype + mpi_type_id(const std::complex *) + { + return MPI_COMPLEX; + } - inline MPI_Datatype - mpi_type_id(const std::complex *) - { - return MPI_DOUBLE_COMPLEX; - } + inline MPI_Datatype + mpi_type_id(const std::complex *) + { + return MPI_DOUBLE_COMPLEX; + } # endif +#endif + } // namespace MPIDataTypes + } // namespace internal + + + + /** + * A template variable that translates from the data type given as + * template argument to the corresponding + * `MPI_Datatype` to be used for MPI communication. + * + * As an example, the value of `mpi_type_id` is `MPI_INT`. A + * common way to use this variable is when sending an object `obj` + * via MPI functions to another process, and using + * `mpi_type_id` to infer the correct MPI type to + * use for the communication. + */ + template + const MPI_Datatype mpi_type_id = internal::MPIDataTypes::mpi_type_id( + static_cast>> *>(nullptr)); +#ifndef DOXYGEN namespace internal { // declaration for an internal function that lives in mpi.templates.h @@ -1808,7 +1834,7 @@ namespace Utilities const int ierr = MPI_Bcast(buffer + total_sent_count, current_count, - mpi_type_id(buffer), + mpi_type_id, root, comm); AssertThrowMPI(ierr); @@ -1846,8 +1872,11 @@ namespace Utilities } // Exchange the size of buffer - int ierr = MPI_Bcast( - &buffer_size, 1, mpi_type_id(&buffer_size), root_process, comm); + int ierr = MPI_Bcast(&buffer_size, + 1, + mpi_type_id, + root_process, + comm); AssertThrowMPI(ierr); // If not on the root process, correctly size the buffer to diff --git a/include/deal.II/base/mpi.templates.h b/include/deal.II/base/mpi.templates.h index ab5b455752..0af2709a1c 100644 --- a/include/deal.II/base/mpi.templates.h +++ b/include/deal.II/base/mpi.templates.h @@ -85,7 +85,7 @@ namespace Utilities MPI_IN_PLACE, static_cast(output.data()), static_cast(values.size()), - mpi_type_id(values.data()), + mpi_type_id, mpi_op, mpi_communicator); AssertThrowMPI(ierr); @@ -125,7 +125,7 @@ namespace Utilities MPI_IN_PLACE, static_cast(output.data()), static_cast(values.size() * 2), - mpi_type_id(static_cast(nullptr)), + mpi_type_id, mpi_op, mpi_communicator); AssertThrowMPI(ierr); diff --git a/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h b/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h index cb4b6f286c..e342968988 100644 --- a/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h +++ b/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h @@ -106,14 +106,13 @@ namespace Utilities AssertIndexRange(recv_ranks.size(), recv_ptr.size()); for (types::global_dof_index i = 0; i < recv_ranks.size(); ++i) { - const int ierr = - MPI_Irecv(buffers.data() + recv_ptr[i], - recv_ptr[i + 1] - recv_ptr[i], - Utilities::MPI::mpi_type_id(buffers.data()), - recv_ranks[i], - tag, - communicator, - &requests[i + send_ranks.size()]); + const int ierr = MPI_Irecv(buffers.data() + recv_ptr[i], + recv_ptr[i + 1] - recv_ptr[i], + Utilities::MPI::mpi_type_id, + recv_ranks[i], + tag, + communicator, + &requests[i + send_ranks.size()]); AssertThrowMPI(ierr); } @@ -135,14 +134,13 @@ namespace Utilities (send_ptr[i] == buffers.size() && send_ptr[i + 1] == send_ptr[i]), ExcMessage("The input buffer doesn't contain enough entries")); - const int ierr = - MPI_Isend(buffers.data() + send_ptr[i], - send_ptr[i + 1] - send_ptr[i], - Utilities::MPI::mpi_type_id(buffers.data()), - send_ranks[i], - tag, - communicator, - &requests[i]); + const int ierr = MPI_Isend(buffers.data() + send_ptr[i], + send_ptr[i + 1] - send_ptr[i], + Utilities::MPI::mpi_type_id, + send_ranks[i], + tag, + communicator, + &requests[i]); AssertThrowMPI(ierr); } #endif diff --git a/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h b/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h index 72b4aa717e..742f8bc05a 100644 --- a/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h +++ b/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h @@ -571,7 +571,7 @@ namespace internal const auto ierr_1 = MPI_Isend( buffer.data(), buffer.size(), - Utilities::MPI::mpi_type_id(buffer.data()), + Utilities::MPI::mpi_type_id, i.first, Utilities::MPI::internal::Tags::fine_dof_handler_view_reinit, communicator, @@ -638,10 +638,10 @@ namespace internal std::vector buffer; int message_length; - const int ierr_2 = - MPI_Get_count(&status, - Utilities::MPI::mpi_type_id(buffer.data()), - &message_length); + const int ierr_2 = MPI_Get_count( + &status, + Utilities::MPI::mpi_type_id, + &message_length); AssertThrowMPI(ierr_2); buffer.resize(message_length); @@ -649,7 +649,7 @@ namespace internal const int ierr_3 = MPI_Recv( buffer.data(), buffer.size(), - Utilities::MPI::mpi_type_id(buffer.data()), + Utilities::MPI::mpi_type_id, status.MPI_SOURCE, Utilities::MPI::internal::Tags::fine_dof_handler_view_reinit, communicator, diff --git a/source/base/mpi.cc b/source/base/mpi.cc index e392471d4d..3a94760614 100644 --- a/source/base/mpi.cc +++ b/source/base/mpi.cc @@ -87,6 +87,26 @@ namespace Utilities namespace MPI { +#ifdef DEAL_II_WITH_MPI + // Provide definitions of template variables for all valid instantiations. + template const MPI_Datatype mpi_type_id; + template const MPI_Datatype mpi_type_id; + template const MPI_Datatype mpi_type_id; + template const MPI_Datatype mpi_type_id; + template const MPI_Datatype mpi_type_id; + template const MPI_Datatype mpi_type_id; + template const MPI_Datatype mpi_type_id; + template const MPI_Datatype mpi_type_id; + template const MPI_Datatype mpi_type_id; + template const MPI_Datatype mpi_type_id; + template const MPI_Datatype mpi_type_id; + template const MPI_Datatype mpi_type_id; + template const MPI_Datatype mpi_type_id; + template const MPI_Datatype mpi_type_id>; + template const MPI_Datatype mpi_type_id>; +#endif + + MinMaxAvg min_max_avg(const double my_value, const MPI_Comm &mpi_communicator) { diff --git a/source/base/partitioner.cc b/source/base/partitioner.cc index f1b87c5d7e..e6b12f3a2e 100644 --- a/source/base/partitioner.cc +++ b/source/base/partitioner.cc @@ -76,12 +76,13 @@ namespace Utilities types::global_dof_index prefix_sum = 0; #ifdef DEAL_II_WITH_MPI - const int ierr = MPI_Exscan(&local_size, - &prefix_sum, - 1, - Utilities::MPI::mpi_type_id(&prefix_sum), - MPI_SUM, - communicator); + const int ierr = + MPI_Exscan(&local_size, + &prefix_sum, + 1, + Utilities::MPI::mpi_type_id, + MPI_SUM, + communicator); AssertThrowMPI(ierr); #endif diff --git a/source/base/process_grid.cc b/source/base/process_grid.cc index 330df8f706..d7e97a616b 100644 --- a/source/base/process_grid.cc +++ b/source/base/process_grid.cc @@ -247,11 +247,12 @@ namespace Utilities Assert(count > 0, ExcInternalError()); if (mpi_communicator_inactive_with_root != MPI_COMM_NULL) { - const int ierr = MPI_Bcast(value, - count, - Utilities::MPI::mpi_type_id(value), - 0 /*from root*/, - mpi_communicator_inactive_with_root); + const int ierr = + MPI_Bcast(value, + count, + Utilities::MPI::mpi_type_id, + 0 /*from root*/, + mpi_communicator_inactive_with_root); AssertThrowMPI(ierr); } } diff --git a/source/distributed/repartitioning_policy_tools.cc b/source/distributed/repartitioning_policy_tools.cc index e4a7e0c8d4..10b56b5dc7 100644 --- a/source/distributed/repartitioning_policy_tools.cc +++ b/source/distributed/repartitioning_policy_tools.cc @@ -89,8 +89,8 @@ namespace RepartitioningPolicyTools const int ierr = MPI_Exscan(&process_has_active_locally_owned_cells, &offset, 1, - Utilities::MPI::mpi_type_id( - &process_has_active_locally_owned_cells), + Utilities::MPI::mpi_type_id, MPI_SUM, comm); AssertThrowMPI(ierr); @@ -308,12 +308,13 @@ namespace RepartitioningPolicyTools // determine partial sum of weights of this process uint64_t process_local_weight_offset = 0; - int ierr = MPI_Exscan(&process_local_weight, - &process_local_weight_offset, - 1, - Utilities::MPI::mpi_type_id(&process_local_weight), - MPI_SUM, - tria->get_communicator()); + int ierr = + MPI_Exscan(&process_local_weight, + &process_local_weight_offset, + 1, + Utilities::MPI::mpi_type_id, + MPI_SUM, + tria->get_communicator()); AssertThrowMPI(ierr); // total weight of all processes @@ -321,7 +322,7 @@ namespace RepartitioningPolicyTools ierr = MPI_Bcast(&total_weight, 1, - Utilities::MPI::mpi_type_id(&total_weight), + Utilities::MPI::mpi_type_id, n_subdomains - 1, mpi_communicator); AssertThrowMPI(ierr); diff --git a/source/distributed/tria_base.cc b/source/distributed/tria_base.cc index 18736bd4f2..9cbaeb6895 100644 --- a/source/distributed/tria_base.cc +++ b/source/distributed/tria_base.cc @@ -427,7 +427,7 @@ namespace parallel MPI_Exscan(&n_locally_owned_cells, &cell_index, 1, - Utilities::MPI::mpi_type_id(&n_locally_owned_cells), + Utilities::MPI::mpi_type_id, MPI_SUM, this->mpi_communicator); AssertThrowMPI(ierr); @@ -493,13 +493,13 @@ namespace parallel std::vector cell_index( this->n_global_levels(), 0); - int ierr = - MPI_Exscan(n_locally_owned_cells.data(), - cell_index.data(), - this->n_global_levels(), - Utilities::MPI::mpi_type_id(n_locally_owned_cells.data()), - MPI_SUM, - this->mpi_communicator); + int ierr = MPI_Exscan( + n_locally_owned_cells.data(), + cell_index.data(), + this->n_global_levels(), + Utilities::MPI::mpi_type_id, + MPI_SUM, + this->mpi_communicator); AssertThrowMPI(ierr); // 3) determine global number of "active" cells on each level @@ -509,11 +509,12 @@ namespace parallel for (unsigned int l = 0; l < this->n_global_levels(); ++l) n_cells_level[l] = n_locally_owned_cells[l] + cell_index[l]; - ierr = MPI_Bcast(n_cells_level.data(), - this->n_global_levels(), - Utilities::MPI::mpi_type_id(n_cells_level.data()), - this->n_subdomains - 1, - this->mpi_communicator); + ierr = MPI_Bcast( + n_cells_level.data(), + this->n_global_levels(), + Utilities::MPI::mpi_type_id, + this->n_subdomains - 1, + this->mpi_communicator); AssertThrowMPI(ierr); // 4) give global indices to locally-owned cells on level and mark diff --git a/source/lac/scalapack.cc b/source/lac/scalapack.cc index f0b53bc05b..666056d90d 100644 --- a/source/lac/scalapack.cc +++ b/source/lac/scalapack.cc @@ -184,14 +184,14 @@ ScaLAPACKMatrix::ScaLAPACKMatrix( } int ierr = MPI_Bcast(&n_rows, 1, - Utilities::MPI::mpi_type_id(&n_rows), + Utilities::MPI::mpi_type_id, 0 /*from root*/, process_grid->mpi_communicator); AssertThrowMPI(ierr); ierr = MPI_Bcast(&n_columns, 1, - Utilities::MPI::mpi_type_id(&n_columns), + Utilities::MPI::mpi_type_id, 0 /*from root*/, process_grid->mpi_communicator); AssertThrowMPI(ierr); diff --git a/source/matrix_free/vector_data_exchange.cc b/source/matrix_free/vector_data_exchange.cc index afb8ab83bb..dd6a1d4618 100644 --- a/source/matrix_free/vector_data_exchange.cc +++ b/source/matrix_free/vector_data_exchange.cc @@ -896,7 +896,7 @@ namespace internal const int ierr = MPI_Irecv(buffer.data() + ghost_targets_data[i][1] + offset, ghost_targets_data[i][2], - Utilities::MPI::mpi_type_id(buffer.data()), + Utilities::MPI::mpi_type_id, ghost_targets_data[i][0], communication_channel + 1, comm, @@ -917,16 +917,15 @@ namespace internal data_this[import_indices_data.second[j].first + l]; // send data away - const int ierr = - MPI_Isend(temporary_storage.data() + import_targets_data[i][1], - import_targets_data[i][2], - Utilities::MPI::mpi_type_id(data_this.data()), - import_targets_data[i][0], - communication_channel + 1, - comm, - requests.data() + sm_import_ranks.size() + - sm_ghost_ranks.size() + ghost_targets_data.size() + - i); + const int ierr = MPI_Isend( + temporary_storage.data() + import_targets_data[i][1], + import_targets_data[i][2], + Utilities::MPI::mpi_type_id, + import_targets_data[i][0], + communication_channel + 1, + comm, + requests.data() + sm_import_ranks.size() + sm_ghost_ranks.size() + + ghost_targets_data.size() + i); AssertThrowMPI(ierr); } #endif @@ -1173,7 +1172,7 @@ namespace internal const int ierr = MPI_Isend(buffer.data() + ghost_targets_data[i][1], ghost_targets_data[i][2], - Utilities::MPI::mpi_type_id(buffer.data()), + Utilities::MPI::mpi_type_id, ghost_targets_data[i][0], communication_channel + 0, comm, @@ -1184,16 +1183,15 @@ namespace internal for (unsigned int i = 0; i < import_targets_data.size(); ++i) { - const int ierr = - MPI_Irecv(temporary_storage.data() + import_targets_data[i][1], - import_targets_data[i][2], - Utilities::MPI::mpi_type_id(temporary_storage.data()), - import_targets_data[i][0], - communication_channel + 0, - comm, - requests.data() + sm_ghost_ranks.size() + - sm_import_ranks.size() + ghost_targets_data.size() + - i); + const int ierr = MPI_Irecv( + temporary_storage.data() + import_targets_data[i][1], + import_targets_data[i][2], + Utilities::MPI::mpi_type_id, + import_targets_data[i][0], + communication_channel + 0, + comm, + requests.data() + sm_ghost_ranks.size() + sm_import_ranks.size() + + ghost_targets_data.size() + i); AssertThrowMPI(ierr); } #endif