*/
namespace GridTools
{
+ template <int dim, int spacedim>
+ class Cache;
+
/**
* @name Information about meshes and cells
*/
const hp::DoFHandler<dim,spacedim> &mesh,
const Point<spacedim> &p);
+ /**
+ * A version of the previous function that exploits an already existing
+ * GridTools::Cache<dim,spacedim> object.
+ *
+ * @author Luca Heltai, 2017
+ */
+ template <int dim, int spacedim>
+ std::pair<typename Triangulation<dim, spacedim>::active_cell_iterator, Point<dim> >
+ find_active_cell_around_point (const Cache<dim,spacedim> &cache,
+ const Point<spacedim> &p,
+ const typename Triangulation<dim, spacedim>::active_cell_iterator &cell_hint=typename Triangulation<dim, spacedim>::active_cell_iterator(),
+ const std::vector<bool> &marked_vertices = std::vector<bool>());
+
/**
* Return a list of all descendants of the given cell that are active. For
* example, if the current cell is once refined but none of its children are
#include <deal.II/grid/tria_boundary.h>
#include <deal.II/grid/grid_generator.h>
#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/grid_tools_cache.h>
#include <deal.II/grid/grid_reordering.h>
#include <deal.II/dofs/dof_handler.h>
#include <deal.II/dofs/dof_accessor.h>
return id_and_v->first;
}
+ template<int dim, int spacedim>
+ std::pair<typename Triangulation<dim,spacedim>::active_cell_iterator, Point<dim> >
+ find_active_cell_around_point(const Cache<dim,spacedim> &cache,
+ const Point<spacedim> &p,
+ const typename Triangulation<dim,spacedim>::active_cell_iterator &cell_hint,
+ const std::vector<bool> &marked_vertices)
+ {
+ const auto &mesh = cache.get_triangulation();
+ const auto &mapping = cache.get_mapping();
+ const auto &vertex_to_cells = cache.get_vertex_to_cell_map();
+ const auto &vertex_to_cell_centers = cache.get_vertex_to_cell_centers_directions();
+
+ return find_active_cell_around_point(mapping, mesh, p,
+ vertex_to_cells,
+ vertex_to_cell_centers,
+ cell_hint,
+ marked_vertices);
+ }
} /* namespace GridTools */
extract_used_vertices(const Triangulation<deal_II_dimension,deal_II_space_dimension>&mesh,
const Mapping<deal_II_dimension,deal_II_space_dimension> &mapping);
- \}
+ template
+ std::pair<typename Triangulation<deal_II_dimension,deal_II_space_dimension>::active_cell_iterator,
+ Point<deal_II_dimension> >
+ find_active_cell_around_point(const Cache<deal_II_dimension,deal_II_space_dimension>&,
+ const Point<deal_II_space_dimension> &,
+ const typename Triangulation<deal_II_dimension,deal_II_space_dimension>::active_cell_iterator &,
+ const std::vector<bool> &);
+
+
+ \}
#endif
}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 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.
+//
+// ---------------------------------------------------------------------
+
+
+/*
+ * Test cached version of find active cell around point
+ */
+#include "../tests.h"
+#include <deal.II/base/utilities.h>
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/grid_tools_cache.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/hp/fe_collection.h>
+#include <deal.II/hp/fe_values.h>
+#include <deal.II/hp/mapping_collection.h>
+#include <deal.II/hp/q_collection.h>
+#include <deal.II/fe/fe_values.h>
+#include <deal.II/numerics/fe_field_function.h>
+#include <deal.II/lac/vector.h>
+#include <deal.II/numerics/vector_tools.h>
+
+
+template <int dim, int spacedim>
+void test(unsigned int n_ref, unsigned int n_points)
+{
+ deallog << "Testing " << dim << ", " << spacedim << std::endl;
+
+ Triangulation<dim,spacedim> tria;
+ GridGenerator::hyper_cube(tria);
+ tria.refine_global(n_ref);
+ DoFHandler<dim,spacedim> dof_handler(tria);
+
+ std::vector<Point<spacedim>> points;
+
+ deallog << "Points in study: " << n_points << std::endl;
+ for (size_t i=0; i<n_points; ++i)
+ {
+ //We need points in the square [0,1]x[0,1]: this is achieved normalizing with RAND_MAX
+ //RAND_MAX for the used algorithm is 2147483647
+ Point<spacedim> p;
+ for (unsigned int d=0; d<spacedim; ++d)
+ p[d] = double(Testing::rand())/RAND_MAX; //Normalizing the value
+ points.push_back(p);
+ }
+
+ auto &mapping = StaticMappingQ1<dim,spacedim>::mapping;
+
+ GridTools::Cache<dim,spacedim> cache(tria, mapping);
+
+ auto cell = tria.begin_active();
+ for (auto &p : points)
+ {
+ auto c_and_p = GridTools::find_active_cell_around_point(cache, p);
+ auto p2 = mapping.transform_unit_to_real_cell(c_and_p.first, c_and_p.second);
+ if (p2.distance(p)>1e-10)
+ deallog << "NOT OK!" << p << ", " << p2
+ << ", " << c_and_p.first << std::endl;
+ cell = c_and_p.first;
+ }
+ deallog << "OK" << std::endl;
+}
+
+int main ()
+{
+ initlog();
+
+ test<1,1> (4,10);
+ test<2,2> (3,20);
+ test<3,3> (2,30);
+}
--- /dev/null
+
+DEAL::Testing 1, 1
+DEAL::Points in study: 10
+DEAL::OK
+DEAL::Testing 2, 2
+DEAL::Points in study: 20
+DEAL::OK
+DEAL::Testing 3, 3
+DEAL::Points in study: 30
+DEAL::OK
projects / dealii.git / commitdiff