From: Wolfgang Bangerth Date: Sun, 11 Sep 2016 10:29:48 +0000 (-0500) Subject: Move 'cached_primitivity' and 'is_primitive()'. X-Git-Tag: v8.5.0-rc1~677^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=7fd8c20e20ee5be2a6b2f0623f567224c180faaf;p=dealii.git Move 'cached_primitivity' and 'is_primitive()'. These fields were part of FiniteElementData, but this class had no way of computing the information during construction, and it needed to be set at a later time by the derived classes' construction time. This is silly. Move the field and function to class FiniteElement, where we know how to compute the information, and where we already store similar information. --- diff --git a/include/deal.II/fe/fe.h b/include/deal.II/fe/fe.h index d0a0a904ca..a79dacc992 100644 --- a/include/deal.II/fe/fe.h +++ b/include/deal.II/fe/fe.h @@ -1426,6 +1426,16 @@ public: unsigned int n_nonzero_components (const unsigned int i) const; + /** + * Return whether the entire finite element is primitive, in the sense that + * all its shape functions are primitive. If the finite element is scalar, + * then this is always the case. + * + * Since this is an extremely common operation, the result is cached and + * returned by this function. + */ + bool is_primitive () const; + /** * Return whether the @p ith shape function is primitive in the sense that * the shape function is non-zero in only one vector component. Non- @@ -1439,12 +1449,6 @@ public: bool is_primitive (const unsigned int i) const; - /** - * Import function that is overloaded by the one above and would otherwise - * be hidden. - */ - using FiniteElementData::is_primitive; - /** * Number of base elements in a mixed discretization. * @@ -2278,6 +2282,13 @@ protected: */ const std::vector n_nonzero_components_table; + /** + * Store whether all shape functions are primitive. Since finding this out + * is a very common operation, we cache the result, i.e. compute the value + * in the constructor for simpler access. + */ + const bool cached_primitivity; + /** * Return the size of interface constraint matrices. Since this is needed in * every derived finite element class when initializing their size, it is @@ -2925,6 +2936,16 @@ FiniteElement::n_nonzero_components (const unsigned int i) const +template +inline +bool +FiniteElement::is_primitive () const +{ + return cached_primitivity; +} + + + template inline bool diff --git a/include/deal.II/fe/fe_base.h b/include/deal.II/fe/fe_base.h index c52019fa30..97dc0f2835 100644 --- a/include/deal.II/fe/fe_base.h +++ b/include/deal.II/fe/fe_base.h @@ -413,17 +413,6 @@ public: */ const BlockIndices &block_indices() const; - /** - * Return whether the entire finite element is primitive, in the sense that - * all its shape functions are primitive. If the finite element is scalar, - * then this is always the case. - * - * Since this is an extremely common operation, the result is cached in the - * #cached_primitivity variable which is computed in the constructor of the - * derived class FiniteElement. - */ - bool is_primitive () const; - /** * Maximal polynomial degree of a shape function in a single coordinate * direction. @@ -444,24 +433,6 @@ public: * Comparison operator. */ bool operator == (const FiniteElementData &) const; - -protected: - - /** - * Set the primitivity of the element. This is usually done by the - * constructor of a derived class. See - * @ref GlossPrimitive "primitive" - * for details. - */ - void set_primitivity(const bool value); - -private: - /** - * Store whether all shape functions are primitive. Since finding this out - * is a very common operation, we cache the result, i.e. compute the value - * in the constructor for simpler access. - */ - bool cached_primitivity; }; @@ -617,24 +588,6 @@ FiniteElementData::n_components () const -template -inline -bool -FiniteElementData::is_primitive () const -{ - return cached_primitivity; -} - - -template -inline -void -FiniteElementData::set_primitivity (const bool value) -{ - cached_primitivity = value; -} - - template inline const BlockIndices & diff --git a/source/fe/fe.cc b/source/fe/fe.cc index c55b400ae6..1a149f8360 100644 --- a/source/fe/fe.cc +++ b/source/fe/fe.cc @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 1998 - 2015 by the deal.II authors +// Copyright (C) 1998 - 2016 by the deal.II authors // // This file is part of the deal.II library. // @@ -85,16 +85,14 @@ FiniteElement (const FiniteElementData &fe_data, std::vector (fe_data.dofs_per_cell, nonzero_c[0]) : nonzero_c), - n_nonzero_components_table (compute_n_nonzero_components(nonzero_components)) + n_nonzero_components_table (compute_n_nonzero_components(nonzero_components)), + cached_primitivity (std::find_if (n_nonzero_components_table.begin(), + n_nonzero_components_table.end(), + std_cxx11::bind (std::not_equal_to(), + std_cxx11::_1, + 1U)) + == n_nonzero_components_table.end()) { - this->set_primitivity(std::find_if (n_nonzero_components_table.begin(), - n_nonzero_components_table.end(), - std_cxx11::bind (std::not_equal_to(), - std_cxx11::_1, - 1U)) - == n_nonzero_components_table.end()); - - Assert (restriction_is_additive_flags.size() == this->dofs_per_cell, ExcDimensionMismatch(restriction_is_additive_flags.size(), this->dofs_per_cell)); diff --git a/source/fe/fe_data.cc b/source/fe/fe_data.cc index 672cf84930..7eae662dda 100644 --- a/source/fe/fe_data.cc +++ b/source/fe/fe_data.cc @@ -61,12 +61,7 @@ FiniteElementData (const std::vector &dofs_per_object, ? BlockIndices(1, dofs_per_cell) : - block_indices), - // the following field is only set in the FiniteElement - // constructor by calling set_primitivity() of this base - // class. however, to ensure that we always have boolean - // values, initialize it with the safe choice - cached_primitivity (false) + block_indices) { Assert(dofs_per_object.size()==dim+1, ExcDimensionMismatch(dofs_per_object.size()-1,dim)); } diff --git a/tests/base/memory_consumption_01.with_64bit_indices=off.output b/tests/base/memory_consumption_01.with_64bit_indices=off.output index 3607c76697..64f9ba3471 100644 --- a/tests/base/memory_consumption_01.with_64bit_indices=off.output +++ b/tests/base/memory_consumption_01.with_64bit_indices=off.output @@ -5,7 +5,7 @@ DEAL:: Number of active cells: 8 DEAL:: Number of degrees of freedom: 17 DEAL::Memory consumption -- Triangulation: 3090 DEAL::Memory consumption -- DoFHandler: 976 -DEAL::Memory consumption -- FE: 1488 +DEAL::Memory consumption -- FE: 1480 DEAL::Memory consumption -- Constraints: 78 DEAL::Memory consumption -- Sparsity: 456 DEAL::Memory consumption -- Matrix: 632 @@ -16,7 +16,7 @@ DEAL:: Number of active cells: 11 DEAL:: Number of degrees of freedom: 23 DEAL::Memory consumption -- Triangulation: 4026 DEAL::Memory consumption -- DoFHandler: 1164 -DEAL::Memory consumption -- FE: 1488 +DEAL::Memory consumption -- FE: 1480 DEAL::Memory consumption -- Constraints: 78 DEAL::Memory consumption -- Sparsity: 576 DEAL::Memory consumption -- Matrix: 824 @@ -27,7 +27,7 @@ DEAL:: Number of active cells: 15 DEAL:: Number of degrees of freedom: 31 DEAL::Memory consumption -- Triangulation: 4946 DEAL::Memory consumption -- DoFHandler: 1356 -DEAL::Memory consumption -- FE: 1488 +DEAL::Memory consumption -- FE: 1480 DEAL::Memory consumption -- Constraints: 78 DEAL::Memory consumption -- Sparsity: 736 DEAL::Memory consumption -- Matrix: 1080 @@ -38,7 +38,7 @@ DEAL:: Number of active cells: 20 DEAL:: Number of degrees of freedom: 41 DEAL::Memory consumption -- Triangulation: 6090 DEAL::Memory consumption -- DoFHandler: 1584 -DEAL::Memory consumption -- FE: 1488 +DEAL::Memory consumption -- FE: 1480 DEAL::Memory consumption -- Constraints: 78 DEAL::Memory consumption -- Sparsity: 936 DEAL::Memory consumption -- Matrix: 1400 diff --git a/tests/base/memory_consumption_01.with_64bit_indices=on.output b/tests/base/memory_consumption_01.with_64bit_indices=on.output index 751cfa8fa2..53472c326e 100644 --- a/tests/base/memory_consumption_01.with_64bit_indices=on.output +++ b/tests/base/memory_consumption_01.with_64bit_indices=on.output @@ -5,7 +5,7 @@ DEAL:: Number of active cells: 8 DEAL:: Number of degrees of freedom: 17 DEAL::Memory consumption -- Triangulation: 3090 DEAL::Memory consumption -- DoFHandler: 1280 -DEAL::Memory consumption -- FE: 1488 +DEAL::Memory consumption -- FE: 1480 DEAL::Memory consumption -- Constraints: 82 DEAL::Memory consumption -- Sparsity: 792 DEAL::Memory consumption -- Matrix: 632 @@ -16,7 +16,7 @@ DEAL:: Number of active cells: 11 DEAL:: Number of degrees of freedom: 23 DEAL::Memory consumption -- Triangulation: 4026 DEAL::Memory consumption -- DoFHandler: 1576 -DEAL::Memory consumption -- FE: 1488 +DEAL::Memory consumption -- FE: 1480 DEAL::Memory consumption -- Constraints: 82 DEAL::Memory consumption -- Sparsity: 1032 DEAL::Memory consumption -- Matrix: 824 @@ -27,7 +27,7 @@ DEAL:: Number of active cells: 15 DEAL:: Number of degrees of freedom: 31 DEAL::Memory consumption -- Triangulation: 4946 DEAL::Memory consumption -- DoFHandler: 1912 -DEAL::Memory consumption -- FE: 1488 +DEAL::Memory consumption -- FE: 1480 DEAL::Memory consumption -- Constraints: 82 DEAL::Memory consumption -- Sparsity: 1352 DEAL::Memory consumption -- Matrix: 1080 @@ -38,7 +38,7 @@ DEAL:: Number of active cells: 20 DEAL:: Number of degrees of freedom: 41 DEAL::Memory consumption -- Triangulation: 6090 DEAL::Memory consumption -- DoFHandler: 2320 -DEAL::Memory consumption -- FE: 1488 +DEAL::Memory consumption -- FE: 1480 DEAL::Memory consumption -- Constraints: 82 DEAL::Memory consumption -- Sparsity: 1752 DEAL::Memory consumption -- Matrix: 1400 @@ -51,7 +51,7 @@ DEAL:: Number of active cells: 20 DEAL:: Number of degrees of freedom: 89 DEAL::Memory consumption -- Triangulation: 5431 DEAL::Memory consumption -- DoFHandler: 3288 -DEAL::Memory consumption -- FE: 3904 +DEAL::Memory consumption -- FE: 3888 DEAL::Memory consumption -- Constraints: 82 DEAL::Memory consumption -- Sparsity: 11352 DEAL::Memory consumption -- Matrix: 10616 @@ -62,7 +62,7 @@ DEAL:: Number of active cells: 44 DEAL:: Number of degrees of freedom: 209 DEAL::Memory consumption -- Triangulation: 11215 DEAL::Memory consumption -- DoFHandler: 6808 -DEAL::Memory consumption -- FE: 3904 +DEAL::Memory consumption -- FE: 3888 DEAL::Memory consumption -- Constraints: 4594 DEAL::Memory consumption -- Sparsity: 27864 DEAL::Memory consumption -- Matrix: 26168 @@ -73,7 +73,7 @@ DEAL:: Number of active cells: 92 DEAL:: Number of degrees of freedom: 449 DEAL::Memory consumption -- Triangulation: 21927 DEAL::Memory consumption -- DoFHandler: 13672 -DEAL::Memory consumption -- FE: 3904 +DEAL::Memory consumption -- FE: 3888 DEAL::Memory consumption -- Constraints: 15250 DEAL::Memory consumption -- Sparsity: 62360 DEAL::Memory consumption -- Matrix: 58744 @@ -84,7 +84,7 @@ DEAL:: Number of active cells: 200 DEAL:: Number of degrees of freedom: 921 DEAL::Memory consumption -- Triangulation: 44055 DEAL::Memory consumption -- DoFHandler: 28648 -DEAL::Memory consumption -- FE: 3904 +DEAL::Memory consumption -- FE: 3888 DEAL::Memory consumption -- Constraints: 28082 DEAL::Memory consumption -- Sparsity: 128216 DEAL::Memory consumption -- Matrix: 120824 @@ -97,7 +97,7 @@ DEAL:: Number of active cells: 56 DEAL:: Number of degrees of freedom: 517 DEAL::Memory consumption -- Triangulation: 25043 DEAL::Memory consumption -- DoFHandler: 18912 -DEAL::Memory consumption -- FE: 12280 +DEAL::Memory consumption -- FE: 12264 DEAL::Memory consumption -- Constraints: 82 DEAL::Memory consumption -- Sparsity: 240440 DEAL::Memory consumption -- Matrix: 236280 @@ -108,7 +108,7 @@ DEAL:: Number of active cells: 217 DEAL:: Number of degrees of freedom: 2217 DEAL::Memory consumption -- Triangulation: 98569 DEAL::Memory consumption -- DoFHandler: 75320 -DEAL::Memory consumption -- FE: 12280 +DEAL::Memory consumption -- FE: 12264 DEAL::Memory consumption -- Constraints: 78450 DEAL::Memory consumption -- Sparsity: 1134184 DEAL::Memory consumption -- Matrix: 1116424 @@ -119,7 +119,7 @@ DEAL:: Number of active cells: 896 DEAL:: Number of degrees of freedom: 9373 DEAL::Memory consumption -- Triangulation: 395091 DEAL::Memory consumption -- DoFHandler: 301840 -DEAL::Memory consumption -- FE: 12280 +DEAL::Memory consumption -- FE: 12264 DEAL::Memory consumption -- Constraints: 487314 DEAL::Memory consumption -- Sparsity: 5034296 DEAL::Memory consumption -- Matrix: 4959288 @@ -130,7 +130,7 @@ DEAL:: Number of active cells: 3248 DEAL:: Number of degrees of freedom: 32433 DEAL::Memory consumption -- Triangulation: 1390359 DEAL::Memory consumption -- DoFHandler: 1083216 -DEAL::Memory consumption -- FE: 12280 +DEAL::Memory consumption -- FE: 12264 DEAL::Memory consumption -- Constraints: 1729778 DEAL::Memory consumption -- Sparsity: 16643416 DEAL::Memory consumption -- Matrix: 16383928