--- /dev/null
+New: A GridTools::extract_used_vertices() has been added
+that returns a map of index-vertex position, using an
+optional Mapping<dim,spacedim> object to compute the
+vertex locations.
+<br>
+(Luca Heltai, 2017/09/18)
*/
/*@{*/
+ /**
+ * Return a map of index:Point<spacedim>, containing the used vertices of the
+ * given `container`. The key of the returned map is the global index in the
+ * triangulation. The used vertices are obtained by looping over all cells,
+ * and querying for each cell where its vertices are through the (optional)
+ * `mapping` argument.
+ *
+ * The size of the returned map equals Triangulation::n_used_vertices(),
+ * (not Triangulation::n_vertices()). If you use the default `mapping`, the
+ * returned map satisfies the following equality:
+ *
+ * @code
+ * used_vertices = extract_used_vertices(tria);
+ * all_vertices = tria.get_vertices();
+ *
+ * for(auto &id_and_v : used_vertices)
+ * all_vertices[id_and_v.first] == id_and_v.second; // true
+ * @endcode
+ *
+ * Notice that the above is not satisfied for mappings that change the
+ * location of vertices, like MappingQEulerian.
+ *
+ * @ref ConceptMeshType "MeshType concept".
+ * @param container The container to extract vertices from
+ * @param mapping The mapping to use to compute the points locations
+ *
+ * @author Luca Heltai, 2017.
+ */
+ template <int dim, int spacedim>
+ std::map<unsigned int,Point<spacedim>> extract_used_vertices (
+ const Triangulation<dim,spacedim> &container,
+ const Mapping<dim,spacedim> &mapping = StaticMappingQ1<dim,spacedim>::mapping);
+
+
/**
* Find and return the number of the used vertex (or marked vertex) in a
* given mesh that is located closest to a given point.
+ template<int dim, int spacedim>
+ std::map<unsigned int, Point<spacedim> >
+ extract_used_vertices(const Triangulation<dim, spacedim> &container,
+ const Mapping<dim,spacedim> &mapping)
+ {
+ std::map<unsigned int, Point<spacedim> > result;
+ for (const auto &cell : container.active_cell_iterators())
+ {
+ auto vs = mapping.get_vertices(cell);
+ for (unsigned int i=0; i<vs.size(); ++i)
+ result[cell->vertex_index(i)]=vs[i];
+ }
+ Assert(result.size() == container.n_used_vertices(),
+ ExcInternalError());
+ return result;
+ }
+
} /* namespace GridTools */
bool
have_same_coarse_mesh (const X &mesh_1,
const X &mesh_2);
-
\}
#endif
std::map<unsigned int,types::global_vertex_index>
compute_local_to_global_vertex_index_map(const parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &triangulation);
+ template
+ std::map<unsigned int, Point<deal_II_space_dimension> >
+ extract_used_vertices(const Triangulation<deal_II_dimension,deal_II_space_dimension>&mesh,
+ const Mapping<deal_II_dimension,deal_II_space_dimension> &mapping);
- \}
+ \}
#endif
}
--- /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.
+
+#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 ()
+{
+ 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;
+};
+
+
+int main ()
+{
+ initlog();
+ test<1,1> ();
+ test<1,2> ();
+ test<1,3> ();
+ test<2,2> ();
+ test<2,3> ();
+ test<3,3> ();
+ return 0;
+}
+
--- /dev/null
+
+DEAL::dim: 1, spacedim: 1
+DEAL::Vertex: 0: 0.00000
+DEAL::Vertex: 1: 1.00000
+DEAL::Vertex: 2: 0.500000
+DEAL::dim: 1, spacedim: 2
+DEAL::Vertex: 0: 0.00000 0.00000
+DEAL::Vertex: 1: 1.00000 0.00000
+DEAL::Vertex: 2: 0.500000 0.00000
+DEAL::dim: 1, spacedim: 3
+DEAL::Vertex: 0: 0.00000 0.00000 0.00000
+DEAL::Vertex: 1: 1.00000 0.00000 0.00000
+DEAL::Vertex: 2: 0.500000 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::dim: 2, spacedim: 3
+DEAL::Vertex: 0: 0.00000 0.00000 0.00000
+DEAL::Vertex: 1: 1.00000 0.00000 0.00000
+DEAL::Vertex: 2: 0.00000 1.00000 0.00000
+DEAL::Vertex: 3: 1.00000 1.00000 0.00000
+DEAL::Vertex: 4: 0.500000 0.00000 0.00000
+DEAL::Vertex: 5: 0.00000 0.500000 0.00000
+DEAL::Vertex: 6: 1.00000 0.500000 0.00000
+DEAL::Vertex: 7: 0.500000 1.00000 0.00000
+DEAL::Vertex: 8: 0.500000 0.500000 0.00000
+DEAL::dim: 3, spacedim: 3
+DEAL::Vertex: 0: 0.00000 0.00000 0.00000
+DEAL::Vertex: 1: 1.00000 0.00000 0.00000
+DEAL::Vertex: 2: 0.00000 1.00000 0.00000
+DEAL::Vertex: 3: 1.00000 1.00000 0.00000
+DEAL::Vertex: 4: 0.00000 0.00000 1.00000
+DEAL::Vertex: 5: 1.00000 0.00000 1.00000
+DEAL::Vertex: 6: 0.00000 1.00000 1.00000
+DEAL::Vertex: 7: 1.00000 1.00000 1.00000
+DEAL::Vertex: 8: 0.500000 0.00000 0.00000
+DEAL::Vertex: 9: 0.00000 0.500000 0.00000
+DEAL::Vertex: 10: 0.00000 0.00000 0.500000
+DEAL::Vertex: 11: 1.00000 0.500000 0.00000
+DEAL::Vertex: 12: 1.00000 0.00000 0.500000
+DEAL::Vertex: 13: 0.500000 1.00000 0.00000
+DEAL::Vertex: 14: 0.00000 1.00000 0.500000
+DEAL::Vertex: 15: 1.00000 1.00000 0.500000
+DEAL::Vertex: 16: 0.500000 0.00000 1.00000
+DEAL::Vertex: 17: 0.00000 0.500000 1.00000
+DEAL::Vertex: 18: 1.00000 0.500000 1.00000
+DEAL::Vertex: 19: 0.500000 1.00000 1.00000
+DEAL::Vertex: 20: 0.500000 0.00000 0.500000
+DEAL::Vertex: 21: 0.500000 0.500000 0.00000
+DEAL::Vertex: 22: 0.00000 0.500000 0.500000
+DEAL::Vertex: 23: 1.00000 0.500000 0.500000
+DEAL::Vertex: 24: 0.500000 1.00000 0.500000
+DEAL::Vertex: 25: 0.500000 0.500000 1.00000
+DEAL::Vertex: 26: 0.500000 0.500000 0.500000