-New: A new TriangulationInfoCache class
+New: A new GridTools::Cache class
allows caching of some expensive data of the
Triangulation class, computed generally using
functions in GridTools.
--- /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.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_grid_grid_tools_cache_h
+#define dealii_grid_grid_tools_cache_h
+
+
+#include <deal.II/base/config.h>
+#include <deal.II/base/subscriptor.h>
+#include <deal.II/base/exceptions.h>
+#include <deal.II/base/point.h>
+
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_tools_cache_flags.h>
+#include <deal.II/fe/mapping_q1.h>
+
+#include <deal.II/numerics/kdtree.h>
+
+#include <boost/signals2.hpp>
+
+#include <cmath>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace GridTools
+{
+
+ /**
+ * A class that caches computationally intensive information about a
+ * Triangulation.
+ *
+ * This class attaches a signal to the Triangulation passed at construction
+ * time, and rebuilds all the structures specified in the list of specified
+ * TriangulationInfoCacheFlags, whenever the Triangulation is changed.
+ *
+ * If you try to access one of the objects that has not been cached (due to a
+ * missing TriangulationInfoCacheFlags flag), an Assertion is thrown.
+ *
+ * Notice that this class only notices if the underlying Triangulation has
+ * changed due to a Triangulation::Signals::any_change() signal being triggered.
+ *
+ * If the triangulation changes due to a MappingQEulerian being passed to this
+ * class, or because you manually change some vertex locations, then some of
+ * the structures in this class become obsolete, and you will have to call the
+ * method update() manually.
+ *
+ * @author Luca Heltai, 2017.
+ */
+ template <int dim, int spacedim = dim>
+ class Cache : public Subscriptor
+ {
+ public:
+ /**
+ * Constructor.
+ *
+ * You can specify what to cache by passing appropriate
+ * TriangulationInfoCacheFlags. If you provide the optional `mapping`
+ * argument, then this is used whenever a mapping is required.
+ *
+ * @param tria The triangulation for which to store information
+ * @param flags TriangulationInfoCacheFlags that specify what to cache
+ * @param mapping The mapping to use when computing cached objects
+ */
+ Cache (const Triangulation<dim,spacedim> &tria,
+ const TriangulationInfoCacheFlags &flags = cache_nothing,
+ const Mapping<dim,spacedim> &mapping=StaticMappingQ1<dim,spacedim>::mapping);
+
+ ~Cache();
+
+ /**
+ * Loop over all TriangulationInfoCacheFlags and rebuilds the specific data
+ * structure.
+ *
+ * The optional parameter can be set to true to update only data that is
+ * dependent on the position of the vertices. This may happen, for example,
+ * if you have a MappingQEulerian class associated to the triangulation. In
+ * this case, even if the Triangulation has not changed, the vector that
+ * specifies the position of the vertices may actually change *independently*
+ * of the triangulation, requiring you to recompute vertex dependent things.
+ * In this case, the topological information is not changed, and does not
+ * need to be updated.
+ *
+ * @param topology_is_unchanged
+ */
+ void update(bool topology_is_unchanged=false);
+
+
+ /**
+ * Return the cached vertex_to_cell_map as computed by GridTools::vertex_to_cell_map().
+ */
+ const std::vector< std::set<typename Triangulation<dim,spacedim>::active_cell_iterator> >
+ &get_vertex_to_cell_map() const;
+
+ /**
+ * Return the cached vertex_to_cell_centers_directions as computed by
+ * GridTools::vertex_to_cell_centers_directions().
+ */
+ const std::vector<std::vector<Tensor<1,spacedim>>>
+ &get_vertex_to_cell_centers_directions() const;
+
+#ifdef DEAL_II_WITH_NANOFLANN
+ /**
+ * Return the cached vertex_kdtree object, constructed with the vertices of
+ * the stored triangulation.
+ */
+ const KDTree<spacedim> &get_vertex_kdtree() const;
+#endif
+
+ /**
+ * Exception for uninitialized cache fields.
+ */
+ DeclException2(ExcAccessToInvalidCacheObject,
+ TriangulationInfoCacheFlags &,
+ TriangulationInfoCacheFlags &,
+ << "You tried to access a cache object that has not "
+ << "been constructed by this class. You specified "
+ << arg2 << " at construction time, but you need also "
+ << arg1 << " to access the field you are requesting now.");
+ private:
+ /**
+ * A pointer to the Triangulation.
+ */
+ SmartPointer<const Triangulation<dim,spacedim>,
+ Cache<dim,spacedim>> tria;
+
+ /**
+ * Specifies what to cache.
+ */
+ const TriangulationInfoCacheFlags flags;
+
+ /**
+ * Mapping to use when computing on the Triangulation.
+ */
+ SmartPointer<const Mapping<dim,spacedim>,
+ Cache<dim,spacedim>> mapping;
+
+
+ /**
+ * Store vertex to cell map information, as generated by
+ * GridTools::vertex_to_cell_map()
+ */
+ std::vector< std::set<typename Triangulation<dim,spacedim>::active_cell_iterator> > vertex_to_cells;
+
+ /**
+ * Store vertex to cell center directions, as generated by
+ * GridTools::vertex_to_cell_centers_directions.
+ */
+ std::vector<std::vector<Tensor<1,spacedim>>> vertex_to_cell_centers;
+
+#ifdef DEAL_II_WITH_NANOFLANN
+ /**
+ * A KDTree object constructed with the vertices of the triangulation.
+ */
+ KDTree<spacedim> vertex_kdtree;
+#endif
+
+ /**
+ * Storage for the status of the triangulation signal.
+ */
+ boost::signals2::connection tria_signal;
+ };
+
+}
+
+DEAL_II_NAMESPACE_CLOSE
+
+#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.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_grid_tria_info_cache_flags_h
+#define dealii_grid_tria_info_cache_flags_h
+
+
+#include <deal.II/base/config.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace GridTools
+{
+
+ /**
+ * The enum type given to the Cache class to select what
+ * information to cache.
+ *
+ * You can select more than one flag by concatenation using the bitwise or
+ * <code>operator|(TriangulationInfoCacheFlags,TriangulationInfoCacheFlags)</code>.
+ *
+ * @author Luca Heltai, 2017.
+ */
+ enum TriangulationInfoCacheFlags
+ {
+ /**
+ * Cache Nothing.
+ */
+ cache_nothing = 0,
+
+ /**
+ * Cache vertex_to_cell_map, as returned by GridTools::vertex_to_cell_map().
+ */
+ cache_vertex_to_cell_map = 0x0001,
+
+ /**
+ * Cache vertex_to_cell_centers_directions, as returned by
+ * GridTools::vertex_to_cell_centers_directions()
+ */
+ cache_vertex_to_cell_centers_directions = cache_vertex_to_cell_map | 0x0002,
+
+#ifdef DEAL_II_WITH_NANOFLANN
+ /**
+ * Cache a KDTree object, initialized with the vertices of the Triangulation.
+ */
+ cache_vertex_kdtree = 0x0004,
+#endif
+ };
+
+
+ /**
+ * Output operator which outputs assemble flags as a set of or'd text values.
+ *
+ * @ref TriangulationInfoCacheFlags
+ */
+ template <class StreamType>
+ inline
+ StreamType &operator << (StreamType &s, TriangulationInfoCacheFlags u)
+ {
+ s << " TriangulationInfoCacheFlags";
+ if (u & cache_vertex_to_cell_map) s << "|vertex_to_cell_map";
+ if (u & cache_vertex_to_cell_centers_directions) s << "|vertex_to_cells_centers_directions";
+#ifdef DEAL_II_WITH_NANOFLANN
+ if (u & cache_vertex_kdtree) s << "|vertex_kdtree";
+#endif
+ return s;
+ }
+
+
+ /**
+ * Global operator which returns an object in which all bits are set which are
+ * either set in the first or the second argument. This operator exists since
+ * if it did not then the result of the bit-or <tt>operator |</tt> would be an
+ * integer which would in turn trigger a compiler warning when we tried to
+ * assign it to an object of type TriangulationInfoCacheFlags.
+ *
+ * @ref TriangulationInfoCacheFlags
+ */
+ inline
+ TriangulationInfoCacheFlags
+ operator | (TriangulationInfoCacheFlags f1, TriangulationInfoCacheFlags f2)
+ {
+ return static_cast<TriangulationInfoCacheFlags> (
+ static_cast<unsigned int> (f1) |
+ static_cast<unsigned int> (f2));
+ }
+
+
+
+
+ /**
+ * Global operator which sets the bits from the second argument also in the
+ * first one.
+ *
+ * @ref TriangulationInfoCacheFlags
+ */
+ inline
+ TriangulationInfoCacheFlags &
+ operator |= (TriangulationInfoCacheFlags &f1, TriangulationInfoCacheFlags f2)
+ {
+ f1 = f1 | f2;
+ return f1;
+ }
+
+
+ /**
+ * Global operator which returns an object in which all bits are set which are
+ * set in the first as well as the second argument. This operator exists since
+ * if it did not then the result of the bit-and <tt>operator &</tt> would be
+ * an integer which would in turn trigger a compiler warning when we tried to
+ * assign it to an object of type TriangulationInfoCacheFlags.
+ *
+ * @ref TriangulationInfoCacheFlags
+ */
+ inline
+ TriangulationInfoCacheFlags
+ operator & (TriangulationInfoCacheFlags f1, TriangulationInfoCacheFlags f2)
+ {
+ return static_cast<TriangulationInfoCacheFlags> (
+ static_cast<unsigned int> (f1) &
+ static_cast<unsigned int> (f2));
+ }
+
+
+ /**
+ * Global operator which clears all the bits in the first argument if they are
+ * not also set in the second argument.
+ *
+ * @ref TriangulationInfoCacheFlags
+ */
+ inline
+ TriangulationInfoCacheFlags &
+ operator &= (TriangulationInfoCacheFlags &f1, TriangulationInfoCacheFlags f2)
+ {
+ f1 = f1 & f2;
+ return f1;
+ }
+
+}
+DEAL_II_NAMESPACE_CLOSE
+
+#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.
-//
-// ---------------------------------------------------------------------
-
-#ifndef dealii_grid_tria_info_cache_h
-#define dealii_grid_tria_info_cache_h
-
-
-#include <deal.II/base/config.h>
-#include <deal.II/base/subscriptor.h>
-#include <deal.II/base/exceptions.h>
-#include <deal.II/base/point.h>
-
-#include <deal.II/grid/tria.h>
-#include <deal.II/grid/tria_iterator.h>
-#include <deal.II/grid/tria_accessor.h>
-#include <deal.II/grid/tria_info_cache_flags.h>
-#include <deal.II/fe/mapping_q1.h>
-
-#include <deal.II/numerics/kdtree.h>
-
-#include <boost/signals2.hpp>
-
-#include <cmath>
-
-DEAL_II_NAMESPACE_OPEN
-
-/**
- * A class that caches computationally intensive information about a
- * Triangulation.
- *
- * This class attaches a signal to the Triangulation passed at construction
- * time, and rebuilds all the structures specified in the list of specified
- * TriangulationInfoCacheFlags, whenever the Triangulation is changed.
- *
- * If you try to access one of the objects that has not been cached (due to a
- * missing TriangulationInfoCacheFlags flag), an Assertion is thrown.
- *
- * Notice that this class only notices if the underlying Triangulation has
- * changed due to a Triangulation::Signals::any_change() signal being triggered.
- *
- * If the triangulation changes due to a MappingQEulerian being passed to this
- * class, or because you manually change some vertex locations, then some of
- * the structures in this class become obsolete, and you will have to call the
- * method update() manually.
- *
- * @author Luca Heltai, 2017.
- */
-template <int dim, int spacedim = dim>
-class TriangulationInfoCache : public Subscriptor
-{
-public:
- /**
- * Constructor.
- *
- * You can specify what to cache by passing appropriate
- * TriangulationInfoCacheFlags. If you provide the optional `mapping`
- * argument, then this is used whenever a mapping is required.
- *
- * @param tria The triangulation for which to store information
- * @param flags TriangulationInfoCacheFlags that specify what to cache
- * @param mapping The mapping to use when computing cached objects
- */
- TriangulationInfoCache (const Triangulation<dim,spacedim> &tria,
- const TriangulationInfoCacheFlags &flags = cache_nothing,
- const Mapping<dim,spacedim> &mapping=StaticMappingQ1<dim,spacedim>::mapping);
-
- ~TriangulationInfoCache();
-
- /**
- * Loop over all TriangulationInfoCacheFlags and rebuilds the specific data
- * structure.
- *
- * The optional parameter can be set to true to update only data that is
- * dependent on the position of the vertices. This may happen, for example,
- * if you have a MappingQEulerian class associated to the triangulation. In
- * this case, even if the Triangulation has not changed, the vector that
- * specifies the position of the vertices may actually change *independently*
- * of the triangulation, requiring you to recompute vertex dependent things.
- * In this case, the topological information is not changed, and does not
- * need to be updated.
- *
- * @param topology_is_unchanged
- */
- void update(bool topology_is_unchanged=false);
-
-
- /**
- * Return the cached vertex_to_cell_map as computed by GridTools::vertex_to_cell_map().
- */
- const std::vector< std::set<typename Triangulation<dim,spacedim>::active_cell_iterator> >
- &get_vertex_to_cell_map() const;
-
- /**
- * Return the cached vertex_to_cell_centers_directions as computed by
- * GridTools::vertex_to_cell_centers_directions().
- */
- const std::vector<std::vector<Tensor<1,spacedim>>>
- &get_vertex_to_cell_centers_directions() const;
-
-#ifdef DEAL_II_WITH_NANOFLANN
- /**
- * Return the cached vertex_kdtree object, constructed with the vertices of
- * the stored triangulation.
- */
- const KDTree<spacedim> &get_vertex_kdtree() const;
-#endif
-
- /**
- * Exception for uninitialized cache fields.
- */
- DeclException2(ExcAccessToInvalidCacheObject,
- TriangulationInfoCacheFlags &,
- TriangulationInfoCacheFlags &,
- << "You tried to access a cache object that has not "
- << "been constructed by this class. You specified "
- << arg2 << " at construction time, but you need also "
- << arg1 << " to access the field you are requesting now.");
-private:
- /**
- * A pointer to the Triangulation.
- */
- SmartPointer<const Triangulation<dim,spacedim>,
- TriangulationInfoCache<dim,spacedim>> tria;
-
- /**
- * Specifies what to cache.
- */
- const TriangulationInfoCacheFlags flags;
-
- /**
- * Mapping to use when computing on the Triangulation.
- */
- SmartPointer<const Mapping<dim,spacedim>,
- TriangulationInfoCache<dim,spacedim>> mapping;
-
-
- /**
- * Store vertex to cell map information, as generated by
- * GridTools::vertex_to_cell_map()
- */
- std::vector< std::set<typename Triangulation<dim,spacedim>::active_cell_iterator> > vertex_to_cells;
-
- /**
- * Store vertex to cell center directions, as generated by
- * GridTools::vertex_to_cell_centers_directions.
- */
- std::vector<std::vector<Tensor<1,spacedim>>> vertex_to_cell_centers;
-
-#ifdef DEAL_II_WITH_NANOFLANN
- /**
- * A KDTree object constructed with the vertices of the triangulation.
- */
- KDTree<spacedim> vertex_kdtree;
-#endif
-
- /**
- * Storage for the status of the triangulation signal.
- */
- boost::signals2::connection tria_signal;
-};
-
-DEAL_II_NAMESPACE_CLOSE
-
-#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.
-//
-// ---------------------------------------------------------------------
-
-#ifndef dealii_grid_tria_info_cache_flags_h
-#define dealii_grid_tria_info_cache_flags_h
-
-
-#include <deal.II/base/config.h>
-
-DEAL_II_NAMESPACE_OPEN
-
-/**
- * The enum type given to the TriangulationInfoCache class to select what
- * information to cache.
- *
- * You can select more than one flag by concatenation using the bitwise or
- * <code>operator|(TriangulationInfoCacheFlags,TriangulationInfoCacheFlags)</code>.
- *
- * @author Luca Heltai, 2017.
- */
-enum TriangulationInfoCacheFlags
-{
- /**
- * Cache Nothing.
- */
- cache_nothing = 0,
-
- /**
- * Cache vertex_to_cell_map, as returned by GridTools::vertex_to_cell_map().
- */
- cache_vertex_to_cell_map = 0x0001,
-
- /**
- * Cache vertex_to_cell_centers_directions, as returned by
- * GridTools::vertex_to_cell_centers_directions()
- */
- cache_vertex_to_cell_centers_directions = cache_vertex_to_cell_map | 0x0002,
-
-#ifdef DEAL_II_WITH_NANOFLANN
- /**
- * Cache a KDTree object, initialized with the vertices of the Triangulation.
- */
- cache_vertex_kdtree = 0x0004,
-#endif
-};
-
-
-/**
- * Output operator which outputs assemble flags as a set of or'd text values.
- *
- * @ref TriangulationInfoCacheFlags
- */
-template <class StreamType>
-inline
-StreamType &operator << (StreamType &s, TriangulationInfoCacheFlags u)
-{
- s << " TriangulationInfoCacheFlags";
- if (u & cache_vertex_to_cell_map) s << "|vertex_to_cell_map";
- if (u & cache_vertex_to_cell_centers_directions) s << "|vertex_to_cells_centers_directions";
-#ifdef DEAL_II_WITH_NANOFLANN
- if (u & cache_vertex_kdtree) s << "|vertex_kdtree";
-#endif
- return s;
-}
-
-
-/**
- * Global operator which returns an object in which all bits are set which are
- * either set in the first or the second argument. This operator exists since
- * if it did not then the result of the bit-or <tt>operator |</tt> would be an
- * integer which would in turn trigger a compiler warning when we tried to
- * assign it to an object of type TriangulationInfoCacheFlags.
- *
- * @ref TriangulationInfoCacheFlags
- */
-inline
-TriangulationInfoCacheFlags
-operator | (TriangulationInfoCacheFlags f1, TriangulationInfoCacheFlags f2)
-{
- return static_cast<TriangulationInfoCacheFlags> (
- static_cast<unsigned int> (f1) |
- static_cast<unsigned int> (f2));
-}
-
-
-
-
-/**
- * Global operator which sets the bits from the second argument also in the
- * first one.
- *
- * @ref TriangulationInfoCacheFlags
- */
-inline
-TriangulationInfoCacheFlags &
-operator |= (TriangulationInfoCacheFlags &f1, TriangulationInfoCacheFlags f2)
-{
- f1 = f1 | f2;
- return f1;
-}
-
-
-/**
- * Global operator which returns an object in which all bits are set which are
- * set in the first as well as the second argument. This operator exists since
- * if it did not then the result of the bit-and <tt>operator &</tt> would be
- * an integer which would in turn trigger a compiler warning when we tried to
- * assign it to an object of type TriangulationInfoCacheFlags.
- *
- * @ref TriangulationInfoCacheFlags
- */
-inline
-TriangulationInfoCacheFlags
-operator & (TriangulationInfoCacheFlags f1, TriangulationInfoCacheFlags f2)
-{
- return static_cast<TriangulationInfoCacheFlags> (
- static_cast<unsigned int> (f1) &
- static_cast<unsigned int> (f2));
-}
-
-
-/**
- * Global operator which clears all the bits in the first argument if they are
- * not also set in the second argument.
- *
- * @ref TriangulationInfoCacheFlags
- */
-inline
-TriangulationInfoCacheFlags &
-operator &= (TriangulationInfoCacheFlags &f1, TriangulationInfoCacheFlags f2)
-{
- f1 = f1 & f2;
- return f1;
-}
-
-DEAL_II_NAMESPACE_CLOSE
-
-#endif
grid_in.cc
grid_out.cc
grid_refinement.cc
+ grid_tools_cache.cc
intergrid_map.cc
manifold.cc
manifold_lib.cc
tria_accessor.cc
tria_boundary.cc
tria_boundary_lib.cc
- tria_info_cache.cc
tria_faces.cc
tria_levels.cc
tria_objects.cc
grid_out.inst.in
grid_refinement.inst.in
grid_tools.inst.in
+ grid_tools_cache.inst.in
intergrid_map.inst.in
manifold.inst.in
manifold_lib.inst.in
tria_accessor.inst.in
tria_boundary.inst.in
tria_boundary_lib.inst.in
- tria_info_cache.inst.in
tria.inst.in
tria_objects.inst.in
)
--- /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.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/grid/grid_tools_cache.h>
+#include <deal.II/grid/grid_tools.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace GridTools
+{
+
+ template<int dim, int spacedim>
+ Cache<dim,spacedim>::Cache(
+ const Triangulation<dim, spacedim> &tria,
+ const TriangulationInfoCacheFlags &flags,
+ const Mapping<dim, spacedim> &mapping) :
+ tria(&tria),
+ flags(flags),
+ mapping(&mapping)
+ {
+ tria_signal = tria.signals.any_change.connect([&]()
+ {
+ update();
+ });
+
+ if (tria.n_active_cells()>0)
+ update();
+ }
+
+ template<int dim, int spacedim>
+ Cache<dim,spacedim>::~Cache()
+ {
+ // Make sure that the signals that was attached to the triangulation
+ // is removed here.
+ if (tria_signal.connected())
+ tria_signal.disconnect();
+ }
+
+
+
+ template<int dim, int spacedim>
+ void Cache<dim,spacedim>::update(bool topology_is_unchanged)
+ {
+ // If the triangulation is empty, just clear everything.
+ if (tria->n_active_cells() == 0)
+ {
+ vertex_to_cells.clear();
+ vertex_to_cell_centers.clear();
+#ifdef DEAL_II_WITH_NANOFLANN
+ vertex_kdtree.set_points(tria->get_vertices());
+#endif
+ return;
+ }
+
+ if (topology_is_unchanged == false)
+ {
+ if (cache_vertex_to_cell_map & flags)
+ vertex_to_cells = GridTools::vertex_to_cell_map(*tria);
+
+#ifdef DEAL_II_WITH_NANOFLANN
+ if (cache_vertex_kdtree & flags)
+ vertex_kdtree.set_points(tria->get_vertices());
+#endif
+ }
+
+ if (cache_vertex_to_cell_centers_directions & flags)
+ vertex_to_cell_centers =
+ GridTools::vertex_to_cell_centers_directions(*tria, vertex_to_cells);
+
+ }
+
+
+
+ template<int dim, int spacedim>
+ const std::vector<std::set<typename Triangulation<dim,spacedim>::active_cell_iterator> > &
+ Cache<dim,spacedim>::get_vertex_to_cell_map() const
+ {
+ Assert(flags & cache_vertex_to_cell_map,
+ ExcAccessToInvalidCacheObject(cache_vertex_to_cell_map, flags));
+ return vertex_to_cells;
+ }
+
+
+
+ template<int dim, int spacedim>
+ const std::vector<std::vector<Tensor<1,spacedim>>> &
+ Cache<dim,spacedim>::get_vertex_to_cell_centers_directions() const
+ {
+ Assert(flags & cache_vertex_to_cell_centers_directions,
+ ExcAccessToInvalidCacheObject(cache_vertex_to_cell_centers_directions, flags));
+ return vertex_to_cell_centers;
+ }
+
+
+
+#ifdef DEAL_II_WITH_NANOFLANN
+ template<int dim, int spacedim>
+ const KDTree<spacedim> &Cache<dim,spacedim>::get_vertex_kdtree() const
+ {
+ Assert(flags & cache_vertex_kdtree,
+ ExcAccessToInvalidCacheObject(cache_vertex_kdtree, flags));
+ return vertex_kdtree;
+ }
+#endif
+
+#include "grid_tools_cache.inst"
+
+} // GridTools
+
+DEAL_II_NAMESPACE_CLOSE
+
for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
{
#if deal_II_dimension <= deal_II_space_dimension
- template class TriangulationInfoCache<deal_II_dimension, deal_II_space_dimension>;
+ template class Cache<deal_II_dimension, deal_II_space_dimension>;
#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.
-//
-// ---------------------------------------------------------------------
-
-#include <deal.II/grid/tria_info_cache.h>
-#include <deal.II/grid/grid_tools.h>
-
-DEAL_II_NAMESPACE_OPEN
-
-template<int dim, int spacedim>
-TriangulationInfoCache<dim,spacedim>::TriangulationInfoCache(
- const Triangulation<dim, spacedim> &tria,
- const TriangulationInfoCacheFlags &flags,
- const Mapping<dim, spacedim> &mapping) :
- tria(&tria),
- flags(flags),
- mapping(&mapping)
-{
- tria_signal = tria.signals.any_change.connect([&]()
- {
- update();
- });
-
- if (tria.n_active_cells()>0)
- update();
-}
-
-template<int dim, int spacedim>
-TriangulationInfoCache<dim,spacedim>::~TriangulationInfoCache()
-{
- // Make sure that the signals that was attached to the triangulation
- // is removed here.
- if (tria_signal.connected())
- tria_signal.disconnect();
-}
-
-
-
-template<int dim, int spacedim>
-void TriangulationInfoCache<dim,spacedim>::update(bool topology_is_unchanged)
-{
- // If the triangulation is empty, just clear everything.
- if (tria->n_active_cells() == 0)
- {
- vertex_to_cells.clear();
- vertex_to_cell_centers.clear();
-#ifdef DEAL_II_WITH_NANOFLANN
- vertex_kdtree.set_points(tria->get_vertices());
-#endif
- return;
- }
-
- if (topology_is_unchanged == false)
- {
- if (cache_vertex_to_cell_map & flags)
- vertex_to_cells = GridTools::vertex_to_cell_map(*tria);
-
-#ifdef DEAL_II_WITH_NANOFLANN
- if (cache_vertex_kdtree & flags)
- vertex_kdtree.set_points(tria->get_vertices());
-#endif
- }
-
- if (cache_vertex_to_cell_centers_directions & flags)
- vertex_to_cell_centers =
- GridTools::vertex_to_cell_centers_directions(*tria, vertex_to_cells);
-
-}
-
-
-
-template<int dim, int spacedim>
-const std::vector<std::set<typename Triangulation<dim,spacedim>::active_cell_iterator> > &
-TriangulationInfoCache<dim,spacedim>::get_vertex_to_cell_map() const
-{
- Assert(flags & cache_vertex_to_cell_map,
- ExcAccessToInvalidCacheObject(cache_vertex_to_cell_map, flags));
- return vertex_to_cells;
-}
-
-
-
-template<int dim, int spacedim>
-const std::vector<std::vector<Tensor<1,spacedim>>> &
-TriangulationInfoCache<dim,spacedim>::get_vertex_to_cell_centers_directions() const
-{
- Assert(flags & cache_vertex_to_cell_centers_directions,
- ExcAccessToInvalidCacheObject(cache_vertex_to_cell_centers_directions, flags));
- return vertex_to_cell_centers;
-}
-
-
-
-#ifdef DEAL_II_WITH_NANOFLANN
-template<int dim, int spacedim>
-const KDTree<spacedim> &TriangulationInfoCache<dim,spacedim>::get_vertex_kdtree() const
-{
- Assert(flags & cache_vertex_kdtree,
- ExcAccessToInvalidCacheObject(cache_vertex_kdtree, flags));
- return vertex_kdtree;
-}
-#endif
-
-#include "tria_info_cache.inst"
-
-DEAL_II_NAMESPACE_CLOSE
-
//
// ---------------------------------------------------------------------
-// Validate tria_info_cache
+// Validate grid_tools_cache
#include "../tests.h"
#include <deal.II/grid/tria.h>
#include <deal.II/grid/grid_generator.h>
-#include <deal.II/grid/tria_info_cache.h>
+#include <deal.II/grid/grid_tools_cache.h>
template <int dim, int spacedim>
GridGenerator::hyper_cube(tria);
tria.refine_global(1);
- TriangulationInfoCache<dim> cache(tria, cache_vertex_to_cell_map);
+ GridTools::Cache<dim> cache(tria, GridTools::cache_vertex_to_cell_map);
auto m = cache.get_vertex_to_cell_map();
//
// ---------------------------------------------------------------------
-// Validate tria_info_cache. Different construction order. Check signal.
+// Validate grid_tools_cache. Different construction order. Check signal.
#include "../tests.h"
#include <deal.II/grid/tria.h>
#include <deal.II/grid/grid_generator.h>
-#include <deal.II/grid/tria_info_cache.h>
+#include <deal.II/grid/grid_tools_cache.h>
template <int dim, int spacedim>
<< ", spacedim = " << spacedim << std::endl;
Triangulation<dim> tria;
- TriangulationInfoCache<dim> cache(tria, cache_vertex_to_cell_map);
+ GridTools::Cache<dim> cache(tria, GridTools::cache_vertex_to_cell_map);
GridGenerator::hyper_cube(tria);