]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add tests for new BVH
authorBruno Turcksin <bruno.turcksin@gmail.com>
Fri, 8 Jan 2021 20:48:43 +0000 (20:48 +0000)
committerBruno Turcksin <bruno.turcksin@gmail.com>
Fri, 12 Feb 2021 14:21:38 +0000 (14:21 +0000)
tests/arborx/CMakeLists.txt [new file with mode: 0644]
tests/arborx/bounding_box_intersect.cc [new file with mode: 0644]
tests/arborx/bounding_box_intersect.output [new file with mode: 0644]
tests/arborx/point_intersect.cc [new file with mode: 0644]
tests/arborx/point_intersect.output [new file with mode: 0644]
tests/run_testsuite.cmake

diff --git a/tests/arborx/CMakeLists.txt b/tests/arborx/CMakeLists.txt
new file mode 100644 (file)
index 0000000..288e920
--- /dev/null
@@ -0,0 +1,6 @@
+CMAKE_MINIMUM_REQUIRED(VERSION 2.8.9)
+INCLUDE(../setup_testsubproject.cmake)
+PROJECT(testsuite CXX)
+IF(DEAL_II_WITH_ARBORX)
+  DEAL_II_PICKUP_TESTS()
+ENDIF()
diff --git a/tests/arborx/bounding_box_intersect.cc b/tests/arborx/bounding_box_intersect.cc
new file mode 100644 (file)
index 0000000..741d727
--- /dev/null
@@ -0,0 +1,186 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 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.
+//
+// ---------------------------------------------------------------------
+
+
+// Check ArborX wrapper: intersection of boundary boxes
+
+
+#include <deal.II/arborx/access_traits.h>
+#include <deal.II/arborx/bvh.h>
+
+#include <deal.II/base/bounding_box.h>
+
+#include "../tests.h"
+
+
+void
+test_2d()
+{
+  std::vector<BoundingBox<2>> bounding_boxes;
+  std::vector<Point<2>>       points;
+
+  unsigned int n_points_1d = 5;
+  for (unsigned int i = 0; i < n_points_1d; ++i)
+    {
+      for (unsigned int j = 0; j < n_points_1d; ++j)
+        {
+          points.emplace_back(j, i);
+        }
+    }
+
+  for (unsigned int i = 0; i < n_points_1d - 1; ++i)
+    {
+      for (unsigned int j = 0; j < n_points_1d - 1; ++j)
+        {
+          unsigned int point_index = j + i * n_points_1d;
+          bounding_boxes.push_back(
+            std::make_pair(points[point_index],
+                           points[point_index + n_points_1d + 1]));
+        }
+    }
+
+  Point<2> point_a(0.5, 0.5);
+  Point<2> point_b(1.5, 1.5);
+  Point<2> point_c(2.2, 2.2);
+  Point<2> point_d(2.6, 2.6);
+
+  std::vector<BoundingBox<2>> query_bounding_boxes;
+  query_bounding_boxes.emplace_back(std::make_pair(point_a, point_b));
+  query_bounding_boxes.emplace_back(std::make_pair(point_c, point_d));
+
+  ArborXWrappers::BVH                  bvh(bounding_boxes);
+  ArborXWrappers::BoundingBoxIntersect bb_intersect(query_bounding_boxes);
+  auto                                 indices_offset = bvh.query(bb_intersect);
+  std::vector<int>                     indices        = indices_offset.first;
+  std::vector<int>                     offset         = indices_offset.second;
+
+  std::vector<int> indices_ref = {0, 1, 4, 5, 10};
+  std::vector<int> offset_ref  = {0, 4, 5};
+
+  AssertThrow(indices.size() == indices_ref.size(), ExcInternalError());
+  AssertThrow(offset.size() == offset_ref.size(), ExcInternalError());
+  for (unsigned int i = 0; i < offset.size() - 1; ++i)
+    {
+      for (unsigned int j = offset[i]; j < offset[i + 1]; ++j)
+        {
+          // The indices associated to each query are not ordered.
+          bool found = false;
+          for (unsigned int k = offset[i]; k < offset[i + 1]; ++k)
+            {
+              if (indices[j] == indices_ref[k])
+                {
+                  found = true;
+                  break;
+                }
+            }
+          AssertThrow(found, ExcInternalError());
+        }
+    }
+  for (unsigned int i = 0; i < offset.size(); ++i)
+    AssertThrow(offset[i] == offset_ref[i], ExcInternalError());
+
+  deallog << "OK" << std::endl;
+}
+
+void
+test_3d()
+{
+  std::vector<BoundingBox<3>> bounding_boxes;
+  std::vector<Point<3>>       points;
+
+  unsigned int n_points_1d = 5;
+  for (unsigned int i = 0; i < n_points_1d; ++i)
+    {
+      for (unsigned int j = 0; j < n_points_1d; ++j)
+        {
+          for (unsigned int k = 0; k < n_points_1d; ++k)
+            {
+              points.emplace_back(k, j, i);
+            }
+        }
+    }
+
+  for (unsigned int i = 0; i < n_points_1d - 1; ++i)
+    {
+      for (unsigned int j = 0; j < n_points_1d - 1; ++j)
+        {
+          for (unsigned int k = 0; k < n_points_1d - 1; ++k)
+            {
+              unsigned int point_index =
+                k + j * n_points_1d + i * n_points_1d * n_points_1d;
+              bounding_boxes.push_back(
+                std::make_pair(points[point_index],
+                               points[point_index + n_points_1d * n_points_1d +
+                                      n_points_1d + 1]));
+            }
+        }
+    }
+
+  Point<3> point_a(0.5, 0.5, 0.5);
+  Point<3> point_b(1.5, 1.5, 1.5);
+  Point<3> point_c(2.2, 2.2, 2.2);
+  Point<3> point_d(2.6, 2.6, 2.6);
+
+  std::vector<BoundingBox<3>> query_bounding_boxes;
+  query_bounding_boxes.emplace_back(std::make_pair(point_a, point_b));
+  query_bounding_boxes.emplace_back(std::make_pair(point_c, point_d));
+
+  ArborXWrappers::BVH                  bvh(bounding_boxes);
+  ArborXWrappers::BoundingBoxIntersect bb_intersect(query_bounding_boxes);
+  auto                                 indices_offset = bvh.query(bb_intersect);
+  std::vector<int>                     indices        = indices_offset.first;
+  std::vector<int>                     offset         = indices_offset.second;
+
+  std::vector<int> indices_ref = {0, 1, 4, 5, 16, 17, 20, 21, 42};
+  std::vector<int> offset_ref  = {0, 8, 9};
+
+  AssertThrow(indices.size() == indices_ref.size(), ExcInternalError());
+  AssertThrow(offset.size() == offset_ref.size(), ExcInternalError());
+  for (unsigned int i = 0; i < offset.size() - 1; ++i)
+    {
+      for (unsigned int j = offset[i]; j < offset[i + 1]; ++j)
+        {
+          // The indices associated to each query are not ordered.
+          bool found = false;
+          for (unsigned int k = offset[i]; k < offset[i + 1]; ++k)
+            {
+              if (indices[j] == indices_ref[k])
+                {
+                  found = true;
+                  break;
+                }
+            }
+          AssertThrow(found, ExcInternalError());
+        }
+    }
+  for (unsigned int i = 0; i < offset.size(); ++i)
+    AssertThrow(offset[i] == offset_ref[i], ExcInternalError());
+  deallog << "OK" << std::endl;
+}
+
+int
+main(int argc, char **argv)
+{
+  // Initialize Kokkos
+  Kokkos::initialize(argc, argv);
+
+  initlog();
+
+  // tests
+  test_2d();
+  test_3d();
+
+  Kokkos::finalize();
+}
diff --git a/tests/arborx/bounding_box_intersect.output b/tests/arborx/bounding_box_intersect.output
new file mode 100644 (file)
index 0000000..8b3b075
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL::OK
+DEAL::OK
diff --git a/tests/arborx/point_intersect.cc b/tests/arborx/point_intersect.cc
new file mode 100644 (file)
index 0000000..6093454
--- /dev/null
@@ -0,0 +1,186 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 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.
+//
+// ---------------------------------------------------------------------
+
+
+// Check ArborX wrapper: intersection of points with bounding boxes
+
+
+#include <deal.II/arborx/access_traits.h>
+#include <deal.II/arborx/bvh.h>
+
+#include <deal.II/base/point.h>
+
+#include "../tests.h"
+
+
+void
+test_2d()
+{
+  std::vector<BoundingBox<2>> bounding_boxes;
+  std::vector<Point<2>>       points;
+
+  unsigned int n_points_1d = 5;
+  for (unsigned int i = 0; i < n_points_1d; ++i)
+    {
+      for (unsigned int j = 0; j < n_points_1d; ++j)
+        {
+          points.emplace_back(j, i);
+        }
+    }
+
+  for (unsigned int i = 0; i < n_points_1d - 1; ++i)
+    {
+      for (unsigned int j = 0; j < n_points_1d - 1; ++j)
+        {
+          unsigned int point_index = j + i * n_points_1d;
+          bounding_boxes.push_back(
+            std::make_pair(points[point_index],
+                           points[point_index + n_points_1d + 1]));
+        }
+    }
+
+  std::vector<Point<2>> query_points;
+  query_points.emplace_back(0.5, 0.5);
+  query_points.emplace_back(1.5, 1.5);
+  query_points.emplace_back(2.2, 2.2);
+  query_points.emplace_back(2.6, 2.6);
+
+
+  ArborXWrappers::BVH            bvh(bounding_boxes);
+  ArborXWrappers::PointIntersect pt_intersect(query_points);
+  auto                           indices_offset = bvh.query(pt_intersect);
+  std::vector<int>               indices        = indices_offset.first;
+  std::vector<int>               offset         = indices_offset.second;
+
+  std::vector<int> indices_ref = {0, 5, 10, 10};
+  std::vector<int> offset_ref  = {0, 1, 2, 3, 4};
+
+  AssertThrow(indices.size() == indices_ref.size(), ExcInternalError());
+  AssertThrow(offset.size() == offset_ref.size(), ExcInternalError());
+  for (unsigned int i = 0; i < offset.size() - 1; ++i)
+    {
+      for (unsigned int j = offset[i]; j < offset[i + 1]; ++j)
+        {
+          // The indices associated to each query are not ordered.
+          bool found = false;
+          for (unsigned int k = offset[i]; k < offset[i + 1]; ++k)
+            {
+              if (indices[j] == indices_ref[k])
+                {
+                  found = true;
+                  break;
+                }
+            }
+          AssertThrow(found, ExcInternalError());
+        }
+    }
+  for (unsigned int i = 0; i < offset.size(); ++i)
+    AssertThrow(offset[i] == offset_ref[i], ExcInternalError());
+
+  deallog << "OK" << std::endl;
+}
+
+void
+test_3d()
+{
+  std::vector<BoundingBox<3>> bounding_boxes;
+  std::vector<Point<3>>       points;
+
+  unsigned int n_points_1d = 5;
+  for (unsigned int i = 0; i < n_points_1d; ++i)
+    {
+      for (unsigned int j = 0; j < n_points_1d; ++j)
+        {
+          for (unsigned int k = 0; k < n_points_1d; ++k)
+            {
+              points.emplace_back(k, j, i);
+            }
+        }
+    }
+
+  for (unsigned int i = 0; i < n_points_1d - 1; ++i)
+    {
+      for (unsigned int j = 0; j < n_points_1d - 1; ++j)
+        {
+          for (unsigned int k = 0; k < n_points_1d - 1; ++k)
+            {
+              unsigned int point_index =
+                k + j * n_points_1d + i * n_points_1d * n_points_1d;
+              bounding_boxes.push_back(
+                std::make_pair(points[point_index],
+                               points[point_index + n_points_1d * n_points_1d +
+                                      n_points_1d + 1]));
+            }
+        }
+    }
+
+  std::vector<Point<3>> query_points;
+  query_points.emplace_back(0.5, 0.5, 0.5);
+  query_points.emplace_back(1.5, 1.5, 1.5);
+  query_points.emplace_back(2.2, 2.2, 2.2);
+  query_points.emplace_back(2.6, 2.6, 2.6);
+
+
+  ArborXWrappers::BVH            bvh(bounding_boxes);
+  ArborXWrappers::PointIntersect pt_intersect(query_points);
+  auto                           indices_offset = bvh.query(pt_intersect);
+  std::vector<int>               indices        = indices_offset.first;
+  std::vector<int>               offset         = indices_offset.second;
+
+
+  std::vector<int> indices_ref = {0, 21, 42, 42};
+  std::vector<int> offset_ref  = {0, 1, 2, 3, 4};
+
+  AssertThrow(indices.size() == indices_ref.size(), ExcInternalError());
+  AssertThrow(offset.size() == offset_ref.size(), ExcInternalError());
+  for (unsigned int i = 0; i < offset.size() - 1; ++i)
+    {
+      for (unsigned int j = offset[i]; j < offset[i + 1]; ++j)
+        {
+          // The indices associated to each query are not ordered.
+          bool found = false;
+          for (unsigned int k = offset[i]; k < offset[i + 1]; ++k)
+            {
+              if (indices[j] == indices_ref[k])
+                {
+                  found = true;
+                  break;
+                }
+            }
+          AssertThrow(found, ExcInternalError());
+        }
+    }
+  for (unsigned int i = 0; i < offset.size(); ++i)
+    AssertThrow(offset[i] == offset_ref[i], ExcInternalError());
+  deallog << "OK" << std::endl;
+}
+
+int
+main(int argc, char **argv)
+{
+  // Initialize Kokkos
+  Kokkos::initialize(argc, argv);
+
+  initlog();
+
+  // Initialize ArborX
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv);
+
+  // tests
+  test_2d();
+  test_3d();
+
+  Kokkos::finalize();
+}
diff --git a/tests/arborx/point_intersect.output b/tests/arborx/point_intersect.output
new file mode 100644 (file)
index 0000000..8b3b075
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL::OK
+DEAL::OK
index 5106248c67763a5d1afba793a2ed11fae51ab931..39f6c57c038f1a2ff63e0f9ce9e49838589936a8 100644 (file)
@@ -272,7 +272,7 @@ GET_CMAKE_PROPERTY(_variables VARIABLES)
 FOREACH(_var ${_variables})
   IF( _var MATCHES "^(TEST|DEAL_II|ALLOW|WITH|FORCE|COMPONENT)_" OR
       _var MATCHES "^(DOCUMENTATION|EXAMPLES)" OR
-      _var MATCHES "^(ADOLC|ARPACK|BOOST|OPENCASCADE|MUPARSER|HDF5|METIS|MPI)_" OR
+      _var MATCHES "^(ADOLC|ARBORX|ARPACK|BOOST|OPENCASCADE|MUPARSER|HDF5|KOKKOS|METIS|MPI)_" OR
       _var MATCHES "^(GINKGO|P4EST|PETSC|SCALAPACK|SLEPC|THREADS|TBB|TRILINOS)_" OR
       _var MATCHES "^(UMFPACK|ZLIB|LAPACK|MUPARSER|CUDA)_" OR
       _var MATCHES "^(CMAKE|DEAL_II)_(C|CXX|Fortran|BUILD)_(COMPILER|FLAGS)" OR

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.