From aefad5f4e265047f507bf89564b80f2da43560cb Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Wed, 17 Nov 2021 11:57:24 -0500 Subject: [PATCH] Utilities::MPI::create_mpi_data_type_n_bytes part of #12752 --- doc/news/changes/minor/20211117tjhei | 3 + include/deal.II/base/mpi.h | 33 +++++ source/base/mpi.cc | 77 ++++++++++++ tests/mpi/create_mpi_datatype_01.cc | 116 ++++++++++++++++++ .../create_mpi_datatype_01.mpirun=2.output | 21 ++++ 5 files changed, 250 insertions(+) create mode 100644 doc/news/changes/minor/20211117tjhei create mode 100644 tests/mpi/create_mpi_datatype_01.cc create mode 100644 tests/mpi/create_mpi_datatype_01.mpirun=2.output diff --git a/doc/news/changes/minor/20211117tjhei b/doc/news/changes/minor/20211117tjhei new file mode 100644 index 0000000000..162100bc31 --- /dev/null +++ b/doc/news/changes/minor/20211117tjhei @@ -0,0 +1,3 @@ +New: Utilities::MPI::create_mpi_data_type_n_bytes() helps with large MPI communication and IO. +
+(Timo Heister, 2021/11/17) diff --git a/include/deal.II/base/mpi.h b/include/deal.II/base/mpi.h index 33eb2c1dfe..406dde250a 100644 --- a/include/deal.II/base/mpi.h +++ b/include/deal.II/base/mpi.h @@ -498,6 +498,39 @@ namespace Utilities const MPI_Comm &comm); #endif + + /** + * Create an MPI_Datatype that consists of @p n_bytes bytes. + * + * The resulting data type can be used in MPI send/recv or MPI IO to process + * messages of sizes larger than 2 GB with MPI_Byte as the underlying data + * type. This helper is required for MPI versions before 4.0 because + * routines like MPI_Send + * use a signed interger for the @p count variable. Instead, you can use this + * data type with the appropriate size set to the size of your message and + * by passing + * 1 as the @p count. + * + * @note You need to free this data type after you are done using it by calling + * MPI_Type_free(&result). + * + * Usage example: + * + * std::vector buffer; + * if (buffer.size()<(1U<<31)) + * MPI_Send(buffer.data(), buffer.size(), MPI_BYTE, dest, tag, comm); + * else + * { + * MPI_Datatype bigtype; + * Utilities::MPI::create_mpi_data_type_n_bytes(bigtype, buffer.size()); + * MPI_Send(buffer.data(), 1, bigtype, dest, tag, comm); + * MPI_Type_free(&bigtype); + * } + * + */ + void + create_mpi_data_type_n_bytes(MPI_Datatype &result, std::size_t n_bytes); + /** * Return the sum over all processors of the value @p t. This function is * collective over all processors given in the diff --git a/source/base/mpi.cc b/source/base/mpi.cc index 861853eed9..b39d298fce 100644 --- a/source/base/mpi.cc +++ b/source/base/mpi.cc @@ -281,6 +281,8 @@ namespace Utilities return res; } + + IndexSet create_evenly_distributed_partitioning(const MPI_Comm & comm, const IndexSet::size_type total_size) @@ -295,6 +297,81 @@ namespace Utilities + void + create_mpi_data_type_n_bytes(MPI_Datatype &result, std::size_t n_bytes) + { +# ifdef DEAL_II_WITH_MPI + // Simplified version from BigMPI repository, see + // https://github.com/jeffhammond/BigMPI/blob/5300b18cc8ec1b2431bf269ee494054ee7bd9f72/src/type_contiguous_x.c#L74 + // (code is MIT licensed) + + // We create an MPI datatype that has the layout A*nB where A is + // max_signed_int bytes repeated n times and B is the remainder. + + const MPI_Count max_signed_int = (1U << 31) - 1; + + const MPI_Count n_chunks = n_bytes / max_signed_int; + const MPI_Count n_bytes_remainder = n_bytes % max_signed_int; + + Assert(static_cast(max_signed_int * n_chunks + + n_bytes_remainder) == n_bytes, + ExcInternalError()); + + MPI_Datatype chunks; + + int ierr = MPI_Type_vector( + n_chunks, max_signed_int, max_signed_int, MPI_BYTE, &chunks); + AssertThrowMPI(ierr); + + MPI_Datatype remainder; + ierr = MPI_Type_contiguous(n_bytes_remainder, MPI_BYTE, &remainder); + AssertThrowMPI(ierr); + + const int blocklengths[2] = {1, 1}; + const MPI_Aint displacements[2] = {0, + static_cast(n_chunks) * + max_signed_int}; + + // This fails if Aint happens to be 32 bits (maybe on some 32bit + // systems as it has type "long" which is usually 64bits) or the + // message is very, very big. + AssertThrow( + displacements[1] == n_chunks * max_signed_int, + ExcMessage( + "ERROR in create_mpi_data_type_n_bytes(), size too big to support.")); + + const MPI_Datatype types[2] = {chunks, remainder}; + ierr = + MPI_Type_create_struct(2, blocklengths, displacements, types, &result); + AssertThrowMPI(ierr); + + ierr = MPI_Type_commit(&result); + AssertThrowMPI(ierr); + + ierr = MPI_Type_free(&chunks); + AssertThrowMPI(ierr); + ierr = MPI_Type_free(&remainder); + AssertThrowMPI(ierr); + +# ifdef DEBUG +# if DEAL_II_MPI_VERSION_GTE(3, 0) + MPI_Count size64; + // this function is only available starting MPI 3.0: + ierr = MPI_Type_size_x(result, &size64); + AssertThrowMPI(ierr); + + Assert(size64 == static_cast(n_bytes), ExcInternalError()); +# endif +# endif + +# else + (void)result; + (void)n_bytes; +# endif + } + + + std::vector compute_point_to_point_communication_pattern( const MPI_Comm & mpi_comm, diff --git a/tests/mpi/create_mpi_datatype_01.cc b/tests/mpi/create_mpi_datatype_01.cc new file mode 100644 index 0000000000..f7564bda40 --- /dev/null +++ b/tests/mpi/create_mpi_datatype_01.cc @@ -0,0 +1,116 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2021 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + +// Test Utilities::MPI::create_mpi_data_type_n_bytes + +// To not require a lot of memory, we run this test with smaller +// size. Change this bool to test the real 64bit size communication: +const bool run_big = false; + +#include + +#include "../tests.h" + +using namespace dealii; + +void +test_data_type(long long n_bytes) +{ + MPI_Datatype bigtype; + Utilities::MPI::create_mpi_data_type_n_bytes(bigtype, n_bytes); + + deallog << "checking size " << n_bytes << ":"; + + int size32; + int ierr = MPI_Type_size(bigtype, &size32); + AssertThrowMPI(ierr); + + if (size32 == MPI_UNDEFINED) + deallog << " size32=" + << "UNDEFINED (too big)"; + else + deallog << " size32=" << size32; + +#if DEAL_II_MPI_VERSION_GTE(3, 0) + MPI_Count size64; + ierr = MPI_Type_size_x(bigtype, &size64); + AssertThrowMPI(ierr); + + deallog << " size64=" << size64; +#endif + + deallog << std::endl; + + MPI_Type_free(&bigtype); +} + +void +test_send_recv(MPI_Comm comm) +{ + unsigned int myid = Utilities::MPI::this_mpi_process(comm); + const long long n_bytes = (run_big) ? ((1ULL << 31) + 37) : (37ULL); + + if (myid == 0) + { + std::vector buffer(n_bytes, 'A'); + buffer[n_bytes - 1] = 'B'; + MPI_Datatype bigtype; + Utilities::MPI::create_mpi_data_type_n_bytes(bigtype, buffer.size()); + int ierr = + MPI_Send(buffer.data(), 1, bigtype, 1 /* dest */, 0 /* tag */, comm); + AssertThrowMPI(ierr); + ierr = MPI_Type_free(&bigtype); + AssertThrowMPI(ierr); + } + else if (myid == 1) + { + std::vector buffer(n_bytes, '?'); + MPI_Datatype bigtype; + Utilities::MPI::create_mpi_data_type_n_bytes(bigtype, buffer.size()); + int ierr = MPI_Recv(buffer.data(), + 1, + bigtype, + 0 /* src */, + 0 /* tag */, + comm, + MPI_STATUS_IGNORE); + AssertThrowMPI(ierr); + ierr = MPI_Type_free(&bigtype); + AssertThrowMPI(ierr); + + AssertThrow(buffer[0] == 'A', ExcInternalError()); + AssertThrow(buffer[n_bytes - 1] == 'B', ExcInternalError()); + } + deallog << "send_recv: OK" << std::endl; +} + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll all; + + test_data_type(0); + test_data_type(1); + test_data_type(1ULL << 30); + test_data_type(1ULL << 31); + test_data_type(1ULL << 32); + test_data_type(1ULL << 33); + test_data_type(1ULL << 55); + test_data_type(58493729485ULL); + + test_send_recv(MPI_COMM_WORLD); +} diff --git a/tests/mpi/create_mpi_datatype_01.mpirun=2.output b/tests/mpi/create_mpi_datatype_01.mpirun=2.output new file mode 100644 index 0000000000..181257c7d2 --- /dev/null +++ b/tests/mpi/create_mpi_datatype_01.mpirun=2.output @@ -0,0 +1,21 @@ + +DEAL:0::checking size 0: size32=0 size64=0 +DEAL:0::checking size 1: size32=1 size64=1 +DEAL:0::checking size 1073741824: size32=1073741824 size64=1073741824 +DEAL:0::checking size 2147483648: size32=UNDEFINED (too big) size64=2147483648 +DEAL:0::checking size 4294967296: size32=UNDEFINED (too big) size64=4294967296 +DEAL:0::checking size 8589934592: size32=UNDEFINED (too big) size64=8589934592 +DEAL:0::checking size 36028797018963968: size32=UNDEFINED (too big) size64=36028797018963968 +DEAL:0::checking size 58493729485: size32=UNDEFINED (too big) size64=58493729485 +DEAL:0::send_recv: OK + +DEAL:1::checking size 0: size32=0 size64=0 +DEAL:1::checking size 1: size32=1 size64=1 +DEAL:1::checking size 1073741824: size32=1073741824 size64=1073741824 +DEAL:1::checking size 2147483648: size32=UNDEFINED (too big) size64=2147483648 +DEAL:1::checking size 4294967296: size32=UNDEFINED (too big) size64=4294967296 +DEAL:1::checking size 8589934592: size32=UNDEFINED (too big) size64=8589934592 +DEAL:1::checking size 36028797018963968: size32=UNDEFINED (too big) size64=36028797018963968 +DEAL:1::checking size 58493729485: size32=UNDEFINED (too big) size64=58493729485 +DEAL:1::send_recv: OK + -- 2.39.5