]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added BoundingBox tree to Cache.
authorLuca Heltai <luca.heltai@sissa.it>
Mon, 3 Dec 2018 13:21:40 +0000 (14:21 +0100)
committerLuca Heltai <luca.heltai@sissa.it>
Mon, 3 Dec 2018 13:24:44 +0000 (14:24 +0100)
include/deal.II/grid/grid_tools_cache.h
include/deal.II/grid/grid_tools_cache_update_flags.h
source/grid/grid_tools_cache.cc
tests/boost/cache_boxes_rtree_01.cc [new file with mode: 0644]
tests/boost/cache_boxes_rtree_01.output [new file with mode: 0644]

index f96b15ccaecb3e6dd2954113d6f4fd912f1843a1..3a10319374f44c53fe322dcbd9dc351077baee9b 100644 (file)
@@ -129,6 +129,15 @@ namespace GridTools
     const RTree<std::pair<Point<spacedim>, unsigned int>> &
     get_used_vertices_rtree() const;
 
+    /**
+     * Returned the cached RTree object of the cell bounging boxes, constructed
+     * using the active cell iterators of the stored triangulation.
+     */
+    const RTree<
+      std::pair<BoundingBox<spacedim>,
+                typename Triangulation<dim, spacedim>::active_cell_iterator>> &
+    get_cell_bounding_boxes_rtree() const;
+
     /**
      * Return a reference to the stored triangulation.
      */
@@ -200,6 +209,15 @@ namespace GridTools
      */
     mutable RTree<std::pair<Point<spacedim>, unsigned int>> used_vertices_rtree;
 
+    /**
+     * Store an RTree object, containing the bounding boxes of the cells of the
+     * triangulation.
+     */
+    mutable RTree<
+      std::pair<BoundingBox<spacedim>,
+                typename Triangulation<dim, spacedim>::active_cell_iterator>>
+      cell_bounding_boxes_rtree;
+
     /**
      * Storage for the status of the triangulation signal.
      */
index 27885d44a4f75a6052de11844df7482486b8aea2..326f15d3ffeb8625aac02da0ba1e6bdf0d081163 100644 (file)
@@ -68,6 +68,11 @@ namespace GridTools
      */
     update_used_vertices_rtree = 0x010,
 
+    /**
+     * Update an RTree of the cell bounding boxes.
+     */
+    update_cell_bounding_boxes_rtree = 0x020,
+
     /**
      * Update all objects.
      */
index 03172def928f09f7071e6cf232d524f07141ab47..a969afe419ddf1159416d964178ed027dd2f3c54 100644 (file)
@@ -114,6 +114,26 @@ namespace GridTools
 
 
 
+  template <int dim, int spacedim>
+  const RTree<
+    std::pair<BoundingBox<spacedim>,
+              typename Triangulation<dim, spacedim>::active_cell_iterator>> &
+  Cache<dim, spacedim>::get_cell_bounding_boxes_rtree() const
+  {
+    if (update_flags & update_cell_bounding_boxes_rtree)
+      {
+        std::vector<std::pair<
+          BoundingBox<spacedim>,
+          typename Triangulation<dim, spacedim>::active_cell_iterator>>
+                     boxes(tria->n_active_cells());
+        unsigned int i = 0;
+        for (auto cell : tria->active_cell_iterators())
+          boxes[i++] = std::make_pair(cell->bounding_box(), cell);
+        cell_bounding_boxes_rtree = pack_rtree(boxes);
+      }
+    return cell_bounding_boxes_rtree;
+  }
+
 #ifdef DEAL_II_WITH_NANOFLANN
   template <int dim, int spacedim>
   const KDTree<spacedim> &
diff --git a/tests/boost/cache_boxes_rtree_01.cc b/tests/boost/cache_boxes_rtree_01.cc
new file mode 100644 (file)
index 0000000..7ec73a3
--- /dev/null
@@ -0,0 +1,76 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+// Extract cell bounding boxes rtree from the cache, and try to use it
+
+#include <deal.II/base/patterns.h>
+
+#include <deal.II/boost_adaptors/point.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools_cache.h>
+#include <deal.II/grid/tria.h>
+
+#include <algorithm>
+
+#include "../tests.h"
+
+using Patterns::Tools::to_string;
+
+namespace bg  = boost::geometry;
+namespace bgi = boost::geometry::index;
+
+template <int dim, int spacedim>
+void
+test(const unsigned int ref = 2, const unsigned int n_points = 10)
+{
+  Triangulation<dim, spacedim> tria;
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(ref);
+  GridTools::Cache cache(tria);
+
+  const auto &b_tree = cache.get_cell_bounding_boxes_rtree();
+
+  std::vector<Point<spacedim>> points(n_points);
+  std::generate(points.begin(), points.end(), []() {
+    return random_point<spacedim>();
+  });
+
+  deallog << "Testing dim = " << dim << ", spacedim = " << spacedim
+          << std::endl;
+
+  for (const auto &p : points)
+    {
+      std::vector<
+        std::pair<BoundingBox<spacedim>,
+                  typename Triangulation<dim, spacedim>::active_cell_iterator>>
+        res;
+      b_tree.query(bgi::nearest(p, 1), std::back_inserter(res));
+      deallog << "Nearest cell to " << p << ":  " << res[0].second << std::endl;
+    }
+}
+
+int
+main(int argc, char **argv)
+{
+  initlog();
+
+  test<1, 1>();
+  test<1, 2>();
+  test<1, 3>();
+  test<2, 2>();
+  test<2, 3>();
+  test<3, 3>();
+}
diff --git a/tests/boost/cache_boxes_rtree_01.output b/tests/boost/cache_boxes_rtree_01.output
new file mode 100644 (file)
index 0000000..32864e7
--- /dev/null
@@ -0,0 +1,67 @@
+
+DEAL::Testing dim = 1, spacedim = 1
+DEAL::Nearest cell to 0.840188:  2.3
+DEAL::Nearest cell to 0.394383:  2.1
+DEAL::Nearest cell to 0.783099:  2.3
+DEAL::Nearest cell to 0.798440:  2.3
+DEAL::Nearest cell to 0.911647:  2.3
+DEAL::Nearest cell to 0.197551:  2.0
+DEAL::Nearest cell to 0.335223:  2.1
+DEAL::Nearest cell to 0.768230:  2.3
+DEAL::Nearest cell to 0.277775:  2.1
+DEAL::Nearest cell to 0.553970:  2.2
+DEAL::Testing dim = 1, spacedim = 2
+DEAL::Nearest cell to 0.477397 0.628871:  2.1
+DEAL::Nearest cell to 0.364784 0.513401:  2.1
+DEAL::Nearest cell to 0.952230 0.916195:  2.3
+DEAL::Nearest cell to 0.635712 0.717297:  2.2
+DEAL::Nearest cell to 0.141603 0.606969:  2.0
+DEAL::Nearest cell to 0.0163006 0.242887:  2.0
+DEAL::Nearest cell to 0.137232 0.804177:  2.0
+DEAL::Nearest cell to 0.156679 0.400944:  2.0
+DEAL::Nearest cell to 0.129790 0.108809:  2.0
+DEAL::Nearest cell to 0.998925 0.218257:  2.3
+DEAL::Testing dim = 1, spacedim = 3
+DEAL::Nearest cell to 0.512932 0.839112 0.612640:  2.2
+DEAL::Nearest cell to 0.296032 0.637552 0.524287:  2.1
+DEAL::Nearest cell to 0.493583 0.972775 0.292517:  2.1
+DEAL::Nearest cell to 0.771358 0.526745 0.769914:  2.3
+DEAL::Nearest cell to 0.400229 0.891529 0.283315:  2.1
+DEAL::Nearest cell to 0.352458 0.807725 0.919026:  2.1
+DEAL::Nearest cell to 0.0697553 0.949327 0.525995:  2.0
+DEAL::Nearest cell to 0.0860558 0.192214 0.663227:  2.0
+DEAL::Nearest cell to 0.890233 0.348893 0.0641713:  2.3
+DEAL::Nearest cell to 0.0200230 0.457702 0.0630958:  2.0
+DEAL::Testing dim = 2, spacedim = 2
+DEAL::Nearest cell to 0.238280 0.970634:  2.10
+DEAL::Nearest cell to 0.902208 0.850920:  2.15
+DEAL::Nearest cell to 0.266666 0.539760:  2.9
+DEAL::Nearest cell to 0.375207 0.760249:  2.11
+DEAL::Nearest cell to 0.512535 0.667724:  2.12
+DEAL::Nearest cell to 0.531606 0.0392803:  2.4
+DEAL::Nearest cell to 0.437638 0.931835:  2.11
+DEAL::Nearest cell to 0.930810 0.720952:  2.13
+DEAL::Nearest cell to 0.284293 0.738534:  2.9
+DEAL::Nearest cell to 0.639979 0.354049:  2.6
+DEAL::Testing dim = 2, spacedim = 3
+DEAL::Nearest cell to 0.687861 0.165974 0.440105:  2.4
+DEAL::Nearest cell to 0.880075 0.829201 0.330337:  2.15
+DEAL::Nearest cell to 0.228968 0.893372 0.350360:  2.10
+DEAL::Nearest cell to 0.686670 0.956468 0.588640:  2.14
+DEAL::Nearest cell to 0.657304 0.858676 0.439560:  2.14
+DEAL::Nearest cell to 0.923970 0.398437 0.814767:  2.7
+DEAL::Nearest cell to 0.684219 0.910972 0.482491:  2.14
+DEAL::Nearest cell to 0.215825 0.950252 0.920128:  2.10
+DEAL::Nearest cell to 0.147660 0.881062 0.641081:  2.10
+DEAL::Nearest cell to 0.431953 0.619596 0.281059:  2.9
+DEAL::Testing dim = 3, spacedim = 3
+DEAL::Nearest cell to 0.786002 0.307458 0.447034:  2.15
+DEAL::Nearest cell to 0.226107 0.187533 0.276235:  2.4
+DEAL::Nearest cell to 0.556444 0.416501 0.169607:  2.10
+DEAL::Nearest cell to 0.906804 0.103171 0.126075:  2.9
+DEAL::Nearest cell to 0.495444 0.760475 0.984752:  2.55
+DEAL::Nearest cell to 0.935004 0.684445 0.383188:  2.29
+DEAL::Nearest cell to 0.749771 0.368664 0.294160:  2.14
+DEAL::Nearest cell to 0.232262 0.584489 0.244413:  2.16
+DEAL::Nearest cell to 0.152390 0.732149 0.125475:  2.16
+DEAL::Nearest cell to 0.793470 0.164102 0.745071:  2.41
\ No newline at end of file

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.