]> https://gitweb.dealii.org/ - dealii.git/commitdiff
TriangulationInfoCache -> GridTools::Cache
authorLuca Heltai <luca.heltai@sissa.it>
Sat, 16 Sep 2017 05:05:59 +0000 (07:05 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Sat, 16 Sep 2017 06:17:15 +0000 (08:17 +0200)
13 files changed:
doc/news/changes/major/20170914LucaHeltai
include/deal.II/grid/grid_tools_cache.h [new file with mode: 0644]
include/deal.II/grid/grid_tools_cache_flags.h [new file with mode: 0644]
include/deal.II/grid/tria_info_cache.h [deleted file]
include/deal.II/grid/tria_info_cache_flags.h [deleted file]
source/grid/CMakeLists.txt
source/grid/grid_tools_cache.cc [new file with mode: 0644]
source/grid/grid_tools_cache.inst.in [moved from source/grid/tria_info_cache.inst.in with 90% similarity]
source/grid/tria_info_cache.cc [deleted file]
tests/grid/grid_tools_cache_01.cc [moved from tests/grid/tria_info_cache_01.cc with 89% similarity]
tests/grid/grid_tools_cache_01.output [moved from tests/grid/tria_info_cache_01.output with 100% similarity]
tests/grid/grid_tools_cache_02.cc [moved from tests/grid/tria_info_cache_02.cc with 88% similarity]
tests/grid/grid_tools_cache_02.output [moved from tests/grid/tria_info_cache_02.output with 100% similarity]

index 21f9823320790c959c9606a1362333a3c7fb58bd..b9acc89891a0b3668cc2344c84eac8a7d773a20f 100644 (file)
@@ -1,4 +1,4 @@
-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.
diff --git a/include/deal.II/grid/grid_tools_cache.h b/include/deal.II/grid/grid_tools_cache.h
new file mode 100644 (file)
index 0000000..398486e
--- /dev/null
@@ -0,0 +1,181 @@
+// ---------------------------------------------------------------------
+//
+// 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
diff --git a/include/deal.II/grid/grid_tools_cache_flags.h b/include/deal.II/grid/grid_tools_cache_flags.h
new file mode 100644 (file)
index 0000000..5071efc
--- /dev/null
@@ -0,0 +1,154 @@
+// ---------------------------------------------------------------------
+//
+// 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
diff --git a/include/deal.II/grid/tria_info_cache.h b/include/deal.II/grid/tria_info_cache.h
deleted file mode 100644 (file)
index 26ddfcb..0000000
+++ /dev/null
@@ -1,176 +0,0 @@
-// ---------------------------------------------------------------------
-//
-// 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
diff --git a/include/deal.II/grid/tria_info_cache_flags.h b/include/deal.II/grid/tria_info_cache_flags.h
deleted file mode 100644 (file)
index 0c8a420..0000000
+++ /dev/null
@@ -1,150 +0,0 @@
-// ---------------------------------------------------------------------
-//
-// 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
index 2e26d5f41cc5412a37effae7ee276bd22a4fb05d..d68c7735d4202241b51766bc3a6bad2b253c9dd9 100644 (file)
@@ -21,6 +21,7 @@ SET(_unity_include_src
   grid_in.cc
   grid_out.cc
   grid_refinement.cc
+  grid_tools_cache.cc
   intergrid_map.cc
   manifold.cc
   manifold_lib.cc
@@ -28,7 +29,6 @@ SET(_unity_include_src
   tria_accessor.cc
   tria_boundary.cc
   tria_boundary_lib.cc
-  tria_info_cache.cc
   tria_faces.cc
   tria_levels.cc
   tria_objects.cc
@@ -56,13 +56,13 @@ SET(_inst
   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
   )
diff --git a/source/grid/grid_tools_cache.cc b/source/grid/grid_tools_cache.cc
new file mode 100644 (file)
index 0000000..673dae6
--- /dev/null
@@ -0,0 +1,123 @@
+// ---------------------------------------------------------------------
+//
+// 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
+
similarity index 90%
rename from source/grid/tria_info_cache.inst.in
rename to source/grid/grid_tools_cache.inst.in
index 7f03011ace9a8fcdb5a2ce45583901358c0588d5..b50399a44efea38a0faac4b01e77612310e5159f 100644 (file)
@@ -18,6 +18,6 @@
 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
 }
diff --git a/source/grid/tria_info_cache.cc b/source/grid/tria_info_cache.cc
deleted file mode 100644 (file)
index e9c8913..0000000
+++ /dev/null
@@ -1,118 +0,0 @@
-// ---------------------------------------------------------------------
-//
-// 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
-
similarity index 89%
rename from tests/grid/tria_info_cache_01.cc
rename to tests/grid/grid_tools_cache_01.cc
index 78d5e2dbf4929c76ea99954dc62e843c69f19efa..a71313f1efe481cae839189b98bfff8e9d8ed1d2 100644 (file)
 //
 // ---------------------------------------------------------------------
 
-// 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>
@@ -31,7 +31,7 @@ void test ()
   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();
 
similarity index 88%
rename from tests/grid/tria_info_cache_02.cc
rename to tests/grid/grid_tools_cache_02.cc
index 4280b6c8e3867f3fde351ae64ddd43bca8e4ef70..3dc0a063c08d2e512aaa013c2acf66c04990106e 100644 (file)
 //
 // ---------------------------------------------------------------------
 
-// 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>
@@ -28,7 +28,7 @@ void test ()
           << ", 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);
 

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.