]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added function parallel::GridTools::exchange_local_bounding_boxes and tests 5358/head
authorGiovanni Alzetta <giovannialzetta@hotmail.it>
Fri, 27 Oct 2017 12:45:10 +0000 (12:45 +0000)
committerGiovanni Alzetta <giovannialzetta@hotmail.it>
Wed, 8 Nov 2017 10:47:39 +0000 (10:47 +0000)
doc/news/changes/minor/20171030GiovanniAlzetta [new file with mode: 0644]
include/deal.II/distributed/grid_tools.h
source/distributed/CMakeLists.txt
source/distributed/grid_tools.cc [new file with mode: 0644]
source/distributed/grid_tools.inst.in [new file with mode: 0644]
tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.cc [new file with mode: 0644]
tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.mpirun=12.output [new file with mode: 0644]
tests/distributed_grids/grid_tools_exchange_bounding_boxes_1.mpirun=4.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20171030GiovanniAlzetta b/doc/news/changes/minor/20171030GiovanniAlzetta
new file mode 100644 (file)
index 0000000..1869b1e
--- /dev/null
@@ -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.
+<br>
+(Giovanni Alzetta, 30/10/2017)
index f8eb3d8bf8680629e73f50a00e0fd55dd3ad93fb..47606000f7998f46737216fbcc4baa4d715fdd11 100644 (file)
@@ -17,6 +17,7 @@
 #define dealii_distributed_grid_tools_h
 
 
+#include <deal.II/base/bounding_box.h>
 #include <deal.II/base/config.h>
 #include <deal.II/base/exceptions.h>
 #include <deal.II/distributed/tria_base.h>
@@ -135,6 +136,22 @@ namespace parallel
                                   const std::function<boost::optional<DataType> (const typename MeshType::active_cell_iterator &)> &pack,
                                   const std::function<void (const typename MeshType::active_cell_iterator &, const DataType &)> &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<int spacedim>
+    std::vector< std::vector< BoundingBox<spacedim> > >
+    exchange_local_bounding_boxes(const std::vector< BoundingBox<spacedim> > &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
index 56b6974b564aa80d4187f71135eaa3c97ae5b400..c821f46aec8eff51426165df14cec44c35e4408e 100644 (file)
@@ -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 (file)
index 0000000..0ab6fdb
--- /dev/null
@@ -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 <deal.II/base/point.h>
+#include <deal.II/base/bounding_box.h>
+#include <deal.II/base/mpi.h>
+#include <deal.II/distributed/grid_tools.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+#ifdef DEAL_II_WITH_MPI
+
+namespace parallel
+{
+  namespace GridTools
+  {
+    template<int spacedim>
+    std::vector< std::vector< BoundingBox<spacedim> > >
+    exchange_local_bounding_boxes(const std::vector< BoundingBox<spacedim> > &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<double> loc_data_array(n_local_data);
+      for (unsigned int i=0; i<n_bboxes; ++i)
+        for (unsigned int d=0; d < spacedim; ++d)
+          {
+            // Extracting the coordinates of each boundary point
+            loc_data_array[2*i*spacedim + d] = local_bboxes[i].get_boundary_points().first[d];
+            loc_data_array[2*i*spacedim + spacedim + d] = local_bboxes[i].get_boundary_points().second[d];
+          }
+
+      // Step 2: exchanging the size of local data
+      unsigned int n_procs = Utilities::MPI::n_mpi_processes(mpi_communicator);
+
+      // Vector to store the size of loc_data_array for every process
+      std::vector<int> 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<int> 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<double> 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<spacedim> > > 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<n_bbox_i; ++bbox)
+            {
+              Point<spacedim> p1,p2; // boundary points for bbox
+              for (unsigned int d=0; d<spacedim; ++d)
+                {
+                  p1[d] = data_array[ begin_idx + 2*bbox*spacedim + d];
+                  p2[d] = data_array[ begin_idx + 2*bbox*spacedim + spacedim + d];
+                }
+              BoundingBox<spacedim> 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 (file)
index 0000000..505e37c
--- /dev/null
@@ -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<deal_II_space_dimension> > >
+    exchange_local_bounding_boxes(const std::vector< BoundingBox<deal_II_space_dimension> >&,
+                                  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 (file)
index 0000000..141be4c
--- /dev/null
@@ -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 <deal.II/base/point.h>
+#include <deal.II/base/bounding_box.h>
+#include <deal.II/base/mpi.h>
+#include <deal.II/distributed/grid_tools.h>
+
+template <int spacedim>
+void test_exchange_bbox()
+{
+  // For process i the number of boxes n_bboxes[i%7] is created
+  std::vector<unsigned int> 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<spacedim> > loc_bboxes(n_bboxes[proc%7]);
+  for (unsigned int b = 0; b < loc_bboxes.size(); ++b)
+    {
+      Point<spacedim> pt1,pt2;
+      for (unsigned int d=0; d<spacedim; ++d)
+        {
+          pt1[d] = (proc + 2*b*spacedim + d)*(proc+1);
+          pt2[d] = (proc + (2*b+1)*spacedim + d)*(proc+1);
+        }
+      BoundingBox<spacedim> bbox(std::make_pair(pt1,pt2));
+      loc_bboxes[b] = bbox;
+    }
+
+  std::vector< std::vector< BoundingBox<spacedim> > > 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<n_procs; ++p)
+        {
+          for (unsigned int b = 0; b < n_bboxes[p%7]; ++b)
+            {
+              std::vector<Point<spacedim>> 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<spacedim; ++d)
+                for (unsigned int pt=0; pt<2; ++pt)
+                  if ( std::abs( boundary_pt[pt][d] - (p + (2*b+pt)*spacedim + d)*(p+1)) > 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 (file)
index 0000000..2278be3
--- /dev/null
@@ -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 (file)
index 0000000..6a35125
--- /dev/null
@@ -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
+

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.