From: hartmann Date: Mon, 13 Jun 2005 11:49:17 +0000 (+0000) Subject: To avoid code duplication, use FETools::lexicographic_to_hierarchic_numbering functio... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=520b024538884838a4e2ffdd0e2b5f35f0e96136;p=dealii-svn.git To avoid code duplication, use FETools::lexicographic_to_hierarchic_numbering function instead of reimplementing this functionality in FE_Q. git-svn-id: https://svn.dealii.org/trunk@10858 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 a18e4cabf4..8ca83af781 100644 --- a/deal.II/deal.II/include/fe/fe_q.h +++ b/deal.II/deal.II/include/fe/fe_q.h @@ -328,85 +328,12 @@ class FE_Q : public FE_Poly,dim> * @p FiniteElementData. */ static std::vector get_dpo_vector(const unsigned int degree); - - /** - * Map tensor product data to - * shape function numbering. This - * function is actually an alike - * replica of the respective - * function in the FETools - * class, but is kept for three - * reasons: - * - * 1. It only operates on a - * 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 FETools class needs - * to make sure that the - * 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 FE_Q() - * directly. - * - * 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. - * - * This function is called from - * the constructor. - */ - static - std::vector - lexicographic_to_hierarchic_numbering (const FiniteElementData &fe_data, - const unsigned int degree); /** * This is an analogon to the - * previous function, but working - * on faces. Called from the + * FETools::lexicographic_to_hierarchic_numbering + * function, but working on + * faces. Called from the * constructor. */ static diff --git a/deal.II/deal.II/source/fe/fe_q.cc b/deal.II/deal.II/source/fe/fe_q.cc index 22e879c1ba..e02fc15f0e 100644 --- a/deal.II/deal.II/source/fe/fe_q.cc +++ b/deal.II/deal.II/source/fe/fe_q.cc @@ -12,6 +12,7 @@ //--------------------------------------------------------------------------- #include +#include #ifdef HAVE_STD_STRINGSTREAM # include @@ -431,9 +432,10 @@ FE_Q::FE_Q (const unsigned int degree) get_dpo_vector(degree),1, degree).dofs_per_cell, std::vector(1,true))), face_index_map(FE_Q_Helper::invert_numbering(face_lexicographic_to_hierarchic_numbering (degree))) { - this->poly_space.set_numbering(FE_Q_Helper::invert_numbering( - lexicographic_to_hierarchic_numbering (*this, degree))); - + std::vector renumber (this->dofs_per_cell); + FETools::lexicographic_to_hierarchic_numbering (*this, renumber); + this->poly_space.set_numbering(FE_Q_Helper::invert_numbering(renumber)); + // finally fill in support points // on cell and face initialize_unit_support_points (); @@ -665,245 +667,14 @@ FE_Q::get_dpo_vector(const unsigned int deg) } - -template -std::vector -FE_Q::lexicographic_to_hierarchic_numbering (const FiniteElementData &fe_data, - const unsigned int degree) -{ - std::vector renumber (fe_data.dofs_per_cell); - - const unsigned int n = degree+1; - - if (degree == 0) - { - Assert ((fe_data.dofs_per_vertex == 0) && - ((fe_data.dofs_per_line == 0) || (dim == 1)) && - ((fe_data.dofs_per_quad == 0) || (dim == 2)) && - ((fe_data.dofs_per_hex == 0) || (dim == 3)), - ExcInternalError()); - renumber[0] = 0; - }; - - if (degree > 0) - { - Assert (fe_data.dofs_per_vertex == 1, ExcInternalError()); - for (unsigned int i=0; i::vertices_per_cell; ++i) - { - unsigned int index = 0; - // Find indices of vertices. - // Unfortunately, somebody - // switched the upper corner - // points of a quad. The same - // person decided to find a very - // creative numbering of the - // vertices of a hexahedron. - // Therefore, this looks quite - // sophisticated. - // - // NB: This same person - // claims to have had good - // reasons then, but seems to - // have forgotten about - // them. At least, the - // numbering was discussed - // with the complaining - // person back then when all - // began :-) - switch (dim) - { - case 1: - { - const unsigned int values[GeometryInfo<1>::vertices_per_cell] - = { 0, degree }; - index = values[i]; - break; - }; - - case 2: - { - const unsigned int values[GeometryInfo<2>::vertices_per_cell] - = { 0, degree, n*degree+degree, n*degree }; - index = values[i]; - break; - }; - - case 3: - { - const unsigned int values[GeometryInfo<3>::vertices_per_cell] - = { 0, degree, - n*n*degree + degree, n*n*degree, - n*degree, n*degree+degree, - n*n*degree + n*degree+degree, n*n*degree + n*degree}; - index = values[i]; - break; - }; - - default: - Assert(false, ExcNotImplemented()); - } - - Assert (index 1) - { - Assert (fe_data.dofs_per_line == degree-1, ExcInternalError()); - Assert ((fe_data.dofs_per_quad == (degree-1)*(degree-1)) || - (dim < 2), ExcInternalError()); - Assert ((fe_data.dofs_per_hex == (degree-1)*(degree-1)*(degree-1)) || - (dim < 3), ExcInternalError()); - - for (int i=0; i(GeometryInfo::lines_per_cell); ++i) - { - unsigned int index = fe_data.first_line_index + i*fe_data.dofs_per_line; - unsigned int incr = 0; - unsigned int tensorstart = 0; - // This again looks quite - // strange because of the odd - // numbering scheme. - switch (i+100*dim) - { - // lines in x-direction - case 100: - case 200: case 202: - case 300: case 302: case 304: case 306: - incr = 1; - break; - // lines in y-direction - case 201: case 203: - case 308: case 309: case 310: case 311: - incr = n; - break; - // lines in z-direction - case 301: case 303: case 305: case 307: - incr = n*n; - break; - default: - Assert(false, ExcNotImplemented()); - } - switch (i+100*dim) - { - // x=y=z=0 - case 100: - case 200: case 203: - case 300: case 303: case 308: - tensorstart = 0; - break; - // x=1 y=z=0 - case 201: - case 301: case 309: - tensorstart = degree; - break; - // y=1 x=z=0 - case 202: - case 304: case 307: - tensorstart = n*degree; - break; - // x=z=1 y=0 - case 310: - tensorstart = n*n*degree+degree; - break; - // z=1 x=y=0 - case 302: case 311: - tensorstart = n*n*degree; - break; - // x=y=1 z=0 - case 305: - tensorstart = n*degree+degree; - break; - // y=z=1 x=0 - case 306: - tensorstart = n*n*n-n; - break; - default: - Assert(false, ExcNotImplemented()); - } - - for (unsigned int jx = 1; jx(GeometryInfo::quads_per_cell); ++i) - { - unsigned int index = fe_data.first_quad_index+i*fe_data.dofs_per_quad; - unsigned int tensorstart = 0; - unsigned int incx = 0; - unsigned int incy = 0; - switch (i) - { - case 0: - tensorstart = 0; incx = 1; - if (dim==2) - incy = n; - else - incy = n*n; - break; - case 1: - tensorstart = n*degree; incx = 1; incy = n*n; - break; - case 2: - tensorstart = 0; incx = 1; incy = n; - break; - case 3: - tensorstart = degree; incx = n; incy = n*n; - break; - case 4: - tensorstart = n*n*degree; incx = 1; incy = n; - break; - case 5: - tensorstart = 0; incx = n; incy = n*n; - break; - default: - Assert(false, ExcNotImplemented()); - } - - for (unsigned int jy = 1; jy::hexes_per_cell > 0) - for (int i=0; i(GeometryInfo::hexes_per_cell); ++i) - { - unsigned int index = fe_data.first_hex_index; - - for (unsigned int jz = 1; jz std::vector FE_Q::face_lexicographic_to_hierarchic_numbering (const unsigned int degree) { - const FiniteElementData fe_data(FE_Q::get_dpo_vector(degree),1); - return FE_Q::lexicographic_to_hierarchic_numbering (fe_data, degree); + const FiniteElementData face_data(FE_Q::get_dpo_vector(degree),1); + std::vector face_renumber (face_data.dofs_per_cell); + FETools::lexicographic_to_hierarchic_numbering (face_data, face_renumber); + return face_renumber; }