]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add class MappingQCache
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Thu, 9 May 2019 16:32:19 +0000 (18:32 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Thu, 13 Jun 2019 09:43:58 +0000 (11:43 +0200)
include/deal.II/fe/mapping_q_cache.h [new file with mode: 0644]
include/deal.II/fe/mapping_q_generic.h
source/fe/CMakeLists.txt
source/fe/mapping_q_cache.cc [new file with mode: 0644]
source/fe/mapping_q_cache.inst.in [new file with mode: 0644]

diff --git a/include/deal.II/fe/mapping_q_cache.h b/include/deal.II/fe/mapping_q_cache.h
new file mode 100644 (file)
index 0000000..9f67981
--- /dev/null
@@ -0,0 +1,120 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2019 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_mapping_q_cache_h
+#define dealii_mapping_q_cache_h
+
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/fe/mapping_q_generic.h>
+
+#include <deal.II/grid/tria.h>
+
+
+DEAL_II_NAMESPACE_OPEN
+
+
+
+/*!@addtogroup mapping */
+/*@{*/
+
+
+/**
+ * This class implements a caching strategy for objects of the MappingQ family
+ * in terms of the MappingQGeneric::compute_mapping_support_points() function,
+ * which is used in all operations of MappingQGeneric. The information of the
+ * mapping is pre-computed by the MappingQCache::initialize() function.
+ *
+ * @author Martin Kronbichler, 2019
+ */
+template <int dim, int spacedim = dim>
+class MappingQCache : public MappingQGeneric<dim, spacedim>
+{
+public:
+  /**
+   * Constructor. @p polynomial_degree denotes the polynomial degree of the
+   * polynomials that are used to map cells from the reference to the real
+   * cell.
+   */
+  explicit MappingQCache(const unsigned int polynomial_degree);
+
+  /**
+   * Copy constructor.
+   */
+  explicit MappingQCache(const MappingQCache<dim, spacedim> &mapping);
+
+  /**
+   * Destructor.
+   */
+  ~MappingQCache();
+
+  /**
+   * clone() functionality. For documentation, see Mapping::clone().
+   */
+  virtual std::unique_ptr<Mapping<dim, spacedim>>
+  clone() const override;
+
+  /**
+   * Returns @p false because the preservation of vertex locations depends on
+   * the mapping handed to the reinit() function.
+   */
+  virtual bool
+  preserves_vertex_locations() const override;
+
+  /**
+   * Initialize the data cache by computing the mapping support points for all
+   * cells (on all levels) of the given triangulation. Note that the cache is
+   * invalidated upon the signal Triangulation::Signals::any_change of the
+   * underlying triangulation.
+   */
+  void
+  initialize(const Triangulation<dim, spacedim> &  triangulation,
+             const MappingQGeneric<dim, spacedim> &mapping);
+
+  /**
+   * Return the memory consumption (in bytes) of the cache.
+   */
+  std::size_t
+  memory_consumption() const;
+
+protected:
+  /**
+   * This is the main function overriden from the base class MappingQGeneric.
+   */
+  virtual std::vector<Point<spacedim>>
+  compute_mapping_support_points(
+    const typename Triangulation<dim, spacedim>::cell_iterator &cell)
+    const override;
+
+private:
+  /**
+   * The point cache filled upon calling initialize(). It is made a shared
+   * pointer to allow several instances (created via clone()) to share this
+   * cache.
+   */
+  std::shared_ptr<std::vector<std::vector<std::vector<Point<spacedim>>>>>
+    support_point_cache;
+
+  /**
+   * The connection to Triangulation::signals::any that must be reset once
+   * this class goes out of scope.
+   */
+  boost::signals2::connection clear_signal;
+};
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index aa2e658eef1e76ad22f7c02aec3b9645881a0bd1..5b1f2f8409597102eec25aab1419edebea50e0c8 100644 (file)
@@ -39,6 +39,9 @@ DEAL_II_NAMESPACE_OPEN
 template <int, int>
 class MappingQ;
 
+template <int, int>
+class MappingQCache;
+
 
 /*!@addtogroup mapping */
 /*@{*/
@@ -735,6 +738,13 @@ protected:
    */
   template <int, int>
   friend class MappingQ;
+
+  /**
+   * Make MappingQCache a friend since it needs to call the
+   * compute_mapping_support_points() function.
+   */
+  template <int, int>
+  friend class MappingQCache;
 };
 
 
index 120df5a823c5ce1f46491b5aa1e68f6396c482a8..2beb0674210fe915b824e70495cb22236930107c 100644 (file)
@@ -58,6 +58,7 @@ SET(_unity_include_src
   mapping.cc
   mapping_q1.cc
   mapping_q.cc
+  mapping_q_cache.cc
   mapping_manifold.cc
   )
 
@@ -133,6 +134,7 @@ SET(_inst
   mapping_q_generic.inst.in
   mapping_q1_eulerian.inst.in
   mapping_q1.inst.in
+  mapping_q_cache.inst.in
   mapping_q_eulerian.inst.in
   mapping_q.inst.in
   mapping_manifold.inst.in
diff --git a/source/fe/mapping_q_cache.cc b/source/fe/mapping_q_cache.cc
new file mode 100644 (file)
index 0000000..55925cb
--- /dev/null
@@ -0,0 +1,135 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2019 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/base/memory_consumption.h>
+#include <deal.II/base/work_stream.h>
+
+#include <deal.II/fe/mapping_q_cache.h>
+
+#include <functional>
+
+DEAL_II_NAMESPACE_OPEN
+
+template <int dim, int spacedim>
+MappingQCache<dim, spacedim>::MappingQCache(
+  const unsigned int polynomial_degree)
+  : MappingQGeneric<dim, spacedim>(polynomial_degree)
+{}
+
+
+
+template <int dim, int spacedim>
+MappingQCache<dim, spacedim>::MappingQCache(
+  const MappingQCache<dim, spacedim> &mapping)
+  : MappingQGeneric<dim, spacedim>(mapping)
+  , support_point_cache(mapping.support_point_cache)
+{}
+
+
+
+template <int dim, int spacedim>
+MappingQCache<dim, spacedim>::~MappingQCache()
+{
+  if (clear_signal.connected())
+    clear_signal.disconnect();
+}
+
+
+
+template <int dim, int spacedim>
+std::unique_ptr<Mapping<dim, spacedim>>
+MappingQCache<dim, spacedim>::clone() const
+{
+  return std_cxx14::make_unique<MappingQCache<dim, spacedim>>(*this);
+}
+
+
+
+template <int dim, int spacedim>
+bool
+MappingQCache<dim, spacedim>::preserves_vertex_locations() const
+{
+  return false;
+}
+
+
+
+template <int dim, int spacedim>
+void
+MappingQCache<dim, spacedim>::initialize(
+  const Triangulation<dim, spacedim> &  triangulation,
+  const MappingQGeneric<dim, spacedim> &mapping)
+{
+  AssertDimension(this->get_degree(), mapping.get_degree());
+
+  clear_signal = triangulation.signals.any_change.connect(
+    [&]() -> void { this->support_point_cache.reset(); });
+
+  support_point_cache =
+    std::make_shared<std::vector<std::vector<std::vector<Point<spacedim>>>>>(
+      triangulation.n_levels());
+  for (unsigned int l = 0; l < triangulation.n_levels(); ++l)
+    (*support_point_cache)[l].resize(triangulation.n_raw_cells(l));
+
+  WorkStream::run(
+    triangulation.begin(),
+    triangulation.end(),
+    [&](const typename Triangulation<dim, spacedim>::cell_iterator &cell,
+        void *,
+        void *) {
+      (*support_point_cache)[cell->level()][cell->index()] =
+        mapping.compute_mapping_support_points(cell);
+    },
+    /* copier */ std::function<void(void *)>(),
+    /* scratch_data */ nullptr,
+    /* copy_data */ nullptr,
+    2 * MultithreadInfo::n_threads(),
+    /* chunk_size = */ 1);
+}
+
+
+
+template <int dim, int spacedim>
+std::size_t
+MappingQCache<dim, spacedim>::memory_consumption() const
+{
+  if (support_point_cache.get() != nullptr)
+    return sizeof(*this) +
+           MemoryConsumption::memory_consumption(*support_point_cache);
+  else
+    return sizeof(*this);
+}
+
+
+
+template <int dim, int spacedim>
+std::vector<Point<spacedim>>
+MappingQCache<dim, spacedim>::compute_mapping_support_points(
+  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
+{
+  Assert(support_point_cache.get() != nullptr, ExcNotInitialized());
+
+  AssertIndexRange(cell->level(), support_point_cache->size());
+  AssertIndexRange(cell->index(), (*support_point_cache)[cell->level()].size());
+  return (*support_point_cache)[cell->level()][cell->index()];
+}
+
+
+
+//--------------------------- Explicit instantiations -----------------------
+#include "mapping_q_cache.inst"
+
+
+DEAL_II_NAMESPACE_CLOSE
diff --git a/source/fe/mapping_q_cache.inst.in b/source/fe/mapping_q_cache.inst.in
new file mode 100644 (file)
index 0000000..d57d9e8
--- /dev/null
@@ -0,0 +1,22 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2019 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
+  {
+#if deal_II_dimension <= deal_II_space_dimension
+    template class MappingQCache<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.