--- /dev/null
+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)
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
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 */
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);
}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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;
+}
+
--- /dev/null
+
+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