From: Daniel Arndt Date: Wed, 1 Aug 2018 18:07:17 +0000 (+0200) Subject: Allow merging multiple Triangulation objects at once X-Git-Tag: v9.1.0-rc1~851^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=fd8a21a6a7a68730e0d41c8ceb6a16f1304fb7f2;p=dealii.git Allow merging multiple Triangulation objects at once --- diff --git a/doc/news/changes/minor/20180801DanielArndt b/doc/news/changes/minor/20180801DanielArndt new file mode 100644 index 0000000000..f29df0e876 --- /dev/null +++ b/doc/news/changes/minor/20180801DanielArndt @@ -0,0 +1,4 @@ +New: GridGenerator::merge_triangulations can merge multiple +Triangulation objects at once. +
+(Daniel Arndt, 2018/08/01) diff --git a/include/deal.II/grid/grid_generator.h b/include/deal.II/grid/grid_generator.h index 56da7be560..a96ad98b0d 100644 --- a/include/deal.II/grid/grid_generator.h +++ b/include/deal.II/grid/grid_generator.h @@ -1164,6 +1164,18 @@ namespace GridGenerator const double duplicated_vertex_tolerance = 1.0e-12, const bool copy_manifold_ids = false); + /** + * Same as above but allows to merge more than two triangulations at once. + */ + template + void + merge_triangulations( + const std::initializer_list *const> + triangulations, + Triangulation &result, + const double duplicated_vertex_tolerance = 1.0e-12, + const bool copy_manifold_ids = false); + /** * Given the two triangulations specified as the first two arguments, create * the triangulation that contains the finest cells of both triangulation diff --git a/source/grid/grid_generator.cc b/source/grid/grid_generator.cc index e868b93980..50d0cb86fc 100644 --- a/source/grid/grid_generator.cc +++ b/source/grid/grid_generator.cc @@ -4139,61 +4139,55 @@ namespace GridGenerator template void - merge_triangulations(const Triangulation &triangulation_1, - const Triangulation &triangulation_2, - Triangulation & result, - const double duplicated_vertex_tolerance, - const bool copy_manifold_ids) + merge_triangulations( + const std::initializer_list *const> + triangulations, + Triangulation &result, + const double duplicated_vertex_tolerance, + const bool copy_manifold_ids) { - // if either Triangulation is empty then merging is just a copy. - if (triangulation_1.n_cells() == 0) - { - result.copy_triangulation(triangulation_2); - return; - } - if (triangulation_2.n_cells() == 0) + for (const auto triangulation : triangulations) { - result.copy_triangulation(triangulation_1); - return; + (void)triangulation; + Assert(triangulation->n_levels() == 1, + ExcMessage("The input triangulations must be non-empty " + "and must not be refined.")); } - Assert(triangulation_1.n_levels() == 1, - ExcMessage("The input triangulations must be non-empty " - "and must not be refined.")); - Assert(triangulation_2.n_levels() == 1, - ExcMessage("The input triangulations must be non-empty " - "and must not be refined.")); - // get the union of the set of vertices - std::vector> vertices = triangulation_1.get_vertices(); - vertices.insert(vertices.end(), - triangulation_2.get_vertices().begin(), - triangulation_2.get_vertices().end()); + unsigned int n_vertices = 0; + for (const auto triangulation : triangulations) + n_vertices += triangulation->n_vertices(); + + std::vector> vertices; + vertices.reserve(n_vertices); + for (const auto triangulation : triangulations) + vertices.insert(vertices.end(), + triangulation->get_vertices().begin(), + triangulation->get_vertices().end()); // now form the union of the set of cells std::vector> cells; - cells.reserve(triangulation_1.n_cells() + triangulation_2.n_cells()); - for (const auto &cell : triangulation_1.cell_iterators()) - { - CellData this_cell; - for (unsigned int v = 0; v < GeometryInfo::vertices_per_cell; ++v) - this_cell.vertices[v] = cell->vertex_index(v); - this_cell.material_id = cell->material_id(); - this_cell.manifold_id = cell->manifold_id(); - cells.push_back(std::move(this_cell)); - } + unsigned int n_cells = 0; + for (const auto triangulation : triangulations) + n_cells += triangulation->n_cells(); - // now do the same for the other other mesh. note that we have to - // translate the vertex indices - for (const auto &cell : triangulation_2.cell_iterators()) + cells.reserve(n_cells); + unsigned int n_accumulated_vertices = 0; + for (const auto triangulation : triangulations) { - CellData this_cell; - for (unsigned int v = 0; v < GeometryInfo::vertices_per_cell; ++v) - this_cell.vertices[v] = - cell->vertex_index(v) + triangulation_1.n_vertices(); - this_cell.material_id = cell->material_id(); - this_cell.manifold_id = cell->manifold_id(); - cells.push_back(std::move(this_cell)); + for (const auto &cell : triangulation->cell_iterators()) + { + CellData this_cell; + for (unsigned int v = 0; v < GeometryInfo::vertices_per_cell; + ++v) + this_cell.vertices[v] = + cell->vertex_index(v) + n_accumulated_vertices; + this_cell.material_id = cell->material_id(); + this_cell.manifold_id = cell->manifold_id(); + cells.push_back(std::move(this_cell)); + } + n_accumulated_vertices += triangulation->n_vertices(); } // Now for the SubCellData @@ -4207,8 +4201,11 @@ namespace GridGenerator case 2: { - subcell_data.boundary_lines.reserve(triangulation_1.n_lines() + - triangulation_2.n_lines()); + unsigned int n_boundary_lines = 0; + for (const auto triangulation : triangulations) + n_boundary_lines += triangulation->n_lines(); + + subcell_data.boundary_lines.reserve(n_boundary_lines); auto copy_line_manifold_ids = [&](const Triangulation &tria, @@ -4230,23 +4227,28 @@ namespace GridGenerator } }; - copy_line_manifold_ids(triangulation_1, 0); - copy_line_manifold_ids(triangulation_2, - triangulation_1.n_vertices()); + unsigned int n_accumulated_vertices = 0; + for (const auto triangulation : triangulations) + { + copy_line_manifold_ids(*triangulation, + n_accumulated_vertices); + n_accumulated_vertices += triangulation->n_vertices(); + } break; } case 3: { - subcell_data.boundary_quads.reserve(triangulation_1.n_quads() + - triangulation_2.n_quads()); + unsigned int n_boundary_quads = 0; + for (const auto triangulation : triangulations) + n_boundary_quads += triangulation->n_quads(); + + subcell_data.boundary_quads.reserve(n_boundary_quads); // we can't do better here than to loop over all the lines - // bounding a face. For regular meshes an (interior) line in 3D - // is part of four cells. So this should be an appropriate + // bounding a face. For regular meshes an (interior) line in + // 3D is part of four cells. So this should be an appropriate // guess. - subcell_data.boundary_lines.reserve( - triangulation_1.n_cells() * 4 + - triangulation_2.n_cells() * 4); + subcell_data.boundary_lines.reserve(4 * n_boundary_quads); auto copy_face_and_line_manifold_ids = [&](const Triangulation &tria, @@ -4287,10 +4289,13 @@ namespace GridGenerator boundary_line); // trivially_copyable } }; - - copy_face_and_line_manifold_ids(triangulation_1, 0); - copy_face_and_line_manifold_ids(triangulation_2, - triangulation_1.n_vertices()); + unsigned int n_accumulated_vertices = 0; + for (const auto triangulation : triangulations) + { + copy_face_and_line_manifold_ids(*triangulation, + n_accumulated_vertices); + n_accumulated_vertices += triangulation->n_vertices(); + } break; } default: @@ -4298,14 +4303,19 @@ namespace GridGenerator } } - // throw out duplicated vertices from the two meshes, reorder vertices as - // necessary and create the triangulation - std::vector considered_vertices; - GridTools::delete_duplicated_vertices(vertices, - cells, - subcell_data, - considered_vertices, - duplicated_vertex_tolerance); + // throw out duplicated vertices + // since GridTools::delete_duplicated_vertices can only consider + // one duplications at once, we have to call the function one time less + // than the number of Triangulation objects to be merged. + for (unsigned int i = 0; i < triangulations.size(); ++i) + { + std::vector considered_vertices; + GridTools::delete_duplicated_vertices(vertices, + cells, + subcell_data, + considered_vertices, + duplicated_vertex_tolerance); + } // reorder the cells to ensure that they satisfy the convention for // edge and face directions @@ -4315,6 +4325,34 @@ namespace GridGenerator } + + template + void + merge_triangulations(const Triangulation &triangulation_1, + const Triangulation &triangulation_2, + Triangulation & result, + const double duplicated_vertex_tolerance, + const bool copy_manifold_ids) + { + // if either Triangulation is empty then merging is just a copy. + if (triangulation_1.n_cells() == 0) + { + result.copy_triangulation(triangulation_2); + return; + } + if (triangulation_2.n_cells() == 0) + { + result.copy_triangulation(triangulation_1); + return; + } + merge_triangulations({&triangulation_1, &triangulation_2}, + result, + duplicated_vertex_tolerance, + copy_manifold_ids); + } + + + template void create_union_triangulation( diff --git a/source/grid/grid_generator.inst.in b/source/grid/grid_generator.inst.in index 6ac4262276..659485982f 100644 --- a/source/grid/grid_generator.inst.in +++ b/source/grid/grid_generator.inst.in @@ -69,6 +69,15 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) Triangulation &, const std::vector &); + template void + merge_triangulations( + const std::initializer_list< + const Triangulation + *const> triangulations, + Triangulation &result, + const double duplicated_vertex_tolerance = 1.0e-12, + const bool copy_manifold_ids = false); + template void merge_triangulations( const Triangulation diff --git a/tests/grid/merge_triangulations_06.cc b/tests/grid/merge_triangulations_06.cc new file mode 100644 index 0000000000..4b046796bb --- /dev/null +++ b/tests/grid/merge_triangulations_06.cc @@ -0,0 +1,171 @@ +// --------------------------------------------------------------------- +// +// 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.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + + +#include +#include +#include +#include +#include +#include + +#include "../tests.h" + + +// Test the version of GridTools::merge_triangulations that takes a +// std::initializer_list. + +void +test_2d() +{ + constexpr int dim = 2; + Triangulation tria_1; + { + GridGenerator::hyper_cube(tria_1); + Tensor<1, dim> shift_vector; + shift_vector[0] = -1.; + shift_vector[1] = -1.; + GridTools::shift(shift_vector, tria_1); + } + + Triangulation tria_2; + { + GridGenerator::hyper_cube(tria_2); + Tensor<1, dim> shift_vector; + shift_vector[1] = -1.; + GridTools::shift(shift_vector, tria_2); + } + + Triangulation tria_3; + { + GridGenerator::hyper_cube(tria_3); + Tensor<1, dim> shift_vector; + shift_vector[0] = -1.; + GridTools::shift(shift_vector, tria_3); + } + + Triangulation tria_4; + { + GridGenerator::hyper_cube(tria_4); + } + + // now merge triangulations + Triangulation result; + GridGenerator::merge_triangulations({&tria_1, &tria_2, &tria_3, &tria_4}, + result); + + GridOut().write_gnuplot(result, deallog.get_file_stream()); + + deallog << " Total number of cells = " << result.n_cells() + << std::endl + << " Total number of vertices = " << result.n_used_vertices() + << std::endl; +} + +void +test_3d() +{ + constexpr int dim = 3; + Triangulation tria_1; + { + GridGenerator::hyper_cube(tria_1); + Tensor<1, dim> shift_vector; + shift_vector[0] = -1.; + shift_vector[1] = -1.; + shift_vector[2] = -1.; + GridTools::shift(shift_vector, tria_1); + } + + Triangulation tria_2; + { + GridGenerator::hyper_cube(tria_2); + Tensor<1, dim> shift_vector; + shift_vector[1] = -1.; + shift_vector[2] = -1.; + GridTools::shift(shift_vector, tria_2); + } + + Triangulation tria_3; + { + GridGenerator::hyper_cube(tria_3); + Tensor<1, dim> shift_vector; + shift_vector[0] = -1.; + shift_vector[2] = -1.; + GridTools::shift(shift_vector, tria_3); + } + + Triangulation tria_4; + { + Tensor<1, dim> shift_vector; + shift_vector[2] = -1.; + GridGenerator::hyper_cube(tria_4); + GridTools::shift(shift_vector, tria_4); + } + + Triangulation tria_5; + { + GridGenerator::hyper_cube(tria_5); + Tensor<1, dim> shift_vector; + shift_vector[0] = -1.; + shift_vector[1] = -1.; + GridTools::shift(shift_vector, tria_5); + } + + Triangulation tria_6; + { + GridGenerator::hyper_cube(tria_6); + Tensor<1, dim> shift_vector; + shift_vector[1] = -1.; + GridTools::shift(shift_vector, tria_6); + } + + Triangulation tria_7; + { + GridGenerator::hyper_cube(tria_7); + Tensor<1, dim> shift_vector; + shift_vector[0] = -1.; + GridTools::shift(shift_vector, tria_7); + } + + Triangulation tria_8; + { + GridGenerator::hyper_cube(tria_8); + } + + // now merge triangulations + Triangulation result; + GridGenerator::merge_triangulations( + {&tria_1, &tria_2, &tria_3, &tria_4, &tria_5, &tria_6, &tria_7, &tria_8}, + result); + + GridOut().write_gnuplot(result, deallog.get_file_stream()); + + deallog << " Total number of cells = " << result.n_cells() + << std::endl + << " Total number of vertices = " << result.n_used_vertices() + << std::endl; +} + +int +main() +{ + initlog(); + + test_2d(); + test_3d(); + + return 0; +} diff --git a/tests/grid/merge_triangulations_06.output b/tests/grid/merge_triangulations_06.output new file mode 100644 index 0000000000..3248035357 --- /dev/null +++ b/tests/grid/merge_triangulations_06.output @@ -0,0 +1,225 @@ + +-1.00000 -1.00000 0 0 +0.00000 -1.00000 0 0 +0.00000 0.00000 0 0 +-1.00000 0.00000 0 0 +-1.00000 -1.00000 0 0 + + +0.00000 -1.00000 0 0 +1.00000 -1.00000 0 0 +1.00000 0.00000 0 0 +0.00000 0.00000 0 0 +0.00000 -1.00000 0 0 + + +-1.00000 0.00000 0 0 +0.00000 0.00000 0 0 +0.00000 1.00000 0 0 +-1.00000 1.00000 0 0 +-1.00000 0.00000 0 0 + + +0.00000 0.00000 0 0 +1.00000 0.00000 0 0 +1.00000 1.00000 0 0 +0.00000 1.00000 0 0 +0.00000 0.00000 0 0 + + +DEAL:: Total number of cells = 4 +DEAL:: Total number of vertices = 9 +-1.00000 -1.00000 -1.00000 0 0 +0.00000 -1.00000 -1.00000 0 0 +0.00000 -1.00000 0.00000 0 0 +-1.00000 -1.00000 0.00000 0 0 +-1.00000 -1.00000 -1.00000 0 0 + +-1.00000 0.00000 -1.00000 0 0 +0.00000 0.00000 -1.00000 0 0 +0.00000 0.00000 0.00000 0 0 +-1.00000 0.00000 0.00000 0 0 +-1.00000 0.00000 -1.00000 0 0 + +-1.00000 -1.00000 -1.00000 0 0 +-1.00000 0.00000 -1.00000 0 0 + +0.00000 -1.00000 -1.00000 0 0 +0.00000 0.00000 -1.00000 0 0 + +0.00000 -1.00000 0.00000 0 0 +0.00000 0.00000 0.00000 0 0 + +-1.00000 -1.00000 0.00000 0 0 +-1.00000 0.00000 0.00000 0 0 + +0.00000 -1.00000 -1.00000 0 0 +1.00000 -1.00000 -1.00000 0 0 +1.00000 -1.00000 0.00000 0 0 +0.00000 -1.00000 0.00000 0 0 +0.00000 -1.00000 -1.00000 0 0 + +0.00000 0.00000 -1.00000 0 0 +1.00000 0.00000 -1.00000 0 0 +1.00000 0.00000 0.00000 0 0 +0.00000 0.00000 0.00000 0 0 +0.00000 0.00000 -1.00000 0 0 + +0.00000 -1.00000 -1.00000 0 0 +0.00000 0.00000 -1.00000 0 0 + +1.00000 -1.00000 -1.00000 0 0 +1.00000 0.00000 -1.00000 0 0 + +1.00000 -1.00000 0.00000 0 0 +1.00000 0.00000 0.00000 0 0 + +0.00000 -1.00000 0.00000 0 0 +0.00000 0.00000 0.00000 0 0 + +-1.00000 0.00000 -1.00000 0 0 +0.00000 0.00000 -1.00000 0 0 +0.00000 0.00000 0.00000 0 0 +-1.00000 0.00000 0.00000 0 0 +-1.00000 0.00000 -1.00000 0 0 + +-1.00000 1.00000 -1.00000 0 0 +0.00000 1.00000 -1.00000 0 0 +0.00000 1.00000 0.00000 0 0 +-1.00000 1.00000 0.00000 0 0 +-1.00000 1.00000 -1.00000 0 0 + +-1.00000 0.00000 -1.00000 0 0 +-1.00000 1.00000 -1.00000 0 0 + +0.00000 0.00000 -1.00000 0 0 +0.00000 1.00000 -1.00000 0 0 + +0.00000 0.00000 0.00000 0 0 +0.00000 1.00000 0.00000 0 0 + +-1.00000 0.00000 0.00000 0 0 +-1.00000 1.00000 0.00000 0 0 + +0.00000 0.00000 -1.00000 0 0 +1.00000 0.00000 -1.00000 0 0 +1.00000 0.00000 0.00000 0 0 +0.00000 0.00000 0.00000 0 0 +0.00000 0.00000 -1.00000 0 0 + +0.00000 1.00000 -1.00000 0 0 +1.00000 1.00000 -1.00000 0 0 +1.00000 1.00000 0.00000 0 0 +0.00000 1.00000 0.00000 0 0 +0.00000 1.00000 -1.00000 0 0 + +0.00000 0.00000 -1.00000 0 0 +0.00000 1.00000 -1.00000 0 0 + +1.00000 0.00000 -1.00000 0 0 +1.00000 1.00000 -1.00000 0 0 + +1.00000 0.00000 0.00000 0 0 +1.00000 1.00000 0.00000 0 0 + +0.00000 0.00000 0.00000 0 0 +0.00000 1.00000 0.00000 0 0 + +-1.00000 -1.00000 0.00000 0 0 +0.00000 -1.00000 0.00000 0 0 +0.00000 -1.00000 1.00000 0 0 +-1.00000 -1.00000 1.00000 0 0 +-1.00000 -1.00000 0.00000 0 0 + +-1.00000 0.00000 0.00000 0 0 +0.00000 0.00000 0.00000 0 0 +0.00000 0.00000 1.00000 0 0 +-1.00000 0.00000 1.00000 0 0 +-1.00000 0.00000 0.00000 0 0 + +-1.00000 -1.00000 0.00000 0 0 +-1.00000 0.00000 0.00000 0 0 + +0.00000 -1.00000 0.00000 0 0 +0.00000 0.00000 0.00000 0 0 + +0.00000 -1.00000 1.00000 0 0 +0.00000 0.00000 1.00000 0 0 + +-1.00000 -1.00000 1.00000 0 0 +-1.00000 0.00000 1.00000 0 0 + +0.00000 -1.00000 0.00000 0 0 +1.00000 -1.00000 0.00000 0 0 +1.00000 -1.00000 1.00000 0 0 +0.00000 -1.00000 1.00000 0 0 +0.00000 -1.00000 0.00000 0 0 + +0.00000 0.00000 0.00000 0 0 +1.00000 0.00000 0.00000 0 0 +1.00000 0.00000 1.00000 0 0 +0.00000 0.00000 1.00000 0 0 +0.00000 0.00000 0.00000 0 0 + +0.00000 -1.00000 0.00000 0 0 +0.00000 0.00000 0.00000 0 0 + +1.00000 -1.00000 0.00000 0 0 +1.00000 0.00000 0.00000 0 0 + +1.00000 -1.00000 1.00000 0 0 +1.00000 0.00000 1.00000 0 0 + +0.00000 -1.00000 1.00000 0 0 +0.00000 0.00000 1.00000 0 0 + +-1.00000 0.00000 0.00000 0 0 +0.00000 0.00000 0.00000 0 0 +0.00000 0.00000 1.00000 0 0 +-1.00000 0.00000 1.00000 0 0 +-1.00000 0.00000 0.00000 0 0 + +-1.00000 1.00000 0.00000 0 0 +0.00000 1.00000 0.00000 0 0 +0.00000 1.00000 1.00000 0 0 +-1.00000 1.00000 1.00000 0 0 +-1.00000 1.00000 0.00000 0 0 + +-1.00000 0.00000 0.00000 0 0 +-1.00000 1.00000 0.00000 0 0 + +0.00000 0.00000 0.00000 0 0 +0.00000 1.00000 0.00000 0 0 + +0.00000 0.00000 1.00000 0 0 +0.00000 1.00000 1.00000 0 0 + +-1.00000 0.00000 1.00000 0 0 +-1.00000 1.00000 1.00000 0 0 + +0.00000 0.00000 0.00000 0 0 +1.00000 0.00000 0.00000 0 0 +1.00000 0.00000 1.00000 0 0 +0.00000 0.00000 1.00000 0 0 +0.00000 0.00000 0.00000 0 0 + +0.00000 1.00000 0.00000 0 0 +1.00000 1.00000 0.00000 0 0 +1.00000 1.00000 1.00000 0 0 +0.00000 1.00000 1.00000 0 0 +0.00000 1.00000 0.00000 0 0 + +0.00000 0.00000 0.00000 0 0 +0.00000 1.00000 0.00000 0 0 + +1.00000 0.00000 0.00000 0 0 +1.00000 1.00000 0.00000 0 0 + +1.00000 0.00000 1.00000 0 0 +1.00000 1.00000 1.00000 0 0 + +0.00000 0.00000 1.00000 0 0 +0.00000 1.00000 1.00000 0 0 + +DEAL:: Total number of cells = 8 +DEAL:: Total number of vertices = 27