--- /dev/null
+New: Given a vector of points, the GridTools::guess_point_owner() function uses bounding boxes describing the mesh to guess the processes that own these points. A test for the function has been added.
+<br>
+(Giovanni Alzetta, 2017/11/13)
const std::function<bool (const typename MeshType::active_cell_iterator &)> &predicate,
const unsigned int &refinement_level = 0, const bool &allow_merge = false, const unsigned int &max_boxes = numbers::invalid_unsigned_int);
+ /**
+ * Given an array of points, use the global bounding box description obtained using
+ * GridTools::compute_mesh_predicate_bounding_box to guess, for each of them,
+ * which process might own it.
+ *
+ * @param[in] global_bboxes Vector of bounding boxes describing the portion of
+ * mesh with a property for each process.
+ * @param[in] points Array of points to test.
+ *
+ * @param[out] A tuple containing the following information:
+ * - A vector indicized with ranks of processes. For each rank it contains
+ * a vector of the indices of points it might own.
+ * - A map from the index <code>unsigned int</code> of the point in @p points
+ * to the rank of the owner.
+ * - A map from the index <code>unsigned int</code> of the point in @p points
+ * to the ranks of the guessed owners.
+ *
+ * @author Giovanni Alzetta, 2017
+ */
+ template <int spacedim>
+ std::tuple< std::vector< std::vector< unsigned int > >,
+ std::map< unsigned int, unsigned int>,
+ std::map< unsigned int, std::vector< unsigned int > > >
+ guess_point_owner (const std::vector< std::vector< BoundingBox<spacedim> > >
+ &global_bboxes,
+ const std::vector< Point<spacedim> > &points);
+
+
/**
* Return the adjacent cells of all the vertices. If a vertex is also a
* hanging node, the associated coarse cell is also returned. The vertices
}
}
+
+
+
+ template <int spacedim>
+ std::tuple< std::vector< std::vector< unsigned int > >,
+ std::map< unsigned int, unsigned int>,
+ std::map< unsigned int, std::vector< unsigned int > > >
+ guess_point_owner (const std::vector< std::vector< BoundingBox<spacedim> > >
+ &global_bboxes,
+ const std::vector< Point<spacedim> > &points)
+ {
+ unsigned int n_procs = global_bboxes.size();
+ std::vector< std::vector< unsigned int > > point_owners(n_procs);
+ std::map< unsigned int, unsigned int> map_owners_found;
+ std::map< unsigned int, std::vector< unsigned int > > map_owners_guessed;
+
+ unsigned int n_points = points.size();
+ for (unsigned int pt=0; pt<n_points; ++pt)
+ {
+ // Keep track of how many processes we guess to own the point
+ std::vector< unsigned int > owners_found;
+ // Check in which other processes the point might be
+ for (unsigned int rk=0; rk<n_procs; ++rk)
+ {
+ for (const BoundingBox<spacedim> &bbox: global_bboxes[rk])
+ if (bbox.point_inside(points[pt]))
+ {
+ point_owners[rk].emplace_back(pt);
+ owners_found.emplace_back(rk);
+ break; // We can check now the next process
+ }
+ }
+ Assert(owners_found.size() > 0,
+ ExcMessage("No owners found for the point " + std::to_string(pt)));
+ if (owners_found.size()==1)
+ map_owners_found[pt] = owners_found[0];
+ else
+ // Multiple owners
+ map_owners_guessed[pt] = owners_found;
+ }
+
+ std::tuple< std::vector< std::vector< unsigned int > >,
+ std::map< unsigned int, unsigned int>,
+ std::map< unsigned int, std::vector< unsigned int > > >
+ output_tuple;
+
+ std::get<0>(output_tuple) = point_owners;
+ std::get<1>(output_tuple) = map_owners_found;
+ std::get<2>(output_tuple) = map_owners_guessed;
+
+ return output_tuple;
+ }
+
+
+
template <int dim, int spacedim>
std::vector<std::set<typename Triangulation<dim,spacedim>::active_cell_iterator> >
vertex_to_cell_map(const Triangulation<dim,spacedim> &triangulation)
template std::vector< std::vector< BoundingBox<deal_II_space_dimension> > >
GridTools::exchange_local_bounding_boxes(const std::vector< BoundingBox<deal_II_space_dimension> >&,
MPI_Comm);
+
+ template
+ std::tuple< std::vector< std::vector< unsigned int > >,
+ std::map< unsigned int, unsigned int>,
+ std::map< unsigned int, std::vector< unsigned int > > >
+ GridTools::guess_point_owner (const std::vector< std::vector< BoundingBox<deal_II_space_dimension> > > &,
+ const std::vector<Point<deal_II_space_dimension> > &);
+
}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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 for GridTools::guess_point_owner
+
+#include "../tests.h"
+
+#include <deal.II/base/point.h>
+#include <deal.II/base/bounding_box.h>
+#include <deal.II/grid/grid_tools.h>
+
+#include <tuple>
+
+template <int spacedim>
+void test_point_owner(unsigned int n_procs)
+{
+ // Step 1: creating the vector of bounding boxes
+ // We want to simulate n_procs bounding boxes:
+ // each process of rank rk has rk + 1 bounding boxes
+ // the bounding boxes describe a mesh which is a
+ // segment (in dim 1), a rectangle (in dim 2),
+ // a parallelogram (in dim 3)
+ // and using the center point of each bounding box
+
+ std::vector< std::vector< BoundingBox<spacedim> > >
+ global_bboxes(n_procs);
+
+ unsigned int tot_bbox = 0;
+
+ for (unsigned int rk=0; rk<n_procs; ++rk)
+ for (unsigned int box=0; box<rk; ++box)
+ {
+ std::pair<Point<spacedim>,Point<spacedim>> boundaries;
+ boundaries.first[0] = tot_bbox;
+ boundaries.second[0] = tot_bbox + 1 ;
+ for (int i=1; i<spacedim; i++)
+ boundaries.second[i] = 1;
+
+ BoundingBox<spacedim> new_box(boundaries);
+ global_bboxes[rk].push_back(new_box);
+ ++tot_bbox;
+ }
+
+ // Creating the vector of center points
+ std::vector< Point< spacedim > > points(tot_bbox);
+ for (unsigned int i=0; i<tot_bbox; ++i)
+ {
+ points[i][0] = i + 0.5;
+ for (unsigned int j=1; j<spacedim; j++)
+ points[i][j] = 0.5;
+ }
+
+ // Step 2: testing the function
+ auto output_tp = GridTools::guess_point_owner(global_bboxes,points);
+
+ bool test_passed = true;
+
+ if (std::get<2>(output_tp).size()!=0)
+ {
+ test_passed = false;
+ deallog << " Test 1 failed: multiple owners" << std::endl;
+ }
+ else
+ {
+ // Check the points are all in the correct rank
+ unsigned int tot_pt = 0;
+ for (unsigned int rk=0; rk<n_procs; ++rk)
+ {
+ const auto &rk_points = std::get<0>(output_tp)[rk];
+ for (unsigned int box=0; box<rk; ++box)
+ {
+ if ( std::find( rk_points.begin(), rk_points.end(),tot_pt) == rk_points.end())
+ {
+ deallog << "Point " << tot_pt << " not found in rank " << rk << std::endl;
+ test_passed = false;
+ }
+ else
+ ++tot_pt;
+ }
+ }
+ }
+
+ // Step 3: creating random points and adding bounding boxes
+
+ // Creating the random points
+ points.clear();
+
+ unsigned int n_points = 2*tot_bbox;
+ for (size_t i=0; i<n_points; ++i)
+ {
+ Point<spacedim> p;
+ p[0] = tot_bbox*double(Testing::rand())/RAND_MAX;
+ for (unsigned int d=1; d<spacedim; ++d)
+ p[d] = double(Testing::rand())/RAND_MAX;
+ points.push_back(p);
+ }
+
+ // Adding a bounding box to some processes which covers
+ // part of the domain already covered by other processes
+
+ for (unsigned int rk=1; rk<n_procs; rk+=2)
+ {
+ std::pair<Point<spacedim>,Point<spacedim>> boundaries
+ = global_bboxes[rk][0].get_boundary_points();
+
+ boundaries.first[0] -= 0.5*double(Testing::rand())/RAND_MAX;
+ boundaries.second[0] -= 0.5*double(Testing::rand())/RAND_MAX;
+
+ BoundingBox<spacedim> new_box(boundaries);
+ global_bboxes[rk].push_back(new_box);
+ }
+
+ // Running the function again and checking the output
+ output_tp = GridTools::guess_point_owner(global_bboxes,points);
+
+ for (unsigned int rk=0; rk<n_procs ; ++rk)
+ {
+ for (const auto &pt: std::get<0>(output_tp)[rk])
+ {
+ bool found = false;
+ for (const auto &bbox: global_bboxes[rk])
+ if (bbox.point_inside(points[pt]))
+ {
+ found = true;
+ break;
+ }
+
+ if ( ! found )
+ {
+ test_passed = false;
+ deallog << "Couldn't find point " << pt << " in process " << rk << std::endl;
+ }
+ else
+ {
+ // Check if the point is contained also in one of the maps
+ if ( std::get<1>(output_tp).find(pt) == std::get<1>(output_tp).end() &&
+ std::get<2>(output_tp).find(pt) == std::get<2>(output_tp).end())
+ {
+ test_passed = false;
+ deallog << "Output maps missing point " << pt << " in process " << rk << std::endl;
+ }
+ }
+ }
+ }
+
+ if (test_passed)
+ deallog << "TEST PASSED" << std::endl;
+ else
+ deallog << "TEST FAILED" << std::endl;
+}
+
+int main()
+{
+ initlog();
+
+ deallog << "Test for GridTools::guess_point_owner " << std::endl;
+ deallog << std::endl << "Test for dimension 1" << std::endl;
+ for (unsigned int d=1; d<4; ++d)
+ {
+ deallog << "Simulating " << d << " processes" << std::endl;
+ test_point_owner<1>(d);
+ }
+
+ deallog << std::endl << "Test for dimension 2" << std::endl;
+ for (unsigned int d=2; d<10; ++d)
+ {
+ deallog << "Simulating " << d << " processes" << std::endl;
+ test_point_owner<2>(d);
+ }
+
+
+ deallog << std::endl << "Test for dimension 3" << std::endl;
+ for (unsigned int d=1; d<17; d+=2)
+ {
+ deallog << "Simulating " << d << " processes" << std::endl;
+ test_point_owner<3>(d);
+ }
+
+}
--- /dev/null
+
+DEAL::Test for GridTools::guess_point_owner
+DEAL::
+DEAL::Test for dimension 1
+DEAL::Simulating 1 processes
+DEAL::TEST PASSED
+DEAL::Simulating 2 processes
+DEAL::TEST PASSED
+DEAL::Simulating 3 processes
+DEAL::TEST PASSED
+DEAL::
+DEAL::Test for dimension 2
+DEAL::Simulating 2 processes
+DEAL::TEST PASSED
+DEAL::Simulating 3 processes
+DEAL::TEST PASSED
+DEAL::Simulating 4 processes
+DEAL::TEST PASSED
+DEAL::Simulating 5 processes
+DEAL::TEST PASSED
+DEAL::Simulating 6 processes
+DEAL::TEST PASSED
+DEAL::Simulating 7 processes
+DEAL::TEST PASSED
+DEAL::Simulating 8 processes
+DEAL::TEST PASSED
+DEAL::Simulating 9 processes
+DEAL::TEST PASSED
+DEAL::
+DEAL::Test for dimension 3
+DEAL::Simulating 1 processes
+DEAL::TEST PASSED
+DEAL::Simulating 3 processes
+DEAL::TEST PASSED
+DEAL::Simulating 5 processes
+DEAL::TEST PASSED
+DEAL::Simulating 7 processes
+DEAL::TEST PASSED
+DEAL::Simulating 9 processes
+DEAL::TEST PASSED
+DEAL::Simulating 11 processes
+DEAL::TEST PASSED
+DEAL::Simulating 13 processes
+DEAL::TEST PASSED
+DEAL::Simulating 15 processes
+DEAL::TEST PASSED