From 00f9888947731922488421ec1257241da14dca99 Mon Sep 17 00:00:00 2001 From: Giovanni Alzetta Date: Mon, 13 Nov 2017 12:26:00 +0000 Subject: [PATCH] Added GridTools::guess_point_owner and test --- .../changes/minor/20171113GiovanniAlzetta | 3 + include/deal.II/grid/grid_tools.h | 28 +++ source/grid/grid_tools.cc | 55 +++++ source/grid/grid_tools.inst.in | 8 + tests/grid/grid_tools_guess_pt_owner_1.cc | 192 ++++++++++++++++++ tests/grid/grid_tools_guess_pt_owner_1.output | 46 +++++ 6 files changed, 332 insertions(+) create mode 100644 doc/news/changes/minor/20171113GiovanniAlzetta create mode 100644 tests/grid/grid_tools_guess_pt_owner_1.cc create mode 100644 tests/grid/grid_tools_guess_pt_owner_1.output diff --git a/doc/news/changes/minor/20171113GiovanniAlzetta b/doc/news/changes/minor/20171113GiovanniAlzetta new file mode 100644 index 0000000000..9af0b3aa46 --- /dev/null +++ b/doc/news/changes/minor/20171113GiovanniAlzetta @@ -0,0 +1,3 @@ +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. +
+(Giovanni Alzetta, 2017/11/13) diff --git a/include/deal.II/grid/grid_tools.h b/include/deal.II/grid/grid_tools.h index d0e80399e2..c0fb24880f 100644 --- a/include/deal.II/grid/grid_tools.h +++ b/include/deal.II/grid/grid_tools.h @@ -1288,6 +1288,34 @@ namespace GridTools const std::function &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 unsigned int of the point in @p points + * to the rank of the owner. + * - A map from the index unsigned int of the point in @p points + * to the ranks of the guessed owners. + * + * @author Giovanni Alzetta, 2017 + */ + template + 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 > > + &global_bboxes, + const std::vector< Point > &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 diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc index 6f64aa1560..cc6c40ff50 100644 --- a/source/grid/grid_tools.cc +++ b/source/grid/grid_tools.cc @@ -2339,6 +2339,61 @@ next_cell: } } + + + + template + 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 > > + &global_bboxes, + const std::vector< Point > &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 owners_found; + // Check in which other processes the point might be + for (unsigned int rk=0; rk &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 std::vector::active_cell_iterator> > vertex_to_cell_map(const Triangulation &triangulation) diff --git a/source/grid/grid_tools.inst.in b/source/grid/grid_tools.inst.in index e5bea7e2ba..7552215724 100644 --- a/source/grid/grid_tools.inst.in +++ b/source/grid/grid_tools.inst.in @@ -184,6 +184,14 @@ for (deal_II_space_dimension : SPACE_DIMENSIONS) template std::vector< std::vector< BoundingBox > > GridTools::exchange_local_bounding_boxes(const std::vector< BoundingBox >&, 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 > > &, + const std::vector > &); + } diff --git a/tests/grid/grid_tools_guess_pt_owner_1.cc b/tests/grid/grid_tools_guess_pt_owner_1.cc new file mode 100644 index 0000000000..1372939cf3 --- /dev/null +++ b/tests/grid/grid_tools_guess_pt_owner_1.cc @@ -0,0 +1,192 @@ +// --------------------------------------------------------------------- +// +// 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 +#include +#include + +#include + +template +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 > > + global_bboxes(n_procs); + + unsigned int tot_bbox = 0; + + for (unsigned int rk=0; rk,Point> boundaries; + boundaries.first[0] = tot_bbox; + boundaries.second[0] = tot_bbox + 1 ; + for (int i=1; i 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(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(output_tp)[rk]; + for (unsigned int box=0; box p; + p[0] = tot_bbox*double(Testing::rand())/RAND_MAX; + for (unsigned int d=1; d,Point> 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 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(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); + } + +} diff --git a/tests/grid/grid_tools_guess_pt_owner_1.output b/tests/grid/grid_tools_guess_pt_owner_1.output new file mode 100644 index 0000000000..06872b5511 --- /dev/null +++ b/tests/grid/grid_tools_guess_pt_owner_1.output @@ -0,0 +1,46 @@ + +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 -- 2.39.5