]> https://gitweb.dealii.org/ - dealii.git/commitdiff
update -> mark_for_update. Do work only on demand.
authorLuca Heltai <luca.heltai@sissa.it>
Sat, 16 Sep 2017 06:12:05 +0000 (08:12 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Sat, 16 Sep 2017 06:17:15 +0000 (08:17 +0200)
include/deal.II/grid/grid_tools_cache.h
include/deal.II/grid/grid_tools_cache_flags.h
source/grid/grid_tools_cache.cc
tests/grid/grid_tools_cache_01.cc
tests/grid/grid_tools_cache_02.cc

index 398486e845e263934deb502cabe6c1ce615b9ec4..ee48e160db6eb20fdea31250dc26d16a415c680d 100644 (file)
@@ -44,19 +44,19 @@ namespace GridTools
    * 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.
+   * time to keep track about refinement changes, and allows the user to query
+   * some of the data structures constructed using functions in the GridTools
+   * namespace which are computed only once, and then cached inside this class
+   * for faster access whenever the triangulation has not changed.
    *
    * Notice that this class only notices if the underlying Triangulation has
-   * changed due to a Triangulation::Signals::any_change() signal being triggered.
+   * 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.
+   * 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.
    */
@@ -67,36 +67,28 @@ namespace GridTools
     /**
      * 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.
+     * 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.
+     * Make sure that the objects marked for update are recomputed during
+     * subsequent calls to the `get_*` functions defined in this class.
      *
-     * 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.
+     * Notice that no work is performed when you call this function. The actual
+     * data structures are computed the next time you call the corresponding
+     * `get_*` method.
      *
-     * @param topology_is_unchanged
+     * @param flags What to mark for update
      */
-    void update(bool topology_is_unchanged=false);
+    void mark_for_update(const CacheFlags &flags=update_all);
 
 
     /**
@@ -120,28 +112,18 @@ namespace GridTools
     const KDTree<spacedim> &get_vertex_kdtree() const;
 #endif
 
+  private:
     /**
-     * Exception for uninitialized cache fields.
+     * Keep track of what needs to be updated next.
      */
-    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:
+    mutable CacheFlags update_flags;
+
     /**
      * 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.
      */
@@ -153,19 +135,20 @@ namespace GridTools
      * 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;
+    mutable 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;
+    mutable 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;
+    mutable KDTree<spacedim> vertex_kdtree;
 #endif
 
     /**
index 5071efc04f234f16ee198bce6a324cb5ea31b35a..d3f6880ca62ed7baa73165180747e85c561585b0 100644 (file)
@@ -26,54 +26,59 @@ namespace GridTools
 
   /**
    * The enum type given to the Cache class to select what
-   * information to cache.
+   * information to update.
    *
    * You can select more than one flag by concatenation using the bitwise or
-   * <code>operator|(TriangulationInfoCacheFlags,TriangulationInfoCacheFlags)</code>.
+   * <code>operator|(CacheFlags,CacheFlags)</code>.
    *
    * @author Luca Heltai, 2017.
    */
-  enum TriangulationInfoCacheFlags
+  enum CacheFlags
   {
     /**
-     * Cache Nothing.
+     * Update Nothing.
      */
-    cache_nothing = 0,
+    update_nothing = 0x00,
 
     /**
-     * Cache vertex_to_cell_map, as returned by GridTools::vertex_to_cell_map().
+     * Update vertex_to_cell_map, as returned by GridTools::vertex_to_cell_map().
      */
-    cache_vertex_to_cell_map = 0x0001,
+    update_vertex_to_cell_map = 0x01,
 
     /**
-     * Cache vertex_to_cell_centers_directions, as returned by
+     * Update 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,
+    update_vertex_to_cell_centers_directions = update_vertex_to_cell_map | 0x0002,
 
 #ifdef DEAL_II_WITH_NANOFLANN
     /**
-     * Cache a KDTree object, initialized with the vertices of the Triangulation.
+     * Update a KDTree object, initialized with the vertices of the Triangulation.
      */
-    cache_vertex_kdtree = 0x0004,
+    update_vertex_kdtree = 0x04,
 #endif
+
+    /**
+     * Update all objects.
+     */
+    update_all = 0xFF,
   };
 
 
   /**
    * Output operator which outputs assemble flags as a set of or'd text values.
    *
-   * @ref TriangulationInfoCacheFlags
+   * @ref CacheFlags
    */
   template <class StreamType>
   inline
-  StreamType &operator << (StreamType &s, TriangulationInfoCacheFlags u)
+  StreamType &operator << (StreamType &s, CacheFlags 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";
+    s << " CacheFlags";
+    if (u & update_vertex_to_cell_map)                 s << "|vertex_to_cell_map";
+    if (u & update_vertex_to_cell_centers_directions)  s << "|vertex_to_cells_centers_directions";
 #ifdef DEAL_II_WITH_NANOFLANN
-    if (u & cache_vertex_kdtree)                      s << "|vertex_kdtree";
+    if (u & update_vertex_kdtree)                      s << "|vertex_kdtree";
 #endif
     return s;
   }
@@ -84,19 +89,36 @@ namespace GridTools
    * 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.
+   * assign it to an object of type CacheFlags.
    *
-   * @ref TriangulationInfoCacheFlags
+   * @ref CacheFlags
    */
   inline
-  TriangulationInfoCacheFlags
-  operator | (TriangulationInfoCacheFlags f1, TriangulationInfoCacheFlags f2)
+  CacheFlags
+  operator | (CacheFlags f1, CacheFlags f2)
   {
-    return static_cast<TriangulationInfoCacheFlags> (
+    return static_cast<CacheFlags> (
              static_cast<unsigned int> (f1) |
              static_cast<unsigned int> (f2));
   }
 
+  /**
+   * Global operator which returns an object in which all bits are set which are
+   * not set in the argument. This operator exists since
+   * if it did not then the result of the bit-negation <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 CacheFlags.
+   *
+   * @ref CacheFlags
+   */
+  inline
+  CacheFlags
+  operator ~ (CacheFlags f1)
+  {
+    return static_cast<CacheFlags> (
+             ~static_cast<unsigned int> (f1)
+           );
+  }
 
 
 
@@ -104,11 +126,11 @@ namespace GridTools
    * Global operator which sets the bits from the second argument also in the
    * first one.
    *
-   * @ref TriangulationInfoCacheFlags
+   * @ref CacheFlags
    */
   inline
-  TriangulationInfoCacheFlags &
-  operator |= (TriangulationInfoCacheFlags &f1, TriangulationInfoCacheFlags f2)
+  CacheFlags &
+  operator |= (CacheFlags &f1, CacheFlags f2)
   {
     f1 = f1 | f2;
     return f1;
@@ -120,15 +142,15 @@ namespace GridTools
    * 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.
+   * assign it to an object of type CacheFlags.
    *
-   * @ref TriangulationInfoCacheFlags
+   * @ref CacheFlags
    */
   inline
-  TriangulationInfoCacheFlags
-  operator & (TriangulationInfoCacheFlags f1, TriangulationInfoCacheFlags f2)
+  CacheFlags
+  operator & (CacheFlags f1, CacheFlags f2)
   {
-    return static_cast<TriangulationInfoCacheFlags> (
+    return static_cast<CacheFlags> (
              static_cast<unsigned int> (f1) &
              static_cast<unsigned int> (f2));
   }
@@ -138,11 +160,11 @@ namespace GridTools
    * Global operator which clears all the bits in the first argument if they are
    * not also set in the second argument.
    *
-   * @ref TriangulationInfoCacheFlags
+   * @ref CacheFlags
    */
   inline
-  TriangulationInfoCacheFlags &
-  operator &= (TriangulationInfoCacheFlags &f1, TriangulationInfoCacheFlags f2)
+  CacheFlags &
+  operator &= (CacheFlags &f1, CacheFlags f2)
   {
     f1 = f1 & f2;
     return f1;
index 673dae6fd06917700a521eee24b35105052c2e20..e1cba45da28fa1b1eb77228ea68cffa092ad6974 100644 (file)
@@ -24,19 +24,15 @@ namespace GridTools
   template<int dim, int spacedim>
   Cache<dim,spacedim>::Cache(
     const Triangulation<dim, spacedim> &tria,
-    const TriangulationInfoCacheFlags &flags,
     const Mapping<dim, spacedim> &mapping) :
+    update_flags(update_all),
     tria(&tria),
-    flags(flags),
     mapping(&mapping)
   {
     tria_signal = tria.signals.any_change.connect([&]()
     {
-      update();
+      mark_for_update(update_all);
     });
-
-    if (tria.n_active_cells()>0)
-      update();
   }
 
   template<int dim, int spacedim>
@@ -51,34 +47,9 @@ namespace GridTools
 
 
   template<int dim, int spacedim>
-  void Cache<dim,spacedim>::update(bool topology_is_unchanged)
+  void Cache<dim,spacedim>::mark_for_update(const CacheFlags &flags)
   {
-    // 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);
-
+    update_flags |= flags;
   }
 
 
@@ -87,8 +58,11 @@ namespace GridTools
   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));
+    if (update_flags & update_vertex_to_cell_map)
+      {
+        vertex_to_cells = GridTools::vertex_to_cell_map(*tria);
+        update_flags = update_flags & ~update_vertex_to_cell_map;
+      }
     return vertex_to_cells;
   }
 
@@ -98,8 +72,12 @@ namespace GridTools
   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));
+    if (update_flags & update_vertex_to_cell_centers_directions)
+      {
+        vertex_to_cell_centers =
+          GridTools::vertex_to_cell_centers_directions(*tria, get_vertex_to_cell_map());
+        update_flags = update_flags & ~update_vertex_to_cell_centers_directions;
+      }
     return vertex_to_cell_centers;
   }
 
@@ -109,8 +87,11 @@ namespace GridTools
   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));
+    if (update_flags & update_vertex_kdtree)
+      {
+        vertex_kdtree.set_points(tria->get_vertices());
+        update_flags  = update_flags & ~update_vertex_kdtree;
+      }
     return vertex_kdtree;
   }
 #endif
index a71313f1efe481cae839189b98bfff8e9d8ed1d2..cca324a610055196616f492613a087d9839bfee4 100644 (file)
@@ -31,7 +31,7 @@ void test ()
   GridGenerator::hyper_cube(tria);
   tria.refine_global(1);
 
-  GridTools::Cache<dim> cache(tria, GridTools::cache_vertex_to_cell_map);
+  GridTools::Cache<dim> cache(tria);
 
   auto m = cache.get_vertex_to_cell_map();
 
index 3dc0a063c08d2e512aaa013c2acf66c04990106e..88283c40e899eed96e4d6badade5600468c7cea3 100644 (file)
@@ -28,7 +28,7 @@ void test ()
           << ", spacedim = " << spacedim << std::endl;
 
   Triangulation<dim> tria;
-  GridTools::Cache<dim> cache(tria, GridTools::cache_vertex_to_cell_map);
+  GridTools::Cache<dim> cache(tria);
 
   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.