From 8d357ba88e5fb96c267a1162d0fce3e7542e5775 Mon Sep 17 00:00:00 2001 From: Bruno Turcksin Date: Wed, 30 Sep 2015 12:16:21 -0500 Subject: [PATCH] Add vertex to cells map. --- include/deal.II/grid/grid_tools.h | 10 ++ source/grid/grid_tools.cc | 31 +++++++ source/grid/grid_tools.inst.in | 4 + tests/grid/grid_tools_vertex_to_cell_map.cc | 93 +++++++++++++++++++ ...d_tools_vertex_to_cell_map.mpirun=2.output | 5 + 5 files changed, 143 insertions(+) create mode 100644 tests/grid/grid_tools_vertex_to_cell_map.cc create mode 100644 tests/grid/grid_tools_vertex_to_cell_map.mpirun=2.output diff --git a/include/deal.II/grid/grid_tools.h b/include/deal.II/grid/grid_tools.h index a3155cf56d..5daf06dd66 100644 --- a/include/deal.II/grid/grid_tools.h +++ b/include/deal.II/grid/grid_tools.h @@ -26,6 +26,7 @@ #include #include +#include DEAL_II_NAMESPACE_OPEN @@ -630,6 +631,15 @@ namespace GridTools compute_ghost_cell_halo_layer (const Container &container); + /** + * Return the adjacent cells of all the vertices. The vertices are ordered by + * their ID. + */ + template + std::vector::active_cell_iterator> > + vertex_to_cell_map(const Triangulation &triangulation); + + /*@}*/ /** * @name Partitions and subdomains of triangulations diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc index dce6a7f6c7..828d0ccd92 100644 --- a/source/grid/grid_tools.cc +++ b/source/grid/grid_tools.cc @@ -1594,6 +1594,37 @@ next_cell: + template + std::vector::active_cell_iterator> > + vertex_to_cell_map(const Triangulation &triangulation) + { + std::vector::active_cell_iterator> > + vertex_to_cell_map(triangulation.n_vertices()); + typename Triangulation::active_cell_iterator cell = triangulation.begin_active(), + endc = triangulation.end(); + for (; cell!=endc; ++cell) + for (unsigned int i=0; i::vertices_per_cell; ++i) + vertex_to_cell_map[cell->vertex_index(i)].insert(cell); + + // Take care of hanging nodes + cell = triangulation.begin_active(); + for (; cell!=endc; ++cell) + for (unsigned int i=0; i::faces_per_cell; ++i) + { + if ((cell->neighbor_level(i)>=0) && (cell->neighbor_level(i)level())) + { + typename Triangulation::active_cell_iterator adjacent_cell = + cell->neighbor(i); + for (unsigned int j=0; j::vertices_per_face; ++j) + vertex_to_cell_map[cell->face(i)->vertex_index(j)].insert(adjacent_cell); + } + } + + return vertex_to_cell_map; + } + + + template void get_face_connectivity_of_cells (const Triangulation &triangulation, diff --git a/source/grid/grid_tools.inst.in b/source/grid/grid_tools.inst.in index e6a85a1827..5e63a31814 100644 --- a/source/grid/grid_tools.inst.in +++ b/source/grid/grid_tools.inst.in @@ -234,6 +234,10 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS std::map > get_all_vertices_at_boundary (const Triangulation &tria); + template + std::vector::active_cell_iterator> > + vertex_to_cell_map(const Triangulation &triangulation); + #if deal_II_dimension == deal_II_space_dimension # if deal_II_dimension > 1 template diff --git a/tests/grid/grid_tools_vertex_to_cell_map.cc b/tests/grid/grid_tools_vertex_to_cell_map.cc new file mode 100644 index 0000000000..3e66bf9902 --- /dev/null +++ b/tests/grid/grid_tools_vertex_to_cell_map.cc @@ -0,0 +1,93 @@ +// +// Copyright (C) 2015 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 "../tests.h" +#include +#include +#include +#include +#include +#include + +#include + +void test() +{ + const MPI_Comm mpi_comm = MPI_COMM_WORLD; + + parallel::distributed::Triangulation<2> tria(mpi_comm); + GridGenerator::hyper_cube(tria); + tria.refine_global(2); + unsigned int rank = Utilities::MPI::this_mpi_process(mpi_comm); + if (rank==0) + { + typename Triangulation<2>::active_cell_iterator cell = tria.begin_active(), + end_cell = tria.end(); + for (; cell!=end_cell; ++cell) + if (cell->is_locally_owned()) + cell->set_refine_flag(); + } + + tria.execute_coarsening_and_refinement(); + + std::vector::active_cell_iterator> > + vertex_to_cell = GridTools::vertex_to_cell_map(tria); + + AssertThrow(tria.n_vertices()==vertex_to_cell.size(), + ExcMessage("Wrong number of vertices")); + std::vector n_cells; + for (unsigned int i=0; i histogram(5,0); + for (unsigned int i=0; i histogram(5,0); + for (unsigned int i=0; i