]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Implement FE_Poly::memory_consumption so that memory_consumption works for FE_Q and...
authorPeter Munch <peterrmuench@gmail.com>
Fri, 20 Dec 2019 10:49:44 +0000 (11:49 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Sat, 21 Dec 2019 23:23:54 +0000 (00:23 +0100)
22 files changed:
doc/news/changes/minor/20191222PeterMunch [new file with mode: 0644]
include/deal.II/base/polynomial.h
include/deal.II/base/polynomials_piecewise.h
include/deal.II/base/scalar_polynomials_base.h
include/deal.II/base/tensor_product_polynomials.h
include/deal.II/fe/fe.h
include/deal.II/fe/fe_dgq.h
include/deal.II/fe/fe_poly.h
include/deal.II/fe/fe_poly.templates.h
include/deal.II/fe/mapping.h
include/deal.II/fe/mapping_cartesian.h
include/deal.II/fe/mapping_fe_field.h
include/deal.II/fe/mapping_manifold.h
include/deal.II/fe/mapping_q.h
include/deal.II/fe/mapping_q_generic.h
source/base/polynomial.cc
source/base/polynomials_piecewise.cc
source/base/scalar_polynomials_base.cc
source/base/tensor_product_polynomials.cc
source/fe/fe_dgq.cc
tests/fe/memory_consumption_01.cc [new file with mode: 0644]
tests/fe/memory_consumption_01.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20191222PeterMunch b/doc/news/changes/minor/20191222PeterMunch
new file mode 100644 (file)
index 0000000..d32db2f
--- /dev/null
@@ -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)
index 3d99c7010b027f8e966b13ef2c21615bc1adc5da..bff0d248a936e0c4902aa6be05d8d1cc3ec3ebdd 100644 (file)
@@ -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.
index ee21f0d6d1abccf711fc7d529eca132fd831fa42..cb1659254cfe10ac3dc82dda8ac242f648c635dd 100644 (file)
@@ -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
index 78e263c8300d8f2359e8e18db95a9f157d08a204..d65ef3451e53534e352bea6a71cc6c4e8b8b61cd 100644 (file)
@@ -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.
index 9423032102faff10bb748647928516de03001ddd..25e23d5e543efc9a99202f21c26ab3bbdf3690c3 100644 (file)
@@ -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.
index 583e233a79ccb865d042ab7e60aaef9b001e326e..73cbda8c1d95bc6e368bcd30ac04884ba5df09be 100644 (file)
@@ -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;
index 5d9123e4483dca36b3bb64c9d250eae10fd6312f..26705d2c3c4691303f0b548c2d063a9d1374b7da 100644 (file)
@@ -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;
 
index b971dcbac4ab3a159fd29ea80d9dc1db63013d35..b8f30add18bed2ca4e9499b16b9b68b8ff6e7d05 100644 (file)
@@ -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
index fb7f7948eb9606253f80c0fde10c19398d0271bb..3e705a2fc74811431797da001a1c083d706f1c57 100644 (file)
@@ -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
index 82a4c1a776be42cee17270ebee9f334889079e38..793376a7d116a7c369b9645bb2cb8669ef13b4bf 100644 (file)
@@ -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;
index 89f115096c49db54709028004778938d4008b490..7f0d58cf95ef2ab98159a4aa408832af203f0426 100644 (file)
@@ -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;
index b4d5d5ead196a636c0727a350fee6d27ddcfdc66..4930988eca3c0abe8cdebd4bfbfd4a60bdb6face 100644 (file)
@@ -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;
index b2911598ac568bfeb83d290e4f748f5e82eb2449..828f74d011bde57299da5331b1f7d5465314b16b 100644 (file)
@@ -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;
index 33dd1de1d899a7bac00f749efcf084f583ce2ce1..e66dfdb137b4ec21d35c1ea3cc70894997c6f19f 100644 (file)
@@ -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;
index f4669f6d3359e877236496fcfc5d897be24cb12e..18dedc8bed7e5a9c0bd7dfe5b986f28f5859ae52 100644 (file)
@@ -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;
index 3eacc200a2a8ea8354e9a376dd91276fb1dd7af8..bff1f5c672841ea78014e6a466ffd35b01c03142 100644 (file)
@@ -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 -------------------------- //
 
index 75bc20b8ac821ea9f4cfb0fb5d20c28dbeaadf3a..842f767ed77e0105a600209c9bad83b61b5d2844 100644 (file)
@@ -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,
index f472f5f46de9bc8f6c4392c840b2910c85c74d38..8fe926612c5c7d23bcef42c923d04439aac71bbd 100644 (file)
@@ -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>;
index b1fb82078a4e4c6a89bf974dc3b44b0e3491940b..a4c344b444be68f72cfc440f343124523838d3f4 100644 (file)
@@ -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 -------------- */
 
 
index 5f4fe48a94e1aba15dd3810d45e953b3e87dbc35..802dfcb9b14a3d90405f9b6034d68f0a0961c1aa 100644 (file)
@@ -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 (file)
index 0000000..e0a3a11
--- /dev/null
@@ -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 (file)
index 0000000..2ba6abd
--- /dev/null
@@ -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

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.