From 773fa71887b2667840917b14f425e1b0e45da76f Mon Sep 17 00:00:00 2001
From: Peter Munch <peterrmuench@gmail.com>
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.  
+<br>
+(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<ScalarPolynomialsBase<dim>>
   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 <tt>pols</tt> 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<Vector<double>> &support_point_values,
     std::vector<double> &              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<FiniteElement<dim, spacedim>>
   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<dim> & 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<PolynomialType, dim, spacedim>::get_poly_space_numbering_inverse() const
 
 
 
+template <class PolynomialType, int dim, int spacedim>
+std::size_t
+FE_Poly<PolynomialType, dim, spacedim>::memory_consumption() const
+{
+  return FiniteElement<dim, spacedim>::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<dim> &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<dim, spacedim>::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 <deal.II/base/exceptions.h>
+#include <deal.II/base/memory_consumption.h>
 #include <deal.II/base/point.h>
 #include <deal.II/base/polynomial.h>
 #include <deal.II/base/quadrature_lib.h>
@@ -660,6 +661,17 @@ namespace Polynomials
   }
 
 
+  template <typename number>
+  std::size_t
+  Polynomial<number>::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 <deal.II/base/memory_consumption.h>
 #include <deal.II/base/polynomials_piecewise.h>
 
 
@@ -117,6 +118,18 @@ namespace Polynomials
 
 
 
+  template <typename number>
+  std::size_t
+  PiecewisePolynomial<number>::memory_consumption() const
+  {
+    return (polynomial.memory_consumption() +
+            MemoryConsumption::memory_consumption(n_intervals) +
+            MemoryConsumption::memory_consumption(interval) +
+            MemoryConsumption::memory_consumption(spans_two_intervals));
+  }
+
+
+
   std::vector<PiecewisePolynomial<double>>
   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<dim>::ScalarPolynomialsBase(
 
 
 
+template <int dim>
+std::size_t
+ScalarPolynomialsBase<dim>::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 <deal.II/base/exceptions.h>
+#include <deal.II/base/memory_consumption.h>
 #include <deal.II/base/polynomials_piecewise.h>
 #include <deal.II/base/std_cxx14/memory.h>
 #include <deal.II/base/tensor_product_polynomials.h>
@@ -435,6 +436,17 @@ TensorProductPolynomials<dim, PolynomialType>::clone() const
 
 
 
+template <int dim, typename PolynomialType>
+std::size_t
+TensorProductPolynomials<dim, PolynomialType>::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<dim, spacedim>::get_constant_modes() const
 
 
 
-template <int dim, int spacedim>
-std::size_t
-FE_DGQ<dim, spacedim>::memory_consumption() const
-{
-  Assert(false, ExcNotImplemented());
-  return 0;
-}
-
-
-
 // ------------------------------ FE_DGQArbitraryNodes -----------------------
 
 template <int dim, int spacedim>
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 <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_q.h>
+
+#include "../tests.h"
+
+using namespace dealii;
+
+template <int dim>
+void
+test(const FiniteElement<dim> &fe)
+{
+  deallog << fe.get_name() << " " << fe.memory_consumption() << std::endl;
+}
+
+template <int dim>
+void
+tests()
+{
+  // FE_Q
+  for (unsigned int i = 1; i < 4; i++)
+    test(FE_Q<dim>(i));
+
+  // FE_DGQ
+  for (unsigned int i = 0; i < 4; i++)
+    test(FE_DGQ<dim>(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