From 2289a4059e4924401001440ae07b6c0a7692c29f Mon Sep 17 00:00:00 2001 From: Giovanni Alzetta Date: Fri, 27 Oct 2017 12:45:10 +0000 Subject: [PATCH] Added function parallel::GridTools::exchange_local_bounding_boxes and tests --- .../changes/minor/20171030GiovanniAlzetta | 4 + include/deal.II/distributed/grid_tools.h | 17 +++ source/distributed/CMakeLists.txt | 2 + source/distributed/grid_tools.cc | 111 ++++++++++++++ source/distributed/grid_tools.inst.in | 26 ++++ .../grid_tools_exchange_bounding_boxes_1.cc | 129 ++++++++++++++++ ...exchange_bounding_boxes_1.mpirun=12.output | 143 ++++++++++++++++++ ..._exchange_bounding_boxes_1.mpirun=4.output | 47 ++++++ 8 files changed, 479 insertions(+) create mode 100644 doc/news/changes/minor/20171030GiovanniAlzetta create mode 100644 source/distributed/grid_tools.cc create mode 100644 source/distributed/grid_tools.inst.in create mode 100644 tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.cc create mode 100644 tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.mpirun=12.output create mode 100644 tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.mpirun=4.output diff --git a/doc/news/changes/minor/20171030GiovanniAlzetta b/doc/news/changes/minor/20171030GiovanniAlzetta new file mode 100644 index 0000000000..1869b1ec99 --- /dev/null +++ b/doc/news/changes/minor/20171030GiovanniAlzetta @@ -0,0 +1,4 @@ +New: Added parallel::GridTools::exchange_local_bounding_boxes function and tests to exchange +vectors of bounding boxes and gather them in every process. +
+(Giovanni Alzetta, 30/10/2017) diff --git a/include/deal.II/distributed/grid_tools.h b/include/deal.II/distributed/grid_tools.h index f8eb3d8bf8..47606000f7 100644 --- a/include/deal.II/distributed/grid_tools.h +++ b/include/deal.II/distributed/grid_tools.h @@ -17,6 +17,7 @@ #define dealii_distributed_grid_tools_h +#include #include #include #include @@ -135,6 +136,22 @@ namespace parallel const std::function (const typename MeshType::active_cell_iterator &)> &pack, const std::function &unpack); + /** + * Exchange with all processors of the MPI communicator @p mpi_communicator the vector of bounding + * boxes @p local_bboxes. + * + * This function is meant to exchange bounding boxes describing the locally owned + * cells in a distributed triangulation obtained with the function + * GridTools::compute_mesh_predicate_bounding_box . + * + * The output vector's size is the number of processes of the MPI communicator: + * its i-th entry contains the vector @p local_bboxes of the i-th process. + */ + template + std::vector< std::vector< BoundingBox > > + exchange_local_bounding_boxes(const std::vector< BoundingBox > &local_bboxes, + MPI_Comm mpi_communicator); + /** * A structure that allows the transfer of cell data of type @p T from one processor * to another. It corresponds to a packed buffer that stores a vector of diff --git a/source/distributed/CMakeLists.txt b/source/distributed/CMakeLists.txt index 56b6974b56..c821f46aec 100644 --- a/source/distributed/CMakeLists.txt +++ b/source/distributed/CMakeLists.txt @@ -16,6 +16,7 @@ INCLUDE_DIRECTORIES(BEFORE ${CMAKE_CURRENT_BINARY_DIR}) SET(_unity_include_src + grid_tools.cc grid_refinement.cc solution_transfer.cc tria.cc @@ -37,6 +38,7 @@ SETUP_SOURCE_LIST("${_unity_include_src}" ) SET(_inst + grid_tools.inst.in grid_refinement.inst.in solution_transfer.inst.in tria.inst.in diff --git a/source/distributed/grid_tools.cc b/source/distributed/grid_tools.cc new file mode 100644 index 0000000000..0ab6fdbd21 --- /dev/null +++ b/source/distributed/grid_tools.cc @@ -0,0 +1,111 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#include +#include +#include +#include + +DEAL_II_NAMESPACE_OPEN + +#ifdef DEAL_II_WITH_MPI + +namespace parallel +{ + namespace GridTools + { + template + std::vector< std::vector< BoundingBox > > + exchange_local_bounding_boxes(const std::vector< BoundingBox > &local_bboxes, + MPI_Comm mpi_communicator) + { +#ifndef DEAL_II_WITH_MPI + (void)local_bboxes; + (void)mpi_communicator; + Assert(false, ExcMessage("parallel::GridTools::exchange_local_bounding_boxes() requires MPI.")); +#else + // Step 1: preparing data to be sent + unsigned int n_bboxes = local_bboxes.size(); + // Dimension of the array to be exchanged (number of double) + int n_local_data = 2*spacedim*n_bboxes; + // data array stores each entry of each point describing the bounding boxes + std::vector loc_data_array(n_local_data); + for (unsigned int i=0; i size_all_data(n_procs); + + // Exchanging the number of bboxes + MPI_Allgather(&n_local_data, 1, MPI_INT, + &(size_all_data[0]), 1, MPI_INT, + mpi_communicator); + + // Now computing the the displacement, relative to recvbuf, + // at which to store the incoming data + std::vector rdispls(n_procs); + rdispls[0] = 0; + for (unsigned int i=1; i < n_procs; ++i) + rdispls[i] = rdispls[i-1] + size_all_data[i-1]; + + // Step 3: exchange the data and bounding boxes: + // Allocating a vector to contain all the received data + std::vector data_array(rdispls.back() + size_all_data.back()); + + MPI_Allgatherv(&(loc_data_array[0]), n_local_data, MPI_DOUBLE, + &(data_array[0]), &(size_all_data[0]), + &(rdispls[0]), MPI_DOUBLE, mpi_communicator); + + // Step 4: create the array of bboxes for output + std::vector< std::vector< BoundingBox > > global_bboxes(n_procs); + unsigned int begin_idx = 0; + for (unsigned int i=0; i < n_procs; ++i) + { + // Number of local bounding boxes + unsigned int n_bbox_i = size_all_data[i]/(spacedim*2); + global_bboxes[i].resize(n_bbox_i); + for (unsigned int bbox=0; bbox p1,p2; // boundary points for bbox + for (unsigned int d=0; d loc_bbox(std::make_pair(p1,p2)); + global_bboxes[i][bbox] = loc_bbox; + } + // Shifting the first index to the start of the next vector + begin_idx += size_all_data[i]; + } + return global_bboxes; +#endif // DEAL_II_WITH_MPI + } + } +} + +// explicit instantiations +#include "grid_tools.inst" + +#endif // DEAL_II_WITH_MPI +DEAL_II_NAMESPACE_CLOSE diff --git a/source/distributed/grid_tools.inst.in b/source/distributed/grid_tools.inst.in new file mode 100644 index 0000000000..505e37cd5f --- /dev/null +++ b/source/distributed/grid_tools.inst.in @@ -0,0 +1,26 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +for (deal_II_space_dimension : SPACE_DIMENSIONS) +{ + namespace parallel \{ + namespace GridTools \{ + template + std::vector< std::vector< BoundingBox > > + exchange_local_bounding_boxes(const std::vector< BoundingBox >&, + MPI_Comm); + \} + \} +} diff --git a/tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.cc b/tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.cc new file mode 100644 index 0000000000..141be4c297 --- /dev/null +++ b/tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.cc @@ -0,0 +1,129 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// test parallel::GridTools::exchange_local_bounding_boxes + +#include "../tests.h" +#include +#include +#include +#include + +template +void test_exchange_bbox() +{ + // For process i the number of boxes n_bboxes[i%7] is created + std::vector n_bboxes = {2,4,3,5,1,3,8}; + + const MPI_Comm &mpi_communicator = MPI_COMM_WORLD; + unsigned int n_procs = Utilities::MPI::n_mpi_processes(mpi_communicator); + unsigned int proc = Utilities::MPI::this_mpi_process(mpi_communicator); + + deallog << "Test for: dimension " << spacedim << std::endl; + deallog << n_procs << " mpi processes" << std::endl; + + // Creating the local bboxes + // The first bbox first coordinate is proc: then use + // the natural numbers increasingly + std::vector< BoundingBox > loc_bboxes(n_bboxes[proc%7]); + for (unsigned int b = 0; b < loc_bboxes.size(); ++b) + { + Point pt1,pt2; + for (unsigned int d=0; d bbox(std::make_pair(pt1,pt2)); + loc_bboxes[b] = bbox; + } + + std::vector< std::vector< BoundingBox > > global_boxes + = parallel::GridTools::exchange_local_bounding_boxes(loc_bboxes,mpi_communicator); + + + bool passed = true; + // First check if the dimensions are correct: + for (unsigned int i=0; i< n_procs; ++i) + { + if ( global_boxes[i].size() != n_bboxes[i%7]) + { + deallog << "Test FAILED: dimension check failed for process " << i << std::endl; + deallog << "Expected number of bboxes: " << n_bboxes[i%7] << std::endl; + deallog << "Received number of bboxes: " << global_boxes[i].size() << std::endl; + passed = false; + } + } + + // Checking if all received data is correct + if (passed) + { + for (unsigned int p=0; p> boundary_pt(2); + boundary_pt[0] = global_boxes[p][b].get_boundary_points().first; + boundary_pt[1] = global_boxes[p][b].get_boundary_points().second; + for (unsigned int d=0; d 1e-10 ) + { + passed = false; + deallog << "Test FAILED, value not corresponding for:" << std::endl; + deallog << "Current proc: " << proc << std::endl; + deallog << "Value for proc: " << p << std::endl; + deallog << "BBox " << b << " Point " << pt << " Dimension " << d << std::endl; + deallog << "Value found: " << (int) boundary_pt[pt][d] << std::endl; + deallog << "Value expected: " << (p + (2*b+pt)*spacedim + d)*(p+1) << std::endl; + } + } + } + } + if (passed) + deallog << "Test passed" << std::endl; + else + { + deallog << "Test failed" << std::endl; + deallog << "Current proc values" << std::endl; + for (auto b: loc_bboxes) + { + deallog << b.get_boundary_points().first << " and " << std::endl; + deallog << b.get_boundary_points().second << std::endl; + } + deallog << "Received values" << std::endl; + for (auto vec_box: global_boxes) + { + for (auto b: vec_box) + { + deallog << b.get_boundary_points().first << " and " << std::endl; + deallog << b.get_boundary_points().second << std::endl; + } + } + } +} + +int main (int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1); + MPILogInitAll log; + + deallog << "Test: parallel::GridTools::exchange_local_bounding_boxes " << std::endl; + + test_exchange_bbox<1>(); + test_exchange_bbox<2>(); + test_exchange_bbox<3>(); +} diff --git a/tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.mpirun=12.output b/tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.mpirun=12.output new file mode 100644 index 0000000000..2278be30b2 --- /dev/null +++ b/tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.mpirun=12.output @@ -0,0 +1,143 @@ + +DEAL:0::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:0::Test for: dimension 1 +DEAL:0::12 mpi processes +DEAL:0::Test passed +DEAL:0::Test for: dimension 2 +DEAL:0::12 mpi processes +DEAL:0::Test passed +DEAL:0::Test for: dimension 3 +DEAL:0::12 mpi processes +DEAL:0::Test passed + +DEAL:1::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:1::Test for: dimension 1 +DEAL:1::12 mpi processes +DEAL:1::Test passed +DEAL:1::Test for: dimension 2 +DEAL:1::12 mpi processes +DEAL:1::Test passed +DEAL:1::Test for: dimension 3 +DEAL:1::12 mpi processes +DEAL:1::Test passed + + +DEAL:2::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:2::Test for: dimension 1 +DEAL:2::12 mpi processes +DEAL:2::Test passed +DEAL:2::Test for: dimension 2 +DEAL:2::12 mpi processes +DEAL:2::Test passed +DEAL:2::Test for: dimension 3 +DEAL:2::12 mpi processes +DEAL:2::Test passed + + +DEAL:3::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:3::Test for: dimension 1 +DEAL:3::12 mpi processes +DEAL:3::Test passed +DEAL:3::Test for: dimension 2 +DEAL:3::12 mpi processes +DEAL:3::Test passed +DEAL:3::Test for: dimension 3 +DEAL:3::12 mpi processes +DEAL:3::Test passed + + +DEAL:4::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:4::Test for: dimension 1 +DEAL:4::12 mpi processes +DEAL:4::Test passed +DEAL:4::Test for: dimension 2 +DEAL:4::12 mpi processes +DEAL:4::Test passed +DEAL:4::Test for: dimension 3 +DEAL:4::12 mpi processes +DEAL:4::Test passed + + +DEAL:5::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:5::Test for: dimension 1 +DEAL:5::12 mpi processes +DEAL:5::Test passed +DEAL:5::Test for: dimension 2 +DEAL:5::12 mpi processes +DEAL:5::Test passed +DEAL:5::Test for: dimension 3 +DEAL:5::12 mpi processes +DEAL:5::Test passed + + +DEAL:6::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:6::Test for: dimension 1 +DEAL:6::12 mpi processes +DEAL:6::Test passed +DEAL:6::Test for: dimension 2 +DEAL:6::12 mpi processes +DEAL:6::Test passed +DEAL:6::Test for: dimension 3 +DEAL:6::12 mpi processes +DEAL:6::Test passed + + +DEAL:7::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:7::Test for: dimension 1 +DEAL:7::12 mpi processes +DEAL:7::Test passed +DEAL:7::Test for: dimension 2 +DEAL:7::12 mpi processes +DEAL:7::Test passed +DEAL:7::Test for: dimension 3 +DEAL:7::12 mpi processes +DEAL:7::Test passed + + +DEAL:8::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:8::Test for: dimension 1 +DEAL:8::12 mpi processes +DEAL:8::Test passed +DEAL:8::Test for: dimension 2 +DEAL:8::12 mpi processes +DEAL:8::Test passed +DEAL:8::Test for: dimension 3 +DEAL:8::12 mpi processes +DEAL:8::Test passed + + +DEAL:9::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:9::Test for: dimension 1 +DEAL:9::12 mpi processes +DEAL:9::Test passed +DEAL:9::Test for: dimension 2 +DEAL:9::12 mpi processes +DEAL:9::Test passed +DEAL:9::Test for: dimension 3 +DEAL:9::12 mpi processes +DEAL:9::Test passed + + +DEAL:10::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:10::Test for: dimension 1 +DEAL:10::12 mpi processes +DEAL:10::Test passed +DEAL:10::Test for: dimension 2 +DEAL:10::12 mpi processes +DEAL:10::Test passed +DEAL:10::Test for: dimension 3 +DEAL:10::12 mpi processes +DEAL:10::Test passed + + +DEAL:11::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:11::Test for: dimension 1 +DEAL:11::12 mpi processes +DEAL:11::Test passed +DEAL:11::Test for: dimension 2 +DEAL:11::12 mpi processes +DEAL:11::Test passed +DEAL:11::Test for: dimension 3 +DEAL:11::12 mpi processes +DEAL:11::Test passed + diff --git a/tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.mpirun=4.output b/tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.mpirun=4.output new file mode 100644 index 0000000000..6a35125161 --- /dev/null +++ b/tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.mpirun=4.output @@ -0,0 +1,47 @@ + +DEAL:0::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:0::Test for: dimension 1 +DEAL:0::4 mpi processes +DEAL:0::Test passed +DEAL:0::Test for: dimension 2 +DEAL:0::4 mpi processes +DEAL:0::Test passed +DEAL:0::Test for: dimension 3 +DEAL:0::4 mpi processes +DEAL:0::Test passed + +DEAL:1::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:1::Test for: dimension 1 +DEAL:1::4 mpi processes +DEAL:1::Test passed +DEAL:1::Test for: dimension 2 +DEAL:1::4 mpi processes +DEAL:1::Test passed +DEAL:1::Test for: dimension 3 +DEAL:1::4 mpi processes +DEAL:1::Test passed + + +DEAL:2::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:2::Test for: dimension 1 +DEAL:2::4 mpi processes +DEAL:2::Test passed +DEAL:2::Test for: dimension 2 +DEAL:2::4 mpi processes +DEAL:2::Test passed +DEAL:2::Test for: dimension 3 +DEAL:2::4 mpi processes +DEAL:2::Test passed + + +DEAL:3::Test: parallel::GridTools::exchange_local_bounding_boxes +DEAL:3::Test for: dimension 1 +DEAL:3::4 mpi processes +DEAL:3::Test passed +DEAL:3::Test for: dimension 2 +DEAL:3::4 mpi processes +DEAL:3::Test passed +DEAL:3::Test for: dimension 3 +DEAL:3::4 mpi processes +DEAL:3::Test passed + -- 2.39.5