From 6032f78ccded16696eb281482cf66cfe08c20258 Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Thu, 9 May 2019 18:32:19 +0200 Subject: [PATCH] Add class MappingQCache --- include/deal.II/fe/mapping_q_cache.h | 120 ++++++++++++++++++++++ include/deal.II/fe/mapping_q_generic.h | 10 ++ source/fe/CMakeLists.txt | 2 + source/fe/mapping_q_cache.cc | 135 +++++++++++++++++++++++++ source/fe/mapping_q_cache.inst.in | 22 ++++ 5 files changed, 289 insertions(+) create mode 100644 include/deal.II/fe/mapping_q_cache.h create mode 100644 source/fe/mapping_q_cache.cc create mode 100644 source/fe/mapping_q_cache.inst.in diff --git a/include/deal.II/fe/mapping_q_cache.h b/include/deal.II/fe/mapping_q_cache.h new file mode 100644 index 0000000000..9f67981425 --- /dev/null +++ b/include/deal.II/fe/mapping_q_cache.h @@ -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 + +#include + +#include + + +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 +class MappingQCache : public MappingQGeneric +{ +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 &mapping); + + /** + * Destructor. + */ + ~MappingQCache(); + + /** + * clone() functionality. For documentation, see Mapping::clone(). + */ + virtual std::unique_ptr> + 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 & triangulation, + const MappingQGeneric &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> + compute_mapping_support_points( + const typename Triangulation::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>>>> + 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 diff --git a/include/deal.II/fe/mapping_q_generic.h b/include/deal.II/fe/mapping_q_generic.h index aa2e658eef..5b1f2f8409 100644 --- a/include/deal.II/fe/mapping_q_generic.h +++ b/include/deal.II/fe/mapping_q_generic.h @@ -39,6 +39,9 @@ DEAL_II_NAMESPACE_OPEN template class MappingQ; +template +class MappingQCache; + /*!@addtogroup mapping */ /*@{*/ @@ -735,6 +738,13 @@ protected: */ template friend class MappingQ; + + /** + * Make MappingQCache a friend since it needs to call the + * compute_mapping_support_points() function. + */ + template + friend class MappingQCache; }; diff --git a/source/fe/CMakeLists.txt b/source/fe/CMakeLists.txt index 120df5a823..2beb067421 100644 --- a/source/fe/CMakeLists.txt +++ b/source/fe/CMakeLists.txt @@ -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 index 0000000000..55925cb245 --- /dev/null +++ b/source/fe/mapping_q_cache.cc @@ -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 +#include + +#include + +#include + +DEAL_II_NAMESPACE_OPEN + +template +MappingQCache::MappingQCache( + const unsigned int polynomial_degree) + : MappingQGeneric(polynomial_degree) +{} + + + +template +MappingQCache::MappingQCache( + const MappingQCache &mapping) + : MappingQGeneric(mapping) + , support_point_cache(mapping.support_point_cache) +{} + + + +template +MappingQCache::~MappingQCache() +{ + if (clear_signal.connected()) + clear_signal.disconnect(); +} + + + +template +std::unique_ptr> +MappingQCache::clone() const +{ + return std_cxx14::make_unique>(*this); +} + + + +template +bool +MappingQCache::preserves_vertex_locations() const +{ + return false; +} + + + +template +void +MappingQCache::initialize( + const Triangulation & triangulation, + const MappingQGeneric &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>>>>( + 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::cell_iterator &cell, + void *, + void *) { + (*support_point_cache)[cell->level()][cell->index()] = + mapping.compute_mapping_support_points(cell); + }, + /* copier */ std::function(), + /* scratch_data */ nullptr, + /* copy_data */ nullptr, + 2 * MultithreadInfo::n_threads(), + /* chunk_size = */ 1); +} + + + +template +std::size_t +MappingQCache::memory_consumption() const +{ + if (support_point_cache.get() != nullptr) + return sizeof(*this) + + MemoryConsumption::memory_consumption(*support_point_cache); + else + return sizeof(*this); +} + + + +template +std::vector> +MappingQCache::compute_mapping_support_points( + const typename Triangulation::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 index 0000000000..d57d9e87b5 --- /dev/null +++ b/source/fe/mapping_q_cache.inst.in @@ -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; +#endif + } -- 2.39.5