From: wolf Date: Mon, 3 Sep 2001 08:42:16 +0000 (+0000) Subject: Remove function FE_Q::get_renumber. Slight changes necessary elsewhere. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=3536bb62ce0f0ce21ecc3ff2aaadb9079827cbfb;p=dealii-svn.git Remove function FE_Q::get_renumber. Slight changes necessary elsewhere. git-svn-id: https://svn.dealii.org/trunk@4929 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe_q.h b/deal.II/deal.II/include/fe/fe_q.h index edbb6a15bf..40bfd075ff 100644 --- a/deal.II/deal.II/include/fe/fe_q.h +++ b/deal.II/deal.II/include/fe/fe_q.h @@ -283,22 +283,6 @@ class FE_Q : public FiniteElement virtual Tensor<2,dim> shape_grad_grad (const unsigned int i, const Point &p) const; -//TODO:[WB,RH] make get_renumber private or move it some other place - /** - * Read-only access to the - * renumber vector. - * - * This function shouldn't be - * used, i.e. the users code - * shouldn't rely on the actual - * numbering of the degrees of - * dreedom on each cell. This, - * because this internal - * numbering might change in - * future. - */ - const std::vector & get_renumber() const; - /** * Return the polynomial degree * of this finite element, @@ -466,56 +450,86 @@ class FE_Q : public FiniteElement static std::vector get_dpo_vector(const unsigned int degree); /** - * Map tensor product data to shape - * function numbering. + * Map tensor product data to + * shape function numbering. This + * function is actually an alike + * replica of the respective + * function in the @ref{FETools} + * class, but is kept for three + * reasons: * - * The node values are ordered such - * that vertices are first, - * followed by lines, - * quadrilaterals and - * hexahedra. Furthermore, the - * ordering inside each group may - * be confused, too. Therefore, - * this function computes a mapping - * from lexicographic ordering - * (x,y,z) to the shape function - * structure. + * 1. It only operates on a + * @ref{FiniteElementData} + * structure. This is ok in the + * present context, since we can + * control which types of + * arguments it is called with + * because this is a private + * function. However, the + * publicly visible function in + * the @ref{FETools} class needs + * to make sure that the + * @ref{FiniteElementData} object + * it works on actually + * represents a continuous finite + * element, which we found too + * difficult if we do not pass an + * object of type @ref{FE_Q} + * directly. * - * This function is made - * static. This allows other - * classes (like e.g. @p{MappingQ}) - * to call this functions without a - * need to create a - * FiniteElement. This function - * needs some data from the base - * class @p{FiniteElementData} of - * this @p{FE_Q} class. But this - * data cannot be accessed to by a - * static function. Hence this - * function needs an additional - * @p{FiniteElementData} argument. + * 2. If we would call the + * publicly available version of + * this function instead of this + * one, we would have to pass a + * finite element + * object. However, since the + * construction of an entire + * finite element object can be + * costly, we rather chose to + * retain this function. + * + * 3. Third reason is that we + * want to call this function for + * faces as well, by just calling + * this function for the finite + * element of one dimension + * less. If we would call the + * global function instead, this + * would require us to construct + * a second finite element object + * of one dimension less, just to + * call this function. Since that + * function does not make use of + * hanging nodes constraints, + * interpolation and restriction + * matrices, etc, this would have + * been a waste. Furthermore, it + * would have posed problems with + * template instantiations. + * + * To sum up, the existence of + * this function is a compromise + * between simplicity and proper + * library design, where we have + * chosen to weigh the simplicity + * aspect a little more than + * proper design. */ - static void build_renumbering (const FiniteElementData &fe_data, - const unsigned int degree, - std::vector &numbering); + static + void + lexicographic_to_hierarchic_numbering (const FiniteElementData &fe_data, + const unsigned int degree, + std::vector &numbering); /** - * Map tensor product data to - * shape function numbering on - * first face. - * - * This function does the same as - * @p{build_renumbering}, only on - * the first face. - * - * It is used to compute the - * order of face support points. - * - * In 1d, this function does - * nothing. + * This is an analogon to the + * previous function, but working + * on faces. */ - static void build_face_renumbering (const unsigned int degree, - std::vector &numbering); + static + void + face_lexicographic_to_hierarchic_numbering (const unsigned int degree, + std::vector &numbering); /** * Initialize the @@ -605,21 +619,14 @@ class FE_Q : public FiniteElement * Allow access from other dimensions. */ template friend class FE_Q; - - /** - * Allows @p{MappingQ} class to - * access to @p{build_renumbering} - * function. - */ - friend class MappingQ; }; /* -------------- declaration of explicit specializations ------------- */ template <> void FE_Q<1>::initialize_unit_face_support_points (); -template <> void FE_Q<1>::build_face_renumbering (const unsigned int, - std::vector&); +template <> void FE_Q<1>::face_lexicographic_to_hierarchic_numbering (const unsigned int, + std::vector&); template <> const double * const @@ -661,17 +668,6 @@ template <> const unsigned int FE_Q<3>::Matrices::n_constraint_matrices; -/* ---------------------------- inline functions --------------------- */ - -template -inline -const std::vector & -FE_Q::get_renumber() const -{ - return renumber; -} - - #endif diff --git a/deal.II/deal.II/source/fe/fe_q.cc b/deal.II/deal.II/source/fe/fe_q.cc index d6c9b9695f..be6f229e2f 100644 --- a/deal.II/deal.II/source/fe/fe_q.cc +++ b/deal.II/deal.II/source/fe/fe_q.cc @@ -25,7 +25,6 @@ -//TODO:[RH] move build_renumbering to FiniteElementData class template FE_Q::FE_Q (const unsigned int degree) @@ -51,8 +50,8 @@ FE_Q::FE_Q (const unsigned int degree) // do some internal book-keeping on // cells and faces. if in 1d, the // face function is empty - build_renumbering (*this, degree, renumber); - build_face_renumbering (degree, face_renumber); + lexicographic_to_hierarchic_numbering (*this, degree, renumber); + face_lexicographic_to_hierarchic_numbering (degree, face_renumber); // build inverse of renumbering // vector @@ -723,9 +722,9 @@ FE_Q::get_dpo_vector(const unsigned int deg) template void -FE_Q::build_renumbering (const FiniteElementData &fe_data, - const unsigned int degree, - std::vector &renumber) +FE_Q::lexicographic_to_hierarchic_numbering (const FiniteElementData &fe_data, + const unsigned int degree, + std::vector &renumber) { const unsigned int n = degree+1; @@ -948,11 +947,11 @@ FE_Q::build_renumbering (const FiniteElementData &fe_data, template void -FE_Q::build_face_renumbering (const unsigned int degree, - std::vector &numbering) +FE_Q::face_lexicographic_to_hierarchic_numbering (const unsigned int degree, + std::vector &numbering) { FiniteElementData fe_data(FE_Q::get_dpo_vector(degree),1); - FE_Q::build_renumbering (fe_data, degree, numbering); + FE_Q::lexicographic_to_hierarchic_numbering (fe_data, degree, numbering); } diff --git a/deal.II/deal.II/source/fe/mapping_q.cc b/deal.II/deal.II/source/fe/mapping_q.cc index 926440c822..d55e0e821f 100644 --- a/deal.II/deal.II/source/fe/mapping_q.cc +++ b/deal.II/deal.II/source/fe/mapping_q.cc @@ -22,6 +22,7 @@ #include #include #include +#include #include @@ -118,11 +119,8 @@ MappingQ::MappingQ (const unsigned int p): // shape functions of the Qp // mapping. renumber.resize(n_shape_functions,0); - std::vector dpo(dim+1, static_cast(1)); - for (unsigned int i=1; i fe_data(dpo, 1); - FE_Q::build_renumbering (fe_data, p, renumber); + FETools::lexicographic_to_hierarchic_numbering (FE_Q(degree), + renumber); // build laplace_on_quad_vector if (degree>1)