]> https://gitweb.dealii.org/ - dealii.git/commitdiff
GridTools::find_closest_vertex(map, point)
authorLuca Heltai <luca.heltai@sissa.it>
Tue, 19 Sep 2017 12:49:28 +0000 (14:49 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Wed, 27 Sep 2017 13:07:33 +0000 (15:07 +0200)
doc/news/changes/minor/20170919LucaHeltai [new file with mode: 0644]
include/deal.II/grid/grid_tools.h
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
tests/grid/grid_tools_extract_used_vertices_02.cc [new file with mode: 0644]
tests/grid/grid_tools_extract_used_vertices_02.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20170919LucaHeltai b/doc/news/changes/minor/20170919LucaHeltai
new file mode 100644 (file)
index 0000000..e406b0a
--- /dev/null
@@ -0,0 +1,6 @@
+New: A GridTools::find_closest_vertex() has been added
+that takes a map of index-vertex position and a target 
+point, and returns the index of the closest point in 
+the map. 
+<br>
+(Luca Heltai, 2017/09/19)
index 4585e16b59a1d729ff4af2e7c5f671af0d4ab340..f6fdc81b454fbc412aa8b4d29cd175f9b889fda1 100644 (file)
@@ -593,6 +593,19 @@ namespace GridTools
                                         const Triangulation<dim,spacedim> &container,
                                         const Mapping<dim,spacedim> &mapping = StaticMappingQ1<dim,spacedim>::mapping);
 
+  /**
+   * Find and return the number of the closest vertex to a given point in the
+   * map of vertices passed as the first argument.
+   *
+   * @param vertices A map of index->vertex, as returned by
+   *        extract_used_vertices()
+   * @param p The target point
+   * @return The index of the vertex that is closest to the target point `p`
+   */
+  template<int spacedim>
+  unsigned int
+  find_closest_vertex (const std::map<unsigned int,Point<spacedim>> &vertices,
+                       const Point<spacedim>         &p);
 
   /**
    * Find and return the number of the used vertex (or marked vertex) in a
index 1229339a78c35a13b85fb4a29a19bd53a3d0f4e7..cac93c628a8d5ab37e2497bf21a4ffc7eee260c5 100644 (file)
@@ -4943,6 +4943,21 @@ next_cell:
     return result;
   }
 
+  template<int spacedim>
+  unsigned int
+  find_closest_vertex(const std::map<unsigned int,Point<spacedim> > &vertices,
+                      const Point<spacedim> &p)
+  {
+    auto id_and_v = std::min_element(vertices.begin(), vertices.end(),
+                                     [&](const std::pair<const unsigned int, Point<spacedim>> &p1,
+                                         const std::pair<const unsigned int, Point<spacedim>> &p2) -> bool
+    {
+      return p1.second.distance(p) < p2.second.distance(p);
+    }
+                                    );
+    return id_and_v->first;
+  }
+
 } /* namespace GridTools */
 
 
index 2fc9db2041c65c7e62cfece2bafff8bd6dcfd74d..3938a71118085483d6b6117653c3ade982490314 100644 (file)
@@ -143,6 +143,10 @@ for (deal_II_space_dimension : SPACE_DIMENSIONS)
     find_active_cell_around_point (const Mapping<deal_II_space_dimension> &,
                                    const parallel::distributed::Triangulation<deal_II_space_dimension> &,
                                    const Point<deal_II_space_dimension> &);
+
+    template unsigned int
+    GridTools::find_closest_vertex(const std::map<unsigned int,Point<deal_II_space_dimension> >& vertices,
+                                   const Point<deal_II_space_dimension>& p);
 }
 
 
diff --git a/tests/grid/grid_tools_extract_used_vertices_02.cc b/tests/grid/grid_tools_extract_used_vertices_02.cc
new file mode 100644 (file)
index 0000000..50f3374
--- /dev/null
@@ -0,0 +1,62 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2001 - 2015 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.
+//
+// ---------------------------------------------------------------------
+
+// Check extract used_vertices and find_closest_vertex
+
+#include "../tests.h"
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/grid_out.h>
+
+
+template <int dim, int spacedim>
+void test (const Point<spacedim> &p)
+{
+  deallog << "dim: " << dim << ", spacedim: "
+          << spacedim << std::endl;
+
+  Triangulation<dim,spacedim> tria;
+  GridGenerator::hyper_cube(tria);
+
+  tria.refine_global(1);
+  tria.begin_active()->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+
+  for (auto cell = tria.begin_active(2); cell != tria.end(2); ++cell)
+    cell->set_coarsen_flag();
+
+  tria.execute_coarsening_and_refinement();
+
+  auto m = GridTools::extract_used_vertices(tria);
+
+  for (auto &e: m)
+    deallog << "Vertex: " << e.first << ": " << e.second << std::endl;
+
+  auto i = GridTools::find_closest_vertex(m,p);
+  deallog << "Closest vertex to " << p
+          << ", v["<< i << "] :"
+          << m[i] << std::endl;
+};
+
+
+int main ()
+{
+  initlog();
+  test<2,2> (Point<2>(.2,.2));
+  test<2,2> (Point<2>(.6,.9));
+  return 0;
+}
+
diff --git a/tests/grid/grid_tools_extract_used_vertices_02.output b/tests/grid/grid_tools_extract_used_vertices_02.output
new file mode 100644 (file)
index 0000000..de11874
--- /dev/null
@@ -0,0 +1,23 @@
+
+DEAL::dim: 2, spacedim: 2
+DEAL::Vertex: 0: 0.00000 0.00000
+DEAL::Vertex: 1: 1.00000 0.00000
+DEAL::Vertex: 2: 0.00000 1.00000
+DEAL::Vertex: 3: 1.00000 1.00000
+DEAL::Vertex: 4: 0.500000 0.00000
+DEAL::Vertex: 5: 0.00000 0.500000
+DEAL::Vertex: 6: 1.00000 0.500000
+DEAL::Vertex: 7: 0.500000 1.00000
+DEAL::Vertex: 8: 0.500000 0.500000
+DEAL::Closest vertex to 0.200000 0.200000, v[0] :0.00000 0.00000
+DEAL::dim: 2, spacedim: 2
+DEAL::Vertex: 0: 0.00000 0.00000
+DEAL::Vertex: 1: 1.00000 0.00000
+DEAL::Vertex: 2: 0.00000 1.00000
+DEAL::Vertex: 3: 1.00000 1.00000
+DEAL::Vertex: 4: 0.500000 0.00000
+DEAL::Vertex: 5: 0.00000 0.500000
+DEAL::Vertex: 6: 1.00000 0.500000
+DEAL::Vertex: 7: 0.500000 1.00000
+DEAL::Vertex: 8: 0.500000 0.500000
+DEAL::Closest vertex to 0.600000 0.900000, v[7] :0.500000 1.00000

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.