From 85782cbd3f29dc97bd2fd6aa006685bcdd1a3a98 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Wed, 17 Sep 2014 08:41:06 -0500 Subject: [PATCH] Minor cleanups. Specifically: Make the FiniteElement::clone() function public, as it is of general use. Also: Move several member variables from private to protected as these are generally variables that derived classes are supposed to initialize; this allows to remove rather random friend declarations for these derived classes. --- include/deal.II/fe/fe.h | 186 ++++++++++++++++++++-------------------- 1 file changed, 94 insertions(+), 92 deletions(-) diff --git a/include/deal.II/fe/fe.h b/include/deal.II/fe/fe.h index c47bcb9b90..61b314dc71 100644 --- a/include/deal.II/fe/fe.h +++ b/include/deal.II/fe/fe.h @@ -30,13 +30,12 @@ DEAL_II_NAMESPACE_OPEN template class FEValuesData; template class FEValuesBase; template class FEValues; -template class FEFaceValues; -template class FESubfaceValues; -template class FESystem; -template class FE_PolyTensor; +template class FEFaceValues; +template class FESubfaceValues; +template class FESystem; namespace hp { - template class FECollection; + template class FECollection; } @@ -408,6 +407,14 @@ public: */ virtual ~FiniteElement (); + /** + * A sort of virtual copy constructor. Some places in the library, for + * example the constructors of FESystem as well as the hp::FECollection + * class, need to make copies of finite elements without knowing their exact + * type. They do so through this function. + */ + virtual FiniteElement *clone() const = 0; + /** * Return a string that uniquely identifies a finite element. The general * convention is that this is the class name, followed by the dimension in @@ -1490,7 +1497,8 @@ public: * For more information, see the documentation for the has_support_points() * function. */ - bool has_generalized_face_support_points () const; + bool + has_generalized_face_support_points () const; /** * Interpolate a set of scalar values, computed in the generalized support @@ -1501,8 +1509,10 @@ public: * just the values in the suport points. All other elements must reimplement * it. */ - virtual void interpolate(std::vector &local_dofs, - const std::vector &values) const; + virtual + void + interpolate(std::vector &local_dofs, + const std::vector &values) const; /** * Interpolate a set of vector values, computed in the generalized support @@ -1513,17 +1523,20 @@ public: * be interpolated. Maybe consider changing your data structures to use the * next function. */ - virtual void interpolate(std::vector &local_dofs, - const std::vector > &values, - unsigned int offset = 0) const; + virtual + void + interpolate(std::vector &local_dofs, + const std::vector > &values, + unsigned int offset = 0) const; /** * Interpolate a set of vector values, computed in the generalized support * points. */ - virtual void interpolate( - std::vector &local_dofs, - const VectorSlice > > &values) const; + virtual + void + interpolate(std::vector &local_dofs, + const VectorSlice > > &values) const; //@} @@ -1536,6 +1549,7 @@ public: * itself. */ virtual std::size_t memory_consumption () const; + /** * Exception * @@ -1753,80 +1767,6 @@ protected: */ std::vector adjust_line_dof_index_for_line_orientation_table; - /** - * Return the size of interface constraint matrices. Since this is needed in - * every derived finite element class when initializing their size, it is - * placed into this function, to avoid having to recompute the - * dimension-dependent size of these matrices each time. - * - * Note that some elements do not implement the interface constraints for - * certain polynomial degrees. In this case, this function still returns the - * size these matrices should have when implemented, but the actual matrices - * are empty. - */ - TableIndices<2> - interface_constraints_size () const; - - /** - * Compute second derivatives by finite differences of gradients. - */ - void compute_2nd (const Mapping &mapping, - const typename Triangulation::cell_iterator &cell, - const unsigned int offset, - typename Mapping::InternalDataBase &mapping_internal, - InternalDataBase &fe_internal, - FEValuesData &data) const; - - /** - * Given the pattern of nonzero components for each shape function, compute - * for each entry how many components are non-zero for each shape - * function. This function is used in the constructor of this class. - */ - static - std::vector - compute_n_nonzero_components (const std::vector &nonzero_components); - - /** - * Determine the values a finite element should compute on initialization of - * data for FEValues. - * - * Given a set of flags indicating what quantities are requested from a - * FEValues object, update_once() and update_each() compute which values - * must really be computed. Then, the fill_*_values functions are - * called with the result of these. - * - * Furthermore, values must be computed either on the unit cell or on the - * physical cell. For instance, the function values of FE_Q do only depend - * on the quadrature points on the unit cell. Therefore, this flags will be - * returned by update_once(). The gradients require computation of the - * covariant transformation matrix. Therefore, @p - * update_covariant_transformation and @p update_gradients will be returned - * by update_each(). - * - * For an example see the same function in the derived class FE_Q. - */ - virtual UpdateFlags update_once (const UpdateFlags flags) const = 0; - - /** - * Complementary function for update_once(). - * - * While update_once() returns the values to be computed on the unit cell - * for yielding the required data, this function determines the values that - * must be recomputed on each cell. - * - * Refer to update_once() for more details. - */ - virtual UpdateFlags update_each (const UpdateFlags flags) const = 0; - - /** - * A sort of virtual copy constructor. Some places in the library, for - * example the constructors of FESystem as well as the hp::FECollection - * class, need to make copies of finite elements without knowing their exact - * type. They do so through this function. - */ - virtual FiniteElement *clone() const = 0; - -private: /** * Store what system_to_component_index() will return. */ @@ -1959,6 +1899,71 @@ private: */ static const double fd_step_length; + /** + * Return the size of interface constraint matrices. Since this is needed in + * every derived finite element class when initializing their size, it is + * placed into this function, to avoid having to recompute the + * dimension-dependent size of these matrices each time. + * + * Note that some elements do not implement the interface constraints for + * certain polynomial degrees. In this case, this function still returns the + * size these matrices should have when implemented, but the actual matrices + * are empty. + */ + TableIndices<2> + interface_constraints_size () const; + + /** + * Compute second derivatives by finite differences of gradients. + */ + void compute_2nd (const Mapping &mapping, + const typename Triangulation::cell_iterator &cell, + const unsigned int offset, + typename Mapping::InternalDataBase &mapping_internal, + InternalDataBase &fe_internal, + FEValuesData &data) const; + + /** + * Given the pattern of nonzero components for each shape function, compute + * for each entry how many components are non-zero for each shape + * function. This function is used in the constructor of this class. + */ + static + std::vector + compute_n_nonzero_components (const std::vector &nonzero_components); + + /** + * Determine the values a finite element should compute on initialization of + * data for FEValues. + * + * Given a set of flags indicating what quantities are requested from a + * FEValues object, update_once() and update_each() compute which values + * must really be computed. Then, the fill_*_values functions are + * called with the result of these. + * + * Furthermore, values must be computed either on the unit cell or on the + * physical cell. For instance, the function values of FE_Q do only depend + * on the quadrature points on the unit cell. Therefore, this flags will be + * returned by update_once(). The gradients require computation of the + * covariant transformation matrix. Therefore, @p + * update_covariant_transformation and @p update_gradients will be returned + * by update_each(). + * + * For an example see the same function in the derived class FE_Q. + */ + virtual UpdateFlags update_once (const UpdateFlags flags) const = 0; + + /** + * Complementary function for update_once(). + * + * While update_once() returns the values to be computed on the unit cell + * for yielding the required data, this function determines the values that + * must be recomputed on each cell. + * + * Refer to update_once() for more details. + */ + virtual UpdateFlags update_each (const UpdateFlags flags) const = 0; + /** * Prepare internal data structures and fill in values independent of the * cell. Returns a pointer to an object of which the caller of this function @@ -2046,10 +2051,7 @@ private: friend class FEValues; friend class FEFaceValues; friend class FESubfaceValues; - template friend class FESystem; - template friend class FE_PolyTensor; - friend class hp::FECollection; - + friend class FESystem; }; -- 2.39.5