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.
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-
bool
is_primitive (const unsigned int i) const;
- /**
- * Import function that is overloaded by the one above and would otherwise
- * be hidden.
- */
- using FiniteElementData<dim>::is_primitive;
-
/**
* Number of base elements in a mixed discretization.
*
*/
const std::vector<unsigned int> 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
+template <int dim, int spacedim>
+inline
+bool
+FiniteElement<dim,spacedim>::is_primitive () const
+{
+ return cached_primitivity;
+}
+
+
+
template <int dim, int spacedim>
inline
bool
*/
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.
* 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;
};
-template <int dim>
-inline
-bool
-FiniteElementData<dim>::is_primitive () const
-{
- return cached_primitivity;
-}
-
-
-template <int dim>
-inline
-void
-FiniteElementData<dim>::set_primitivity (const bool value)
-{
- cached_primitivity = value;
-}
-
-
template <int dim>
inline
const BlockIndices &
// ---------------------------------------------------------------------
//
-// 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.
//
std::vector<ComponentMask> (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<unsigned int>(),
+ 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<unsigned int>(),
- 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));
?
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));
}
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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