From 773fa71887b2667840917b14f425e1b0e45da76f Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Fri, 20 Dec 2019 11:49:44 +0100 Subject: [PATCH] Implement FE_Poly::memory_consumption so that memory_consumption works for FE_Q and FE_DGQ --- doc/news/changes/minor/20191222PeterMunch | 4 ++ include/deal.II/base/polynomial.h | 6 ++ include/deal.II/base/polynomials_piecewise.h | 6 ++ .../deal.II/base/scalar_polynomials_base.h | 6 ++ .../deal.II/base/tensor_product_polynomials.h | 6 ++ include/deal.II/fe/fe.h | 2 +- include/deal.II/fe/fe_dgq.h | 11 --- include/deal.II/fe/fe_poly.h | 6 ++ include/deal.II/fe/fe_poly.templates.h | 10 +++ include/deal.II/fe/mapping.h | 2 +- include/deal.II/fe/mapping_cartesian.h | 2 +- include/deal.II/fe/mapping_fe_field.h | 2 +- include/deal.II/fe/mapping_manifold.h | 2 +- include/deal.II/fe/mapping_q.h | 2 +- include/deal.II/fe/mapping_q_generic.h | 2 +- source/base/polynomial.cc | 12 ++++ source/base/polynomials_piecewise.cc | 13 ++++ source/base/scalar_polynomials_base.cc | 10 +++ source/base/tensor_product_polynomials.cc | 12 ++++ source/fe/fe_dgq.cc | 10 --- tests/fe/memory_consumption_01.cc | 67 +++++++++++++++++++ tests/fe/memory_consumption_01.output | 22 ++++++ 22 files changed, 187 insertions(+), 28 deletions(-) create mode 100644 doc/news/changes/minor/20191222PeterMunch create mode 100644 tests/fe/memory_consumption_01.cc create mode 100644 tests/fe/memory_consumption_01.output diff --git a/doc/news/changes/minor/20191222PeterMunch b/doc/news/changes/minor/20191222PeterMunch new file mode 100644 index 0000000000..d32db2f09d --- /dev/null +++ b/doc/news/changes/minor/20191222PeterMunch @@ -0,0 +1,4 @@ +New: The method memory_consumption has been implemented in and added +to a variety of classes, e.g. FE_DGQ. +
+(Peter Munch, 2019/12/22) diff --git a/include/deal.II/base/polynomial.h b/include/deal.II/base/polynomial.h index 3d99c7010b..bff0d248a9 100644 --- a/include/deal.II/base/polynomial.h +++ b/include/deal.II/base/polynomial.h @@ -225,6 +225,12 @@ namespace Polynomials void serialize(Archive &ar, const unsigned int version); + /** + * Return an estimate (in bytes) for the memory consumption of this object. + */ + virtual std::size_t + memory_consumption() const; + protected: /** * This function performs the actual scaling. diff --git a/include/deal.II/base/polynomials_piecewise.h b/include/deal.II/base/polynomials_piecewise.h index ee21f0d6d1..cb1659254c 100644 --- a/include/deal.II/base/polynomials_piecewise.h +++ b/include/deal.II/base/polynomials_piecewise.h @@ -132,6 +132,12 @@ namespace Polynomials void serialize(Archive &ar, const unsigned int version); + /** + * Return an estimate (in bytes) for the memory consumption of this object. + */ + virtual std::size_t + memory_consumption() const; + protected: /** * Underlying polynomial object that is scaled to a subinterval and diff --git a/include/deal.II/base/scalar_polynomials_base.h b/include/deal.II/base/scalar_polynomials_base.h index 78e263c830..d65ef3451e 100644 --- a/include/deal.II/base/scalar_polynomials_base.h +++ b/include/deal.II/base/scalar_polynomials_base.h @@ -138,6 +138,12 @@ public: virtual std::string name() const = 0; + /** + * Return an estimate (in bytes) for the memory consumption of this object. + */ + virtual std::size_t + memory_consumption() const; + private: /** * The highest polynomial degree of this functions represented by this object. diff --git a/include/deal.II/base/tensor_product_polynomials.h b/include/deal.II/base/tensor_product_polynomials.h index 9423032102..25e23d5e54 100644 --- a/include/deal.II/base/tensor_product_polynomials.h +++ b/include/deal.II/base/tensor_product_polynomials.h @@ -209,6 +209,12 @@ public: virtual std::unique_ptr> clone() const override; + /** + * Return an estimate (in bytes) for the memory consumption of this object. + */ + virtual std::size_t + memory_consumption() const override; + protected: /** * Copy of the vector pols of polynomials given to the constructor. diff --git a/include/deal.II/fe/fe.h b/include/deal.II/fe/fe.h index 583e233a79..73cbda8c1d 100644 --- a/include/deal.II/fe/fe.h +++ b/include/deal.II/fe/fe.h @@ -715,7 +715,7 @@ public: UpdateFlags update_each; /** - * Return an estimate (in bytes) or the memory consumption of this object. + * Return an estimate (in bytes) for the memory consumption of this object. */ virtual std::size_t memory_consumption() const; diff --git a/include/deal.II/fe/fe_dgq.h b/include/deal.II/fe/fe_dgq.h index 5d9123e448..26705d2c3c 100644 --- a/include/deal.II/fe/fe_dgq.h +++ b/include/deal.II/fe/fe_dgq.h @@ -320,17 +320,6 @@ public: const std::vector> &support_point_values, std::vector & nodal_values) const override; - /** - * Determine an estimate for the memory consumption (in bytes) of this - * object. - * - * This function is made virtual, since finite element objects are usually - * accessed through pointers to their base class, rather than the class - * itself. - */ - virtual std::size_t - memory_consumption() const override; - virtual std::unique_ptr> clone() const override; diff --git a/include/deal.II/fe/fe_poly.h b/include/deal.II/fe/fe_poly.h index b971dcbac4..b8f30add18 100644 --- a/include/deal.II/fe/fe_poly.h +++ b/include/deal.II/fe/fe_poly.h @@ -227,6 +227,12 @@ public: const Point & p, const unsigned int component) const override; + /** + * Return an estimate (in bytes) for the memory consumption of this object. + */ + virtual std::size_t + memory_consumption() const override; + protected: /* * NOTE: The following function has its definition inlined into the class diff --git a/include/deal.II/fe/fe_poly.templates.h b/include/deal.II/fe/fe_poly.templates.h index fb7f7948eb..3e705a2fc7 100644 --- a/include/deal.II/fe/fe_poly.templates.h +++ b/include/deal.II/fe/fe_poly.templates.h @@ -581,6 +581,16 @@ FE_Poly::get_poly_space_numbering_inverse() const +template +std::size_t +FE_Poly::memory_consumption() const +{ + return FiniteElement::memory_consumption() + + poly_space.memory_consumption(); +} + + + DEAL_II_NAMESPACE_CLOSE #endif diff --git a/include/deal.II/fe/mapping.h b/include/deal.II/fe/mapping.h index 82a4c1a776..793376a7d1 100644 --- a/include/deal.II/fe/mapping.h +++ b/include/deal.II/fe/mapping.h @@ -630,7 +630,7 @@ public: UpdateFlags update_each; /** - * Return an estimate (in bytes) or the memory consumption of this object. + * Return an estimate (in bytes) for the memory consumption of this object. */ virtual std::size_t memory_consumption() const; diff --git a/include/deal.II/fe/mapping_cartesian.h b/include/deal.II/fe/mapping_cartesian.h index 89f115096c..7f0d58cf95 100644 --- a/include/deal.II/fe/mapping_cartesian.h +++ b/include/deal.II/fe/mapping_cartesian.h @@ -182,7 +182,7 @@ private: InternalData(const Quadrature &quadrature); /** - * Return an estimate (in bytes) or the memory consumption of this object. + * Return an estimate (in bytes) for the memory consumption of this object. */ virtual std::size_t memory_consumption() const override; diff --git a/include/deal.II/fe/mapping_fe_field.h b/include/deal.II/fe/mapping_fe_field.h index b4d5d5ead1..4930988eca 100644 --- a/include/deal.II/fe/mapping_fe_field.h +++ b/include/deal.II/fe/mapping_fe_field.h @@ -355,7 +355,7 @@ public: fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr); /** - * Return an estimate (in bytes) or the memory consumption of this object. + * Return an estimate (in bytes) for the memory consumption of this object. */ virtual std::size_t memory_consumption() const override; diff --git a/include/deal.II/fe/mapping_manifold.h b/include/deal.II/fe/mapping_manifold.h index b2911598ac..828f74d011 100644 --- a/include/deal.II/fe/mapping_manifold.h +++ b/include/deal.II/fe/mapping_manifold.h @@ -211,7 +211,7 @@ public: const typename Triangulation::cell_iterator &cell) const; /** - * Return an estimate (in bytes) or the memory consumption of this object. + * Return an estimate (in bytes) for the memory consumption of this object. */ virtual std::size_t memory_consumption() const override; diff --git a/include/deal.II/fe/mapping_q.h b/include/deal.II/fe/mapping_q.h index 33dd1de1d8..e66dfdb137 100644 --- a/include/deal.II/fe/mapping_q.h +++ b/include/deal.II/fe/mapping_q.h @@ -237,7 +237,7 @@ protected: /** - * Return an estimate (in bytes) or the memory consumption of this object. + * Return an estimate (in bytes) for the memory consumption of this object. */ virtual std::size_t memory_consumption() const override; diff --git a/include/deal.II/fe/mapping_q_generic.h b/include/deal.II/fe/mapping_q_generic.h index f4669f6d33..18dedc8bed 100644 --- a/include/deal.II/fe/mapping_q_generic.h +++ b/include/deal.II/fe/mapping_q_generic.h @@ -364,7 +364,7 @@ public: fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr); /** - * Return an estimate (in bytes) or the memory consumption of this object. + * Return an estimate (in bytes) for the memory consumption of this object. */ virtual std::size_t memory_consumption() const override; diff --git a/source/base/polynomial.cc b/source/base/polynomial.cc index 3eacc200a2..bff1f5c672 100644 --- a/source/base/polynomial.cc +++ b/source/base/polynomial.cc @@ -14,6 +14,7 @@ // --------------------------------------------------------------------- #include +#include #include #include #include @@ -660,6 +661,17 @@ namespace Polynomials } + template + std::size_t + Polynomial::memory_consumption() const + { + return (MemoryConsumption::memory_consumption(coefficients) + + MemoryConsumption::memory_consumption(in_lagrange_product_form) + + MemoryConsumption::memory_consumption(lagrange_support_points) + + MemoryConsumption::memory_consumption(lagrange_weight)); + } + + // ------------------ class Monomial -------------------------- // diff --git a/source/base/polynomials_piecewise.cc b/source/base/polynomials_piecewise.cc index 75bc20b8ac..842f767ed7 100644 --- a/source/base/polynomials_piecewise.cc +++ b/source/base/polynomials_piecewise.cc @@ -13,6 +13,7 @@ // // --------------------------------------------------------------------- +#include #include @@ -117,6 +118,18 @@ namespace Polynomials + template + std::size_t + PiecewisePolynomial::memory_consumption() const + { + return (polynomial.memory_consumption() + + MemoryConsumption::memory_consumption(n_intervals) + + MemoryConsumption::memory_consumption(interval) + + MemoryConsumption::memory_consumption(spans_two_intervals)); + } + + + std::vector> generate_complete_Lagrange_basis_on_subdivisions( const unsigned int n_subdivisions, diff --git a/source/base/scalar_polynomials_base.cc b/source/base/scalar_polynomials_base.cc index f472f5f46d..8fe926612c 100644 --- a/source/base/scalar_polynomials_base.cc +++ b/source/base/scalar_polynomials_base.cc @@ -34,6 +34,16 @@ ScalarPolynomialsBase::ScalarPolynomialsBase( +template +std::size_t +ScalarPolynomialsBase::memory_consumption() const +{ + Assert(false, ExcNotImplemented()); + return 0; +} + + + template class ScalarPolynomialsBase<0>; template class ScalarPolynomialsBase<1>; template class ScalarPolynomialsBase<2>; diff --git a/source/base/tensor_product_polynomials.cc b/source/base/tensor_product_polynomials.cc index b1fb82078a..a4c344b444 100644 --- a/source/base/tensor_product_polynomials.cc +++ b/source/base/tensor_product_polynomials.cc @@ -14,6 +14,7 @@ // --------------------------------------------------------------------- #include +#include #include #include #include @@ -435,6 +436,17 @@ TensorProductPolynomials::clone() const +template +std::size_t +TensorProductPolynomials::memory_consumption() const +{ + return (MemoryConsumption::memory_consumption(polynomials) + + MemoryConsumption::memory_consumption(index_map) + + MemoryConsumption::memory_consumption(index_map_inverse)); +} + + + /* ------------------- AnisotropicPolynomials -------------- */ diff --git a/source/fe/fe_dgq.cc b/source/fe/fe_dgq.cc index 5f4fe48a94..802dfcb9b1 100644 --- a/source/fe/fe_dgq.cc +++ b/source/fe/fe_dgq.cc @@ -819,16 +819,6 @@ FE_DGQ::get_constant_modes() const -template -std::size_t -FE_DGQ::memory_consumption() const -{ - Assert(false, ExcNotImplemented()); - return 0; -} - - - // ------------------------------ FE_DGQArbitraryNodes ----------------------- template diff --git a/tests/fe/memory_consumption_01.cc b/tests/fe/memory_consumption_01.cc new file mode 100644 index 0000000000..e0a3a11ca3 --- /dev/null +++ b/tests/fe/memory_consumption_01.cc @@ -0,0 +1,67 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + +// Test memory_consumption of FE_Q and FE_DGQ. + + +#include +#include + +#include "../tests.h" + +using namespace dealii; + +template +void +test(const FiniteElement &fe) +{ + deallog << fe.get_name() << " " << fe.memory_consumption() << std::endl; +} + +template +void +tests() +{ + // FE_Q + for (unsigned int i = 1; i < 4; i++) + test(FE_Q(i)); + + // FE_DGQ + for (unsigned int i = 0; i < 4; i++) + test(FE_DGQ(i)); +} + +int +main() +{ + initlog(); + + { + deallog.push("1d"); + tests<1>(); + deallog.pop(); + } + { + deallog.push("2d"); + tests<2>(); + deallog.pop(); + } + { + deallog.push("3d"); + tests<3>(); + deallog.pop(); + } +} \ No newline at end of file diff --git a/tests/fe/memory_consumption_01.output b/tests/fe/memory_consumption_01.output new file mode 100644 index 0000000000..2ba6abdecd --- /dev/null +++ b/tests/fe/memory_consumption_01.output @@ -0,0 +1,22 @@ + +DEAL:1d::FE_Q<1>(1) 1770 +DEAL:1d::FE_Q<1>(2) 1979 +DEAL:1d::FE_Q<1>(3) 2204 +DEAL:1d::FE_DGQ<1>(0) 1557 +DEAL:1d::FE_DGQ<1>(1) 1750 +DEAL:1d::FE_DGQ<1>(2) 1959 +DEAL:1d::FE_DGQ<1>(3) 2184 +DEAL:2d::FE_Q<2>(1) 4062 +DEAL:2d::FE_Q<2>(2) 4827 +DEAL:2d::FE_Q<2>(3) 5880 +DEAL:2d::FE_DGQ<2>(0) 3573 +DEAL:2d::FE_DGQ<2>(1) 4006 +DEAL:2d::FE_DGQ<2>(2) 4695 +DEAL:2d::FE_DGQ<2>(3) 5640 +DEAL:3d::FE_Q<3>(1) 10678 +DEAL:3d::FE_Q<3>(2) 14499 +DEAL:3d::FE_Q<3>(3) 23432 +DEAL:3d::FE_DGQ<3>(0) 9525 +DEAL:3d::FE_DGQ<3>(1) 10438 +DEAL:3d::FE_DGQ<3>(2) 12807 +DEAL:3d::FE_DGQ<3>(3) 17352 -- 2.39.5