]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Created class Bounding Box and its routine point_inside 5078/head
authorGiovanni Alzetta <giovannialzetta@hotmail.it>
Sun, 17 Sep 2017 23:25:10 +0000 (23:25 +0000)
committerGiovanni Alzetta <giovannialzetta@hotmail.it>
Sun, 17 Sep 2017 23:25:10 +0000 (23:25 +0000)
doc/news/changes/minor/20170915GiovanniAlzetta [new file with mode: 0644]
include/deal.II/base/bounding_box.h [new file with mode: 0644]
source/base/CMakeLists.txt
source/base/bounding_box.cc [new file with mode: 0644]
source/base/bounding_box.inst.in [new file with mode: 0644]
tests/base/bounding_box_1.cc [new file with mode: 0644]
tests/base/bounding_box_1.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20170915GiovanniAlzetta b/doc/news/changes/minor/20170915GiovanniAlzetta
new file mode 100644 (file)
index 0000000..c3bc50a
--- /dev/null
@@ -0,0 +1,4 @@
+New: Added the base class BoundingBox. The BoundingBox is constructed from a std::pair of points in real space, ordered following the convention bottom-left, top-right.
+The class has the method BoundingBox::point_inside(Point<spacedim> p) which checks for the presence of a point inside it.
+<br>
+(Giovanni Alzetta, 2017/09/15)
diff --git a/include/deal.II/base/bounding_box.h b/include/deal.II/base/bounding_box.h
new file mode 100644 (file)
index 0000000..bec903d
--- /dev/null
@@ -0,0 +1,108 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_base_bounding_box_h
+#define dealii_base_bounding_box_h
+
+
+#include <deal.II/base/config.h>
+#include <deal.II/base/exceptions.h>
+#include <deal.II/base/point.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+/**
+ * A class that represents a bounding box in a space with arbitrary dimension
+ * <tt>spacedim</tt>.
+ *
+ * Objects of this class are used to represent bounding boxes. They are,
+ * among other uses, useful in parallel distributed meshes to give a general
+ * description the owners of each portion of the mesh.
+ *
+ * Bounding boxes are represented by two vertices (bottom left and top right).
+ * Geometrically, a bounding box is:
+ * - 1 d: a segment (represented by its vertices in the proper order)
+ * - 2 d: a rectangle (represented by the vertices V at bottom left, top right)
+ * @code
+ * .--------V
+ * |        |
+ * V--------.
+ * @endcode
+ *
+ * - 3 d: a cuboid (in which case the two vertices V follow the convention and
+ * are not owned by the same face)
+ * @code
+ *   .------V
+ *  /      /|
+ * .------. |
+ * |      | /
+ * |      |/
+ * V------.
+ * @endcode
+ * Notice the sides are always parallel to the respective axis.
+ *
+ */
+template <int spacedim, typename Number=double>
+class BoundingBox
+{
+public:
+  /**
+   * Standard constructor. Creates an object that corresponds to an empty box,
+   * i.e. a degenerate box with both points being the origin.
+   */
+  BoundingBox () = default;
+
+  /**
+   * Standard constructor for non-empty boxes: it uses a pair of points
+   * which describe the box: one for the bottom and one for the top
+   * corner.
+   */
+  BoundingBox (const std::pair<Point<spacedim,Number>,Point<spacedim,Number>> &boundary_points);
+
+  /**
+   * Returns the boundary_points
+   */
+  const std::pair<Point<spacedim,Number>,Point<spacedim,Number>> &get_boundary_points () const;
+
+  /**
+   * Returns true if the point is inside the Bounding Box, false otherwise
+   */
+  bool point_inside (const Point<spacedim, Number> &p) const;
+
+private:
+  std::pair<Point<spacedim, Number>,Point<spacedim,Number>> boundary_points;
+};
+
+/*------------------------------- Inline functions: Point ---------------------------*/
+
+#ifndef DOXYGEN
+
+
+template <int spacedim, typename Number>
+inline
+BoundingBox<spacedim, Number>::BoundingBox (const std::pair<Point<spacedim,Number>,Point<spacedim,Number>> &boundary_points)
+{
+  //We check the Bounding Box is not degenerate
+  for (unsigned int i=0; i<spacedim; ++i)
+    Assert (boundary_points.first[i] < boundary_points.second[i],
+            ExcMessage ("Bounding Box can't be created: the point's order should be bottom left, top right!"));
+
+  this->boundary_points = boundary_points;
+}
+
+#endif // DOXYGEN
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index af9cf9853677a78b163ff5d33888f2690d227d16..172607dcafa2df5d0ba47b4916493fd28c88fb4f 100644 (file)
@@ -22,6 +22,7 @@ INCLUDE_DIRECTORIES(BEFORE ${CMAKE_CURRENT_BINARY_DIR})
 #
 SET(_unity_include_src
   auto_derivative_function.cc
+  bounding_box.cc
   conditional_ostream.cc
   config.cc
   convergence_table.cc
@@ -93,6 +94,7 @@ SETUP_SOURCE_LIST("${_unity_include_src}"
   )
 
 SET(_inst
+  bounding_box.inst.in
   data_out_base.inst.in
   function.inst.in
   function_time.inst.in
diff --git a/source/base/bounding_box.cc b/source/base/bounding_box.cc
new file mode 100644 (file)
index 0000000..f2ef5b4
--- /dev/null
@@ -0,0 +1,41 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/base/bounding_box.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+template <int spacedim, typename Number>
+bool BoundingBox<spacedim,Number>::point_inside (const Point<spacedim, Number> &p) const
+{
+  for (unsigned int i=0; i < spacedim; ++i)
+    {
+      //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] || this->boundary_points.second[i] < p[i])
+        return false;
+    }
+  return true;
+}
+
+template <int spacedim, typename Number>
+const std::pair<Point<spacedim,Number>,Point<spacedim,Number>> &BoundingBox<spacedim,Number>::get_boundary_points () const
+{
+  return this->boundary_points;
+}
+
+#include "bounding_box.inst"
+DEAL_II_NAMESPACE_CLOSE
diff --git a/source/base/bounding_box.inst.in b/source/base/bounding_box.inst.in
new file mode 100644 (file)
index 0000000..20af0a2
--- /dev/null
@@ -0,0 +1,21 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+
+for (deal_II_dimension : SPACE_DIMENSIONS; number : REAL_SCALARS)
+{
+    template
+    class BoundingBox<deal_II_dimension,number>;
+}
diff --git a/tests/base/bounding_box_1.cc b/tests/base/bounding_box_1.cc
new file mode 100644 (file)
index 0000000..3507294
--- /dev/null
@@ -0,0 +1,137 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+
+// test for BoundingBox<unsigned int spacedim> which test basic stuff:
+// creation in various dimensions and checking for points inside
+
+#include "../tests.h"
+
+#include <deal.II/base/point.h>
+#include <deal.II/base/bounding_box.h>
+
+template <int spacedim>
+void test_bounding_box()
+{
+  BoundingBox<spacedim> a;
+  deallog << "Empty constructor: " << std::endl;
+  deallog << a.get_boundary_points().first << std::endl;
+  deallog << a.get_boundary_points().second << std::endl;
+
+  std::pair<Point<spacedim>,Point<spacedim>> boundaries;
+
+  for (int i=0; i<spacedim; i++)
+    {
+      boundaries.first[i] = 0.2 -i*0.2;
+      boundaries.second[i] = 0.8 + i*0.8;
+    }
+
+  BoundingBox<spacedim> b(boundaries);
+  deallog << "Boundary points: " << std::endl;
+  deallog << b.get_boundary_points().first << std::endl;
+  deallog << b.get_boundary_points().second << std::endl;
+
+  deallog << "Boundary points are inside: "
+          << b.point_inside(boundaries.first) << " "
+          << b.point_inside(boundaries.second) << std::endl;
+
+  std::vector<Point<spacedim>> test_points;
+
+  //To guarantee points are inside we take a convex combination
+  double c;
+  for (int t=0; t<5; ++t)
+    {
+      Point<spacedim> test_pt;
+      c = (1+t)/5.0;
+      for (unsigned int i=0; i<spacedim; ++i)
+        test_pt[i] = boundaries.first[i]*c+boundaries.second[i]*(1-c);
+
+      test_points.push_back(test_pt);
+    }
+
+  deallog << "Points inside: " << std::endl;
+  for (unsigned int i=0; i<test_points.size(); ++i)
+    deallog << test_points[i] << " is inside: " <<
+            b.point_inside(test_points[i]) << std::endl;
+
+  deallog << std::endl;
+  test_points.clear();
+
+  //To create outside points we take a non-convex combination
+  for (int t=0; t<5; ++t)
+    {
+      Point<spacedim> test_pt;
+      c = (1+t)*2.5;
+      if (t%2==0) c=-c; //Changing the sign sometimes..
+      for (unsigned int i=0; i<spacedim; ++i)
+        test_pt[i] = boundaries.first[i]*c+boundaries.second[i]*c;
+
+      test_points.push_back(test_pt);
+    }
+
+  deallog << "Points outside:" << std::endl;
+  for (unsigned int i=0; i<test_points.size(); ++i)
+    deallog << test_points[i] << " is inside: " <<
+            b.point_inside(test_points[i]) << std::endl;
+  deallog << std::endl;
+}
+
+void test_unitary()
+{
+  std::pair<Point<3>,Point<3>> boundaries;
+
+  for (int i=0; i<3; i++)
+    {
+      boundaries.second[i] = 1.0;
+    }
+
+  BoundingBox<3> b(boundaries);
+  deallog << "Boundary points:" << std::endl;
+  deallog << b.get_boundary_points().first << std::endl;
+  deallog << b.get_boundary_points().second << std::endl;
+
+  Point<3> p1(1.0,0,0);
+  Point<3> p2(0,1.0,0);
+  Point<3> p3(0,0,1.0);
+
+  deallog << "Checking if all vertices are inside: "
+          << b.point_inside(boundaries.first) << " "
+          << b.point_inside(boundaries.second) << std::endl;
+
+  deallog << b.point_inside(p1) << " "
+          << b.point_inside(p2) << " "
+          << b.point_inside(p3) << " "
+          << b.point_inside(p1+p2) << " "
+          << b.point_inside(p2+p3) << " "
+          << b.point_inside(p1+p3) << " " << std::endl;
+}
+
+int main()
+{
+  initlog();
+
+  deallog << "Test: Bounding Box class " << std::endl;
+  deallog << "Test for dimension 1" << std::endl;
+  test_bounding_box<1>();
+
+  deallog << std::endl << "Test for dimension 2" << std::endl;
+  test_bounding_box<2>();
+
+  deallog << std::endl << "Test for dimension 3" << std::endl;
+  test_bounding_box<3>();
+
+  deallog << std::endl << "Test for dimension 3, unitary box" << std::endl;
+  test_unitary();
+}
diff --git a/tests/base/bounding_box_1.output b/tests/base/bounding_box_1.output
new file mode 100644 (file)
index 0000000..b44911a
--- /dev/null
@@ -0,0 +1,77 @@
+
+DEAL::Test: Bounding Box class 
+DEAL::Test for dimension 1
+DEAL::Empty constructor: 
+DEAL::0.00000
+DEAL::0.00000
+DEAL::Boundary points: 
+DEAL::0.200000
+DEAL::0.800000
+DEAL::Boundary points are inside: 1 1
+DEAL::Points inside: 
+DEAL::0.680000 is inside: 1
+DEAL::0.560000 is inside: 1
+DEAL::0.440000 is inside: 1
+DEAL::0.320000 is inside: 1
+DEAL::0.200000 is inside: 1
+DEAL::
+DEAL::Points outside:
+DEAL::-2.50000 is inside: 0
+DEAL::5.00000 is inside: 0
+DEAL::-7.50000 is inside: 0
+DEAL::10.0000 is inside: 0
+DEAL::-12.5000 is inside: 0
+DEAL::
+DEAL::
+DEAL::Test for dimension 2
+DEAL::Empty constructor: 
+DEAL::0.00000 0.00000
+DEAL::0.00000 0.00000
+DEAL::Boundary points: 
+DEAL::0.200000 0.00000
+DEAL::0.800000 1.60000
+DEAL::Boundary points are inside: 1 1
+DEAL::Points inside: 
+DEAL::0.680000 1.28000 is inside: 1
+DEAL::0.560000 0.960000 is inside: 1
+DEAL::0.440000 0.640000 is inside: 1
+DEAL::0.320000 0.320000 is inside: 1
+DEAL::0.200000 0.00000 is inside: 1
+DEAL::
+DEAL::Points outside:
+DEAL::-2.50000 -4.00000 is inside: 0
+DEAL::5.00000 8.00000 is inside: 0
+DEAL::-7.50000 -12.0000 is inside: 0
+DEAL::10.0000 16.0000 is inside: 0
+DEAL::-12.5000 -20.0000 is inside: 0
+DEAL::
+DEAL::
+DEAL::Test for dimension 3
+DEAL::Empty constructor: 
+DEAL::0.00000 0.00000 0.00000
+DEAL::0.00000 0.00000 0.00000
+DEAL::Boundary points: 
+DEAL::0.200000 0.00000 -0.200000
+DEAL::0.800000 1.60000 2.40000
+DEAL::Boundary points are inside: 1 1
+DEAL::Points inside: 
+DEAL::0.680000 1.28000 1.88000 is inside: 1
+DEAL::0.560000 0.960000 1.36000 is inside: 1
+DEAL::0.440000 0.640000 0.840000 is inside: 1
+DEAL::0.320000 0.320000 0.320000 is inside: 1
+DEAL::0.200000 0.00000 -0.200000 is inside: 1
+DEAL::
+DEAL::Points outside:
+DEAL::-2.50000 -4.00000 -5.50000 is inside: 0
+DEAL::5.00000 8.00000 11.0000 is inside: 0
+DEAL::-7.50000 -12.0000 -16.5000 is inside: 0
+DEAL::10.0000 16.0000 22.0000 is inside: 0
+DEAL::-12.5000 -20.0000 -27.5000 is inside: 0
+DEAL::
+DEAL::
+DEAL::Test for dimension 3, unitary box
+DEAL::Boundary points:
+DEAL::0.00000 0.00000 0.00000
+DEAL::1.00000 1.00000 1.00000
+DEAL::Checking if all vertices are inside: 1 1
+DEAL::1 1 1 1 1 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.