--- /dev/null
+Fixed: GridTools::distributed_compute_point_locations now returns the correct
+maps values (third component of the output tuple). Added a test to check it.
+<br>
+(Giovanni Alzetta, 2018/04/23)
+
{
// Point inside locally owned cell: storing all its data
std::vector < Point<spacedim> > cell_points(indices_loc.size());
+ std::vector < unsigned int > cell_points_idx(indices_loc.size());
for (unsigned int i=0; i< indices_loc.size(); ++i)
{
// Adding the point to the cell points
cell_points[i] = local_points[indices_loc[i]];
+
// Storing the index: notice indices loc refer to the local points
// vector, but we need to return the index with respect of
// the points owned by the current process
+ cell_points_idx[i] = local_points_idx[indices_loc[i]];
classified_pts.emplace_back(local_points_idx[indices_loc[i]]);
}
output_unmap.emplace(std::make_pair(cell_loc,
std::make_tuple(q_loc,
- indices_loc,
+ cell_points_idx,
cell_points,
std::vector<unsigned int>
(indices_loc.size(),cell_loc->subdomain_id()))));
// Point inside ghost cell: storing all its information and preparing
// it to be sent
std::vector < Point<spacedim> > cell_points(indices_loc.size());
+ std::vector < unsigned int > cell_points_idx(indices_loc.size());
for (unsigned int i=0; i< indices_loc.size(); ++i)
{
cell_points[i] = local_points[indices_loc[i]];
+ cell_points_idx[i] = local_points_idx[indices_loc[i]];
classified_pts.emplace_back(local_points_idx[indices_loc[i]]);
}
// Each key of the following map represent a process,
// To identify the cell on the other process we use the cell id
std::get<0>(map_tuple_owner).emplace_back(cell_loc->id());
std::get<1>(map_tuple_owner).emplace_back(q_loc);
- std::get<2>(map_tuple_owner).emplace_back(indices_loc);
+ std::get<2>(map_tuple_owner).emplace_back(cell_points_idx);
std::get<3>(map_tuple_owner).emplace_back(cell_points);
}
// else: the cell is artificial, nothing to do
// Rewriting the contents of the map in human readable format
const auto &received_process = rank_and_points.first;
const auto &received_points = rank_and_points.second.first;
- const auto &received_ranks = rank_and_points.second.second;
+ const auto &received_map = rank_and_points.second.second;
// Initializing the vectors needed to store the result of compute point location
std::vector< typename Triangulation<dim, spacedim>::active_cell_iterator > in_cell;
std::vector< Point<spacedim> > cell_points(loc_size);
for (unsigned int pt=0; pt<loc_size; ++pt)
{
- cell_maps[pt] = received_ranks[proc_maps[pt]];
+ cell_maps[pt] = received_map[proc_maps[pt]];
cell_points[pt] = received_points[proc_maps[pt]];
}
in_maps.emplace_back(cell_maps);
{
// Finding/adding in the map the current process
auto ¤t_pts = other_owned_pts[indices.second];
+ // Indices.first is the index of the considered point in local points
current_pts.first.emplace_back(local_points[indices.first]);
current_pts.second.emplace_back(indices.first);
}
// Step 3: construct vectors containing uncertain points i.e. those whose owner
// is known among few guesses
+ // The maps goes from rank of the probable owner to a pair of vectors: the first
+ // containing the points, the second containing the ranks in the current process
std::map<
unsigned int,
std::pair< std::vector < Point<spacedim> >,
- std::vector<unsigned int > > >
+ std::vector< unsigned int > > >
other_check_pts;
- const auto &other_check_idx = std::get<2>(guessed_points);
+ // This map goes from the point index to a vector of
+ // ranks probable owners
+ const std::map< unsigned int, std::vector< unsigned int > >
+ &other_check_idx = std::get<2>(guessed_points);
// Points in classified pts need not to be communicated;
// sorting the array classified pts in order to use
for (const auto &pt_to_guesses: other_check_idx)
{
+ const auto &point_idx = pt_to_guesses.first;
+ const auto &probable_owners_rks = pt_to_guesses.second;
if ( !std::binary_search(
- classified_pts.begin(), classified_pts.end(),pt_to_guesses.first) )
+ classified_pts.begin(), classified_pts.end(), point_idx) )
// The point wasn't found in ghost or locally owned cells: adding it to the map
- for (unsigned int rank=0; rank<pt_to_guesses.second.size(); ++rank)
- if (pt_to_guesses.second[rank] != my_rank)
+ for (unsigned int i=0; i<probable_owners_rks.size(); ++i)
+ if (probable_owners_rks[i] != my_rank)
{
- auto ¤t_pts = other_check_pts[pt_to_guesses.second[rank]];
- current_pts.first.emplace_back(local_points[pt_to_guesses.first]);
- current_pts.second.emplace_back(pt_to_guesses.second[rank]);
+ // add to the data for process probable_owners_rks[i]
+ auto ¤t_pts = other_check_pts[probable_owners_rks[i]];
+ // The point local_points[point_idx]
+ current_pts.first.emplace_back(local_points[point_idx]);
+ // and its index in the current process
+ current_pts.second.emplace_back(point_idx);
}
}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 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 GridTools::compute_point_locations, in particular checks if the
+// maps which are sent/received are correct
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/mpi.h>
+#include <deal.II/base/function_lib.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/grid_tools_cache.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/grid/manifold_lib.h>
+#include <deal.II/grid/filtered_iterator.h>
+#include <deal.II/fe/mapping_q.h>
+
+
+using namespace dealii;
+
+
+
+template <int dim>
+void test_distributed_cpt(unsigned int ref_cube)
+{
+ MPI_Comm mpi_communicator = MPI_COMM_WORLD;
+ unsigned int n_procs = Utilities::MPI::n_mpi_processes(mpi_communicator);
+ unsigned int my_rank = Utilities::MPI::this_mpi_process(mpi_communicator);
+
+ deallog << "Testing for dim = " << dim << " on " << n_procs << " processes" << std::endl;
+ deallog << "Cube refinements: " << ref_cube << std::endl;
+
+ // Creeating the cube on which to run distributed cpt loc
+ parallel::distributed::Triangulation<dim> cube_d(mpi_communicator);
+ GridGenerator::hyper_cube(cube_d);
+ cube_d.refine_global(ref_cube);
+
+ // We shall use the points from a shared grid so that each index is known
+ std::vector<Point< dim > > test_points;
+ Triangulation<dim> cube;
+ GridGenerator::hyper_cube(cube);
+ cube.refine_global(ref_cube);
+ for (auto cell: cube.active_cell_iterators())
+ test_points.emplace_back(cell->center());
+
+ deallog << " Testing on " << test_points.size() << " points" << std::endl;
+
+ // Computing bounding boxes describing the locally owned part of the mesh
+ IteratorFilters::LocallyOwnedCell locally_owned_cell_predicate;
+ std::vector< BoundingBox<dim> > local_bbox = GridTools::compute_mesh_predicate_bounding_box
+ (cube_d, locally_owned_cell_predicate,
+ 1, false, 4);
+
+ // Obtaining the global mesh description through an all to all communication
+ std::vector< std::vector< BoundingBox<dim> > > global_bboxes;
+ global_bboxes = Utilities::MPI::all_gather(mpi_communicator,local_bbox);
+
+ // Initializing the cache
+ GridTools::Cache<dim,dim> cache_d(cube_d);
+ auto output_tuple = distributed_compute_point_locations
+ (cache_d,test_points,global_bboxes);
+ const auto &maps = std::get<2>(output_tuple);
+ const auto &points = std::get<3>(output_tuple);
+ const auto &ranks = std::get<4>(output_tuple);
+
+ // Testing the results: if the map in maps is correct test_points[maps[i][j]] == points[i][j]
+ bool test_passed = true;
+
+ for (unsigned int i=0; i< points.size(); ++i)
+ {
+ for (unsigned int j=0; j< maps[i].size(); ++j)
+ if ( (test_points[maps[i][j]] - points[i][j]).norm() > 1e-10 )
+ {
+ deallog << " Error in cell " << i << " with position " << j << std::endl;
+ deallog << " Received map was: " << maps[i][j] << std::endl;
+ deallog << " From rank: " << ranks[i][j] << std::endl;
+ test_passed = false;
+ }
+ }
+
+ if (test_passed)
+ deallog << "Test passed" << std::endl;
+ else
+ deallog << "Test FAILED" << std::endl;
+}
+
+int main (int argc, char *argv[])
+{
+ Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1);
+ MPILogInitAll log;
+
+ deallog << "Deal.II GridTools::distributed_compute_point_locations" << std::endl;
+ deallog << "Test 3: maps values" << std::endl;
+ deallog << "2D tests:" << std::endl;
+ test_distributed_cpt<2>(3);
+ deallog << "3D tests" << std::endl;
+ test_distributed_cpt<3>(3);
+}
--- /dev/null
+
+DEAL:0::Deal.II GridTools::distributed_compute_point_locations
+DEAL:0::Test 3: maps values
+DEAL:0::2D tests:
+DEAL:0::Testing for dim = 2 on 2 processes
+DEAL:0::Cube refinements: 3
+DEAL:0:: Testing on 64 points
+DEAL:0::Test passed
+DEAL:0::3D tests
+DEAL:0::Testing for dim = 3 on 2 processes
+DEAL:0::Cube refinements: 3
+DEAL:0:: Testing on 512 points
+DEAL:0::Test passed
+
+DEAL:1::Deal.II GridTools::distributed_compute_point_locations
+DEAL:1::Test 3: maps values
+DEAL:1::2D tests:
+DEAL:1::Testing for dim = 2 on 2 processes
+DEAL:1::Cube refinements: 3
+DEAL:1:: Testing on 64 points
+DEAL:1::Test passed
+DEAL:1::3D tests
+DEAL:1::Testing for dim = 3 on 2 processes
+DEAL:1::Cube refinements: 3
+DEAL:1:: Testing on 512 points
+DEAL:1::Test passed
+
--- /dev/null
+
+DEAL:0::Deal.II GridTools::distributed_compute_point_locations
+DEAL:0::Test 3: maps values
+DEAL:0::2D tests:
+DEAL:0::Testing for dim = 2 on 3 processes
+DEAL:0::Cube refinements: 3
+DEAL:0:: Testing on 64 points
+DEAL:0::Test passed
+DEAL:0::3D tests
+DEAL:0::Testing for dim = 3 on 3 processes
+DEAL:0::Cube refinements: 3
+DEAL:0:: Testing on 512 points
+DEAL:0::Test passed
+
+DEAL:1::Deal.II GridTools::distributed_compute_point_locations
+DEAL:1::Test 3: maps values
+DEAL:1::2D tests:
+DEAL:1::Testing for dim = 2 on 3 processes
+DEAL:1::Cube refinements: 3
+DEAL:1:: Testing on 64 points
+DEAL:1::Test passed
+DEAL:1::3D tests
+DEAL:1::Testing for dim = 3 on 3 processes
+DEAL:1::Cube refinements: 3
+DEAL:1:: Testing on 512 points
+DEAL:1::Test passed
+
+
+DEAL:2::Deal.II GridTools::distributed_compute_point_locations
+DEAL:2::Test 3: maps values
+DEAL:2::2D tests:
+DEAL:2::Testing for dim = 2 on 3 processes
+DEAL:2::Cube refinements: 3
+DEAL:2:: Testing on 64 points
+DEAL:2::Test passed
+DEAL:2::3D tests
+DEAL:2::Testing for dim = 3 on 3 processes
+DEAL:2::Cube refinements: 3
+DEAL:2:: Testing on 512 points
+DEAL:2::Test passed
+