From 282df41d84c793f3da0b0a7071fcc07854dd4f7f Mon Sep 17 00:00:00 2001 From: Luca Heltai Date: Mon, 3 Dec 2018 14:14:43 +0100 Subject: [PATCH] Added rtree of vertices to Cache. --- include/deal.II/grid/grid_tools_cache.h | 13 ++++ .../grid/grid_tools_cache_update_flags.h | 17 +++-- source/grid/grid_tools_cache.cc | 19 +++++ tests/boost/cache_vertices_rtree_01.cc | 74 +++++++++++++++++++ tests/boost/cache_vertices_rtree_01.output | 67 +++++++++++++++++ 5 files changed, 184 insertions(+), 6 deletions(-) create mode 100644 tests/boost/cache_vertices_rtree_01.cc create mode 100644 tests/boost/cache_vertices_rtree_01.output diff --git a/include/deal.II/grid/grid_tools_cache.h b/include/deal.II/grid/grid_tools_cache.h index 6b64b1e866..f96b15ccae 100644 --- a/include/deal.II/grid/grid_tools_cache.h +++ b/include/deal.II/grid/grid_tools_cache.h @@ -31,6 +31,7 @@ #include #include +#include #include @@ -121,6 +122,13 @@ namespace GridTools const std::map> & get_used_vertices() const; + /** + * Returned the cached RTree object for the vertices, constructed using the + * used vertices of the triangulation. + */ + const RTree, unsigned int>> & + get_used_vertices_rtree() const; + /** * Return a reference to the stored triangulation. */ @@ -187,6 +195,11 @@ namespace GridTools */ mutable std::map> used_vertices; + /** + * Store an RTree object, containing the used vertices of the triangulation. + */ + mutable RTree, unsigned int>> used_vertices_rtree; + /** * Storage for the status of the triangulation signal. */ diff --git a/include/deal.II/grid/grid_tools_cache_update_flags.h b/include/deal.II/grid/grid_tools_cache_update_flags.h index a3c17c49f3..27885d44a4 100644 --- a/include/deal.II/grid/grid_tools_cache_update_flags.h +++ b/include/deal.II/grid/grid_tools_cache_update_flags.h @@ -37,36 +37,41 @@ namespace GridTools /** * Update Nothing. */ - update_nothing = 0x00, + update_nothing = 0x000, /** * Update vertex_to_cell_map, as returned by * GridTools::vertex_to_cell_map(). */ - update_vertex_to_cell_map = 0x01, + update_vertex_to_cell_map = 0x001, /** * Update vertex_to_cell_centers_directions, as returned by * GridTools::vertex_to_cell_centers_directions() */ update_vertex_to_cell_centers_directions = - update_vertex_to_cell_map | 0x0002, + update_vertex_to_cell_map | 0x002, /** * Update a KDTree object, initialized with the vertices of the * Triangulation. */ - update_vertex_kdtree = 0x04, + update_vertex_kdtree = 0x004, /** * Update a mapping of used vertices. */ - update_used_vertices = 0x08, + update_used_vertices = 0x008, + + /** + * Update an RTree of the used vertices. + */ + update_used_vertices_rtree = 0x010, /** * Update all objects. */ - update_all = 0xFF, + update_all = 0xFFF, }; diff --git a/source/grid/grid_tools_cache.cc b/source/grid/grid_tools_cache.cc index 05698a7cce..03172def92 100644 --- a/source/grid/grid_tools_cache.cc +++ b/source/grid/grid_tools_cache.cc @@ -95,6 +95,25 @@ namespace GridTools + template + const RTree, unsigned int>> & + Cache::get_used_vertices_rtree() const + { + if (update_flags & update_used_vertices_rtree) + { + const auto &used_vertices = get_used_vertices(); + std::vector, unsigned int>> vertices( + used_vertices.size()); + unsigned int i = 0; + for (auto it : used_vertices) + vertices[i++] = std::make_pair(it.second, it.first); + used_vertices_rtree = pack_rtree(vertices); + } + return used_vertices_rtree; + } + + + #ifdef DEAL_II_WITH_NANOFLANN template const KDTree & diff --git a/tests/boost/cache_vertices_rtree_01.cc b/tests/boost/cache_vertices_rtree_01.cc new file mode 100644 index 0000000000..c6ae8d56ae --- /dev/null +++ b/tests/boost/cache_vertices_rtree_01.cc @@ -0,0 +1,74 @@ +// --------------------------------------------------------------------- +// +// 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 vertex rtree from the cache, and try to use it + +#include + +#include + +#include +#include +#include + +#include + +#include "../tests.h" + +using Patterns::Tools::to_string; + +namespace bg = boost::geometry; +namespace bgi = boost::geometry::index; + +template +void +test(const unsigned int ref = 2, const unsigned int n_points = 10) +{ + Triangulation tria; + GridGenerator::hyper_cube(tria); + tria.refine_global(ref); + GridTools::Cache cache(tria); + + const auto &v_tree = cache.get_used_vertices_rtree(); + + std::vector> points(n_points); + std::generate(points.begin(), points.end(), []() { + return random_point(); + }); + + deallog << "Testing dim = " << dim << ", spacedim = " << spacedim + << std::endl; + + for (const auto &p : points) + { + std::vector, unsigned int>> res; + v_tree.query(bgi::nearest(p, 1), std::back_inserter(res)); + deallog << "Nearest vertex to " << p << ": v[" << res[0].second + << "] = " << res[0].first << 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_vertices_rtree_01.output b/tests/boost/cache_vertices_rtree_01.output new file mode 100644 index 0000000000..6e3b4833eb --- /dev/null +++ b/tests/boost/cache_vertices_rtree_01.output @@ -0,0 +1,67 @@ + +DEAL::Testing dim = 1, spacedim = 1 +DEAL::Nearest vertex to 0.840188: v[4] = 0.750000 +DEAL::Nearest vertex to 0.394383: v[2] = 0.500000 +DEAL::Nearest vertex to 0.783099: v[4] = 0.750000 +DEAL::Nearest vertex to 0.798440: v[4] = 0.750000 +DEAL::Nearest vertex to 0.911647: v[1] = 1.00000 +DEAL::Nearest vertex to 0.197551: v[3] = 0.250000 +DEAL::Nearest vertex to 0.335223: v[3] = 0.250000 +DEAL::Nearest vertex to 0.768230: v[4] = 0.750000 +DEAL::Nearest vertex to 0.277775: v[3] = 0.250000 +DEAL::Nearest vertex to 0.553970: v[2] = 0.500000 +DEAL::Testing dim = 1, spacedim = 2 +DEAL::Nearest vertex to 0.477397 0.628871: v[2] = 0.500000 0.00000 +DEAL::Nearest vertex to 0.364784 0.513401: v[3] = 0.250000 0.00000 +DEAL::Nearest vertex to 0.952230 0.916195: v[1] = 1.00000 0.00000 +DEAL::Nearest vertex to 0.635712 0.717297: v[4] = 0.750000 0.00000 +DEAL::Nearest vertex to 0.141603 0.606969: v[3] = 0.250000 0.00000 +DEAL::Nearest vertex to 0.0163006 0.242887: v[0] = 0.00000 0.00000 +DEAL::Nearest vertex to 0.137232 0.804177: v[3] = 0.250000 0.00000 +DEAL::Nearest vertex to 0.156679 0.400944: v[3] = 0.250000 0.00000 +DEAL::Nearest vertex to 0.129790 0.108809: v[3] = 0.250000 0.00000 +DEAL::Nearest vertex to 0.998925 0.218257: v[1] = 1.00000 0.00000 +DEAL::Testing dim = 1, spacedim = 3 +DEAL::Nearest vertex to 0.512932 0.839112 0.612640: v[2] = 0.500000 0.00000 0.00000 +DEAL::Nearest vertex to 0.296032 0.637552 0.524287: v[3] = 0.250000 0.00000 0.00000 +DEAL::Nearest vertex to 0.493583 0.972775 0.292517: v[2] = 0.500000 0.00000 0.00000 +DEAL::Nearest vertex to 0.771358 0.526745 0.769914: v[4] = 0.750000 0.00000 0.00000 +DEAL::Nearest vertex to 0.400229 0.891529 0.283315: v[2] = 0.500000 0.00000 0.00000 +DEAL::Nearest vertex to 0.352458 0.807725 0.919026: v[3] = 0.250000 0.00000 0.00000 +DEAL::Nearest vertex to 0.0697553 0.949327 0.525995: v[0] = 0.00000 0.00000 0.00000 +DEAL::Nearest vertex to 0.0860558 0.192214 0.663227: v[0] = 0.00000 0.00000 0.00000 +DEAL::Nearest vertex to 0.890233 0.348893 0.0641713: v[1] = 1.00000 0.00000 0.00000 +DEAL::Nearest vertex to 0.0200230 0.457702 0.0630958: v[0] = 0.00000 0.00000 0.00000 +DEAL::Testing dim = 2, spacedim = 2 +DEAL::Nearest vertex to 0.238280 0.970634: v[15] = 0.250000 1.00000 +DEAL::Nearest vertex to 0.902208 0.850920: v[14] = 1.00000 0.750000 +DEAL::Nearest vertex to 0.266666 0.539760: v[19] = 0.250000 0.500000 +DEAL::Nearest vertex to 0.375207 0.760249: v[18] = 0.500000 0.750000 +DEAL::Nearest vertex to 0.512535 0.667724: v[18] = 0.500000 0.750000 +DEAL::Nearest vertex to 0.531606 0.0392803: v[4] = 0.500000 0.00000 +DEAL::Nearest vertex to 0.437638 0.931835: v[7] = 0.500000 1.00000 +DEAL::Nearest vertex to 0.930810 0.720952: v[14] = 1.00000 0.750000 +DEAL::Nearest vertex to 0.284293 0.738534: v[23] = 0.250000 0.750000 +DEAL::Nearest vertex to 0.639979 0.354049: v[22] = 0.750000 0.250000 +DEAL::Testing dim = 2, spacedim = 3 +DEAL::Nearest vertex to 0.687861 0.165974 0.440105: v[22] = 0.750000 0.250000 0.00000 +DEAL::Nearest vertex to 0.880075 0.829201 0.330337: v[14] = 1.00000 0.750000 0.00000 +DEAL::Nearest vertex to 0.228968 0.893372 0.350360: v[15] = 0.250000 1.00000 0.00000 +DEAL::Nearest vertex to 0.686670 0.956468 0.588640: v[16] = 0.750000 1.00000 0.00000 +DEAL::Nearest vertex to 0.657304 0.858676 0.439560: v[24] = 0.750000 0.750000 0.00000 +DEAL::Nearest vertex to 0.923970 0.398437 0.814767: v[6] = 1.00000 0.500000 0.00000 +DEAL::Nearest vertex to 0.684219 0.910972 0.482491: v[16] = 0.750000 1.00000 0.00000 +DEAL::Nearest vertex to 0.215825 0.950252 0.920128: v[15] = 0.250000 1.00000 0.00000 +DEAL::Nearest vertex to 0.147660 0.881062 0.641081: v[15] = 0.250000 1.00000 0.00000 +DEAL::Nearest vertex to 0.431953 0.619596 0.281059: v[8] = 0.500000 0.500000 0.00000 +DEAL::Testing dim = 3, spacedim = 3 +DEAL::Nearest vertex to 0.786002 0.307458 0.447034: v[107] = 0.750000 0.250000 0.500000 +DEAL::Nearest vertex to 0.226107 0.187533 0.276235: v[117] = 0.250000 0.250000 0.250000 +DEAL::Nearest vertex to 0.556444 0.416501 0.169607: v[76] = 0.500000 0.500000 0.250000 +DEAL::Nearest vertex to 0.906804 0.103171 0.126075: v[35] = 1.00000 0.00000 0.250000 +DEAL::Nearest vertex to 0.495444 0.760475 0.984752: v[72] = 0.500000 0.750000 1.00000 +DEAL::Nearest vertex to 0.935004 0.684445 0.383188: v[66] = 1.00000 0.750000 0.500000 +DEAL::Nearest vertex to 0.749771 0.368664 0.294160: v[118] = 0.750000 0.250000 0.250000 +DEAL::Nearest vertex to 0.232262 0.584489 0.244413: v[112] = 0.250000 0.500000 0.250000 +DEAL::Nearest vertex to 0.152390 0.732149 0.125475: v[119] = 0.250000 0.750000 0.250000 +DEAL::Nearest vertex to 0.793470 0.164102 0.745071: v[122] = 0.750000 0.250000 0.750000 \ No newline at end of file -- 2.39.5