From 6219117d71ee761da222df2ed3c1b7e489282ec8 Mon Sep 17 00:00:00 2001 From: danshapero Date: Tue, 26 Sep 2017 14:49:06 -0700 Subject: [PATCH] Added functions to compute bounding box for triangulation, accessors --- include/deal.II/grid/grid_tools.h | 17 ++++++ include/deal.II/grid/tria_accessor.h | 6 +++ source/grid/grid_tools.cc | 16 ++++++ source/grid/grid_tools.inst.in | 5 ++ source/grid/tria_accessor.cc | 23 ++++++++ tests/grid/grid_tools_bounding_box_03.cc | 56 ++++++++++++++++++++ tests/grid/grid_tools_bounding_box_03.output | 7 +++ tests/grid/measure_et_al_02.cc | 5 ++ tests/grid/measure_et_al_02.output | 6 +++ 9 files changed, 141 insertions(+) create mode 100644 tests/grid/grid_tools_bounding_box_03.cc create mode 100644 tests/grid/grid_tools_bounding_box_03.output diff --git a/include/deal.II/grid/grid_tools.h b/include/deal.II/grid/grid_tools.h index c03c8eb75b..aea8291e4b 100644 --- a/include/deal.II/grid/grid_tools.h +++ b/include/deal.II/grid/grid_tools.h @@ -18,6 +18,7 @@ #include +#include #include #include #include @@ -170,6 +171,22 @@ namespace GridTools template double cell_measure (const T &, ...); + /** + * Compute the smallest box containing the entire triangulation. + * + * If the input triangulation is a `parallel::distributed::Triangulation`, + * then each processor will compute a bounding box enclosing all locally + * owned, ghost, and artificial cells. In the case of a domain without curved + * boundaries, these bounding boxes will all agree between processors because + * the union of the areas occupied by artificial and ghost cells equals the + * union of the areas occupied by the cells that other processors own. + * However, if the domain has curved boundaries, this is no longer the case. + * The bounding box returned may be appropriate for the current processor, + * but different from the bounding boxes computed on other processors. + */ + template + BoundingBox compute_bounding_box(const Triangulation &triangulation); + /*@}*/ /** * @name Functions supporting the creation of meshes diff --git a/include/deal.II/grid/tria_accessor.h b/include/deal.II/grid/tria_accessor.h index d45031260f..81828d2a3b 100644 --- a/include/deal.II/grid/tria_accessor.h +++ b/include/deal.II/grid/tria_accessor.h @@ -17,6 +17,7 @@ #define dealii_tria_accessor_h +#include #include #include #include @@ -1300,6 +1301,11 @@ public: */ std::pair,double> enclosing_ball () const; + /** + * Return the smallest bounding box that encloses the object. + */ + BoundingBox bounding_box () const; + /** * Length of an object in the direction of the given axis, specified in the * local coordinate system. See the documentation of GeometryInfo for the diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc index e492b2fef4..6fd3160354 100644 --- a/source/grid/grid_tools.cc +++ b/source/grid/grid_tools.cc @@ -364,6 +364,22 @@ namespace GridTools } + + template + BoundingBox + compute_bounding_box(const Triangulation &tria) + { + using iterator = typename Triangulation::active_cell_iterator; + const auto predicate = [](const iterator &) + { + return true; + }; + + return compute_bounding_box(tria, std::function(predicate)); + } + + + template void delete_unused_vertices (std::vector > &vertices, diff --git a/source/grid/grid_tools.inst.in b/source/grid/grid_tools.inst.in index 2b3f29b8e1..df32b3a31c 100644 --- a/source/grid/grid_tools.inst.in +++ b/source/grid/grid_tools.inst.in @@ -159,6 +159,11 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS (const Triangulation &, const Mapping &); + template + BoundingBox + compute_bounding_box + (const Triangulation &); + template void delete_unused_vertices (std::vector > &, std::vector > &, diff --git a/source/grid/tria_accessor.cc b/source/grid/tria_accessor.cc index 106b9615a9..74cf9b5fa9 100644 --- a/source/grid/tria_accessor.cc +++ b/source/grid/tria_accessor.cc @@ -1110,6 +1110,29 @@ measure () const +template +BoundingBox +TriaAccessor:: +bounding_box () const +{ + std::pair, Point > boundary_points = + std::make_pair(this->vertex(0), this->vertex(0)); + + for (unsigned int v=1; v::vertices_per_cell; ++v) + { + const Point &x = this->vertex(v); + for (unsigned int k=0; k(boundary_points); +} + + + template <> double TriaAccessor<1,1,1>::extent_in_direction(const unsigned int axis) const { diff --git a/tests/grid/grid_tools_bounding_box_03.cc b/tests/grid/grid_tools_bounding_box_03.cc new file mode 100644 index 0000000000..a5939dce2d --- /dev/null +++ b/tests/grid/grid_tools_bounding_box_03.cc @@ -0,0 +1,56 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + +// Check computing bounding boxes of entire triangulations and individual cells + +#include "../tests.h" + +#include +#include + +template +void test_tria_bounding_box() +{ + dealii::Point p1, p2; + for (unsigned int k = 0; k < dim; ++k) + { + p1[k] = k; + p2[k] = 2 * k + 1; + } + + Triangulation tria; + GridGenerator::hyper_rectangle(tria, p1, p2); + tria.refine_global(1); + const BoundingBox bounding_box = GridTools::compute_bounding_box(tria); + const std::pair, Point > &boundary_points = + bounding_box.get_boundary_points(); + + deallog << boundary_points.first << ", " << boundary_points.second << std::endl; +} + +int main(void) +{ + initlog(); + + test_tria_bounding_box<1, 1>(); + test_tria_bounding_box<1, 2>(); + test_tria_bounding_box<2, 2>(); + test_tria_bounding_box<1, 3>(); + test_tria_bounding_box<2, 3>(); + test_tria_bounding_box<3, 3>(); + + return 0; +} + diff --git a/tests/grid/grid_tools_bounding_box_03.output b/tests/grid/grid_tools_bounding_box_03.output new file mode 100644 index 0000000000..e7c415a7d1 --- /dev/null +++ b/tests/grid/grid_tools_bounding_box_03.output @@ -0,0 +1,7 @@ + +DEAL::0.00000, 1.00000 +DEAL::0.00000 0.00000, 1.00000 0.00000 +DEAL::0.00000 1.00000, 1.00000 3.00000 +DEAL::0.00000 0.00000 0.00000, 1.00000 0.00000 0.00000 +DEAL::0.00000 1.00000 0.00000, 1.00000 3.00000 0.00000 +DEAL::0.00000 1.00000 2.00000, 1.00000 3.00000 5.00000 diff --git a/tests/grid/measure_et_al_02.cc b/tests/grid/measure_et_al_02.cc index 464871e186..31ce5f90bf 100644 --- a/tests/grid/measure_et_al_02.cc +++ b/tests/grid/measure_et_al_02.cc @@ -100,6 +100,11 @@ void test() << tria.begin_active()->extent_in_direction(0) << std::endl; deallog << "dim" << dim << ":case" << case_no << ":minimum_vertex_distance=" << tria.begin_active()->minimum_vertex_distance() << std::endl; + + const BoundingBox box = tria.begin_active()->bounding_box(); + const std::pair, Point > &pts = box.get_boundary_points(); + deallog << "dim" << dim << ":case" << case_no << ":bounding_box=" + << pts.first << ", " << pts.second << std::endl; tria.clear(); } } diff --git a/tests/grid/measure_et_al_02.output b/tests/grid/measure_et_al_02.output index e1db82bb4a..d7ec0015ed 100644 --- a/tests/grid/measure_et_al_02.output +++ b/tests/grid/measure_et_al_02.output @@ -2,18 +2,24 @@ DEAL::dim2:case0:diameter=2.8284 DEAL::dim2:case0:extent_in_direction=2.0000 DEAL::dim2:case0:minimum_vertex_distance=2.0000 +DEAL::dim2:case0:bounding_box=1.0000 1.0000, 3.0000 3.0000 DEAL::dim2:case1:diameter=3.6056 DEAL::dim2:case1:extent_in_direction=3.0000 DEAL::dim2:case1:minimum_vertex_distance=2.0000 +DEAL::dim2:case1:bounding_box=0.0000 1.0000, 3.0000 3.0000 DEAL::dim2:case2:diameter=4.4721 DEAL::dim2:case2:extent_in_direction=3.0000 DEAL::dim2:case2:minimum_vertex_distance=2.2361 +DEAL::dim2:case2:bounding_box=0.0000 1.0000, 4.0000 3.0000 DEAL::dim3:case0:diameter=3.4641 DEAL::dim3:case0:extent_in_direction=2.0000 DEAL::dim3:case0:minimum_vertex_distance=2.0000 +DEAL::dim3:case0:bounding_box=1.0000 1.0000 1.0000, 3.0000 3.0000 3.0000 DEAL::dim3:case1:diameter=4.1231 DEAL::dim3:case1:extent_in_direction=3.0000 DEAL::dim3:case1:minimum_vertex_distance=2.0000 +DEAL::dim3:case1:bounding_box=0.0000 1.0000 1.0000, 3.0000 3.0000 3.0000 DEAL::dim3:case2:diameter=4.1231 DEAL::dim3:case2:extent_in_direction=3.0000 DEAL::dim3:case2:minimum_vertex_distance=2.0000 +DEAL::dim3:case2:bounding_box=0.0000 1.0000 1.0000, 3.0000 3.0000 3.0000 -- 2.39.5