]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Created base structure for TriaInfoCache.
authorLuca Heltai <luca.heltai@sissa.it>
Wed, 13 Sep 2017 16:10:19 +0000 (18:10 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Wed, 13 Sep 2017 16:13:17 +0000 (18:13 +0200)
include/deal.II/grid/tria_info_cache.h [new file with mode: 0644]
include/deal.II/grid/tria_info_cache_flags.h [new file with mode: 0644]
source/grid/CMakeLists.txt
source/grid/tria_info_cache.cc [new file with mode: 0644]
source/grid/tria_info_cache.inst.in [new file with mode: 0644]

diff --git a/include/deal.II/grid/tria_info_cache.h b/include/deal.II/grid/tria_info_cache.h
new file mode 100644 (file)
index 0000000..9a9febf
--- /dev/null
@@ -0,0 +1,137 @@
+// ---------------------------------------------------------------------
+//
+// 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 <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 beeing 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>
+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 to work on
+   * @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);
+
+  /**
+   * 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;
+
+  /**
+   * 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<Triangulation<dim,spacedim>,
+               TriangulationInfoCache<dim,spacedim>> tria;
+
+  /**
+   * Specifies what to cache.
+   */
+  const TriangulationInfoCacheFlags flags;
+
+  /**
+   * Mapping to use when computing on the Triangulation.
+   */
+  SmartPointer<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;
+};
+
+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
new file mode 100644 (file)
index 0000000..ddc0234
--- /dev/null
@@ -0,0 +1,133 @@
+// ---------------------------------------------------------------------
+//
+// 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,
+};
+
+
+/**
+   * 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"        ;
+  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 aef9879d8da45cf04bbb26452b1a574b7a7ed2f6..2e26d5f41cc5412a37effae7ee276bd22a4fb05d 100644 (file)
@@ -28,6 +28,7 @@ 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
@@ -61,6 +62,7 @@ SET(_inst
   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/tria_info_cache.cc b/source/grid/tria_info_cache.cc
new file mode 100644 (file)
index 0000000..3f4c759
--- /dev/null
@@ -0,0 +1,56 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 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->signals.any_change.connect([&]()
+  {
+    update();
+  });
+}
+
+template<int dim, int spacedim>
+void TriangulationInfoCache<dim,spacedim>::update(bool topology_is_unchanged)
+{
+  if (cache_vertex_to_cell_map & flags)
+    {
+      vertex_to_cells = GridTools::vertex_to_cell_map(*tria);
+    }
+}
+
+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;
+}
+
+
+
+
+DEAL_II_NAMESPACE_CLOSE
diff --git a/source/grid/tria_info_cache.inst.in b/source/grid/tria_info_cache.inst.in
new file mode 100644 (file)
index 0000000..b55e5ef
--- /dev/null
@@ -0,0 +1,23 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2010 - 2016 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+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>;
+#endif
+}

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.