]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add BBox::has_overlap_with() and fix NeighborType for 1D overlapping BBoxes 14361/head
authorJohannes Heinz <johannes.heinz@tuwien.ac.at>
Tue, 18 Oct 2022 15:52:20 +0000 (17:52 +0200)
committerJohannes Heinz <johannes.heinz@tuwien.ac.at>
Wed, 18 Jan 2023 16:58:31 +0000 (17:58 +0100)
doc/news/changes/minor/20230118Heinz [new file with mode: 0644]
include/deal.II/base/bounding_box.h
source/base/bounding_box.cc
tests/base/bounding_box_3.output
tests/base/bounding_box_8.cc [new file with mode: 0644]
tests/base/bounding_box_8.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20230118Heinz b/doc/news/changes/minor/20230118Heinz
new file mode 100644 (file)
index 0000000..011d2a0
--- /dev/null
@@ -0,0 +1,4 @@
+Improved: Introduce a function that only checks if Bounding Boxes are overlapping.
+Fixed: If a 1D BBox is completely contained in another one, it was not recognized as mergable neighbor.
+<br>
+(Johannes Heinz, 2022/07/07)
index dcca720b1707e9d4d955bbe814f75a8ddfef5079..9c2ba06c7a08ef03e95d7d51612a0822683ee2ff 100644 (file)
@@ -187,8 +187,14 @@ public:
   /**
    * Check if the current object and @p other_bbox are neighbors, i.e. if the boxes
    * have dimension spacedim, check if their intersection is non empty.
-   *
-   * Return an enumerator of type NeighborType.
+   */
+  bool
+  has_overlap_with(
+    const BoundingBox<spacedim, Number> &other_bbox,
+    const double tolerance = std::numeric_limits<Number>::epsilon()) const;
+
+  /**
+   * Check which NeighborType @p other_bbox is to the current object.
    */
   NeighborType
   get_neighbor_type(
index 0420499b78e22705bff2d577233f4a8c152117b7..79eda17a7af91ca8fadc660e1491696e42958f48 100644 (file)
@@ -30,12 +30,9 @@ BoundingBox<spacedim, Number>::point_inside(const Point<spacedim, Number> &p,
       // Bottom left-top right convention: the point is outside if it's smaller
       // than the first or bigger than the second boundary point The bounding
       // box is defined as a closed set
-      if ((p[i] < this->boundary_points.first[i] -
-                    tolerance * std::abs(this->boundary_points.second[i] -
-                                         this->boundary_points.first[i])) ||
-          (p[i] > this->boundary_points.second[i] +
-                    tolerance * std::abs(this->boundary_points.second[i] -
-                                         this->boundary_points.first[i])))
+      if ((p[i] <
+           this->boundary_points.first[i] - tolerance * side_length(i)) ||
+          (p[i] > this->boundary_points.second[i] + tolerance * side_length(i)))
         return false;
     }
   return true;
@@ -59,7 +56,23 @@ BoundingBox<spacedim, Number>::merge_with(
     }
 }
 
-
+template <int spacedim, typename Number>
+bool
+BoundingBox<spacedim, Number>::has_overlap_with(
+  const BoundingBox<spacedim, Number> &other_bbox,
+  const double                         tolerance) const
+{
+  for (unsigned int i = 0; i < spacedim; ++i)
+    {
+      // testing if the boxes are close enough to intersect
+      if ((other_bbox.boundary_points.second[i] <
+           this->boundary_points.first[i] - tolerance * side_length(i)) ||
+          (other_bbox.boundary_points.first[i] >
+           this->boundary_points.second[i] + tolerance * side_length(i)))
+        return false;
+    }
+  return true;
+}
 
 template <int spacedim, typename Number>
 NeighborType
@@ -67,14 +80,14 @@ BoundingBox<spacedim, Number>::get_neighbor_type(
   const BoundingBox<spacedim, Number> &other_bbox,
   const double                         tolerance) const
 {
+  if (!has_overlap_with(other_bbox, tolerance))
+    return NeighborType::not_neighbors;
+
   if (spacedim == 1)
     {
       // In dimension 1 if the two bounding box are neighbors
       // we can merge them
-      if (this->point_inside(other_bbox.boundary_points.first, tolerance) ||
-          this->point_inside(other_bbox.boundary_points.second, tolerance))
-        return NeighborType::mergeable_neighbors;
-      return NeighborType::not_neighbors;
+      return NeighborType::mergeable_neighbors;
     }
   else
     {
@@ -85,12 +98,6 @@ BoundingBox<spacedim, Number>::get_neighbor_type(
         {other_bbox.get_boundary_points().first,
          other_bbox.get_boundary_points().second}};
 
-      // Step 1: testing if the boxes are close enough to intersect
-      for (unsigned int d = 0; d < spacedim; ++d)
-        if (bbox1[0][d] * (1 - tolerance) > bbox2[1][d] ||
-            bbox2[0][d] * (1 - tolerance) > bbox1[1][d])
-          return NeighborType::not_neighbors;
-
       // The boxes intersect: we need to understand now how they intersect.
       // We begin by computing the intersection:
       std::array<double, spacedim> intersect_bbox_min;
index d12f6c43b74fe86cad2c753318b5d5d67a604e9a..74b9c830e74a18a00aa31d292a77e5ed11792e0a 100644 (file)
@@ -65,7 +65,7 @@ DEAL::Is neighbor G with A: 3
 DEAL::Vice-versa: 3
 DEAL::Is neighbor G with B: 3
 DEAL::Vice-versa: 3
-DEAL::Is neighbor G with C: 0
+DEAL::Is neighbor G with C: 3
 DEAL::Vice-versa: 3
 DEAL::Is neighbor G with D: 3
 DEAL::Vice-versa: 3
diff --git a/tests/base/bounding_box_8.cc b/tests/base/bounding_box_8.cc
new file mode 100644 (file)
index 0000000..f260bfc
--- /dev/null
@@ -0,0 +1,77 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 - 2021 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.
+//
+// ---------------------------------------------------------------------
+
+
+// test for BoundingBox<int dim> which tests the function
+// has_overlap_with()
+
+#include <deal.II/base/bounding_box.h>
+#include <deal.II/base/point.h>
+
+#include "../tests.h"
+
+BoundingBox<1>
+generate_bbox(const double left, const double right)
+{
+  Point<1> p1;
+  Point<1> p2;
+
+  p1[0] = left;
+  p2[0] = right;
+
+  BoundingBox<1> bbox(std::make_pair(p1, p2));
+
+  return bbox;
+}
+
+void
+test_bounding_box(const double left_a,
+                  const double right_a,
+                  const double left_b,
+                  const double right_b,
+                  const double tolerance)
+{
+  const auto bbox_a = generate_bbox(left_a, right_a);
+  const auto bbox_b = generate_bbox(left_b, right_b);
+
+  deallog << "Bounding box A: ";
+  deallog << "[" << bbox_a.get_boundary_points().first << ", "
+          << bbox_a.get_boundary_points().second << "]" << std::endl;
+  deallog << "Bounding box B: ";
+  deallog << "[" << bbox_b.get_boundary_points().first << ", "
+          << bbox_b.get_boundary_points().second << "]" << std::endl;
+  deallog << "Has overlap with: " << bbox_a.has_overlap_with(bbox_b, tolerance)
+          << std::endl;
+}
+
+int
+main()
+{
+  initlog();
+
+  test_bounding_box(1.0, 2.0, 2.0, 3.0, 1e-12);
+  test_bounding_box(1.0, 2.0, 2.0 + 1e-11, 3.0, 1e-12);
+  test_bounding_box(1.0, 2.0, 2.0 + 1e-11, 3.0, 1e-10);
+
+  test_bounding_box(-1.0, 0.0, 0.0, 1.0, 1e-12);
+  test_bounding_box(-1.0, 0.0, 0.0 + 1e-11, 1.0, 1e-12);
+  test_bounding_box(-1.0, 0.0, 0.0 + 1e-11, 1.0, 1e-10);
+
+  test_bounding_box(-2.0, -1.0, -1.0, 0.0, 1e-12);
+  test_bounding_box(-2.0, -1.0, -1.0 + 1e-11, 0.0, 1e-12);
+  test_bounding_box(-2.0, -1.0, -1.0 + 1e-11, 0.0, 1e-10);
+
+  return 0;
+}
diff --git a/tests/base/bounding_box_8.output b/tests/base/bounding_box_8.output
new file mode 100644 (file)
index 0000000..08517ba
--- /dev/null
@@ -0,0 +1,28 @@
+
+DEAL::Bounding box A: [1.00000, 2.00000]
+DEAL::Bounding box B: [2.00000, 3.00000]
+DEAL::Has overlap with: 1
+DEAL::Bounding box A: [1.00000, 2.00000]
+DEAL::Bounding box B: [2.00000, 3.00000]
+DEAL::Has overlap with: 0
+DEAL::Bounding box A: [1.00000, 2.00000]
+DEAL::Bounding box B: [2.00000, 3.00000]
+DEAL::Has overlap with: 1
+DEAL::Bounding box A: [-1.00000, 0.00000]
+DEAL::Bounding box B: [0.00000, 1.00000]
+DEAL::Has overlap with: 1
+DEAL::Bounding box A: [-1.00000, 0.00000]
+DEAL::Bounding box B: [1.00000e-11, 1.00000]
+DEAL::Has overlap with: 0
+DEAL::Bounding box A: [-1.00000, 0.00000]
+DEAL::Bounding box B: [1.00000e-11, 1.00000]
+DEAL::Has overlap with: 1
+DEAL::Bounding box A: [-2.00000, -1.00000]
+DEAL::Bounding box B: [-1.00000, 0.00000]
+DEAL::Has overlap with: 1
+DEAL::Bounding box A: [-2.00000, -1.00000]
+DEAL::Bounding box B: [-1.00000, 0.00000]
+DEAL::Has overlap with: 0
+DEAL::Bounding box A: [-2.00000, -1.00000]
+DEAL::Bounding box B: [-1.00000, 0.00000]
+DEAL::Has overlap with: 1

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.