From: Martin Kronbichler Date: Tue, 7 Apr 2020 07:17:25 +0000 (+0200) Subject: Deprecate old functions and use the new one X-Git-Tag: v9.2.0-rc1~238^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b15e787fe7ddddfc6db3d5dd11b51b8ac5e94d5a;p=dealii.git Deprecate old functions and use the new one --- diff --git a/include/deal.II/fe/fe_tools.h b/include/deal.II/fe/fe_tools.h index 860f58817c..1872d6f2ae 100644 --- a/include/deal.II/fe/fe_tools.h +++ b/include/deal.II/fe/fe_tools.h @@ -895,75 +895,72 @@ namespace FETools * in y-direction and finally in z-direction. Discontinuous elements of * class FE_DGQ() are numbered in this way, for example. * - * This function constructs a table which lexicographic index each degree of - * freedom in the hierarchic numbering would have. It operates on the - * continuous finite element given as first argument, and outputs the - * lexicographic indices in the second. - * - * Note that since this function uses specifics of the continuous finite - * elements, it can only operate on FiniteElementData objects inherent - * in FE_Q(). However, this function does not take a FE_Q object as it is - * also invoked by the FE_Q() constructor. - * - * It is assumed that the size of the output argument already matches the - * correct size, which is equal to the number of degrees of freedom in the - * finite element. + * This function returns a vector containing information about the + * lexicographic index each degree of freedom in the hierarchic numbering + * would have to a given degree of a continuous finite element. */ template - void - hierarchic_to_lexicographic_numbering(unsigned int degree, - std::vector &h2l); + std::vector + hierarchic_to_lexicographic_numbering(unsigned int degree); /** - * Like the previous function but instead of returning its result through - * the last argument return it as a value. + * Like the previous function but instead of returning its result as a value + * return it through the last argument. + * + * @deprecated Use the function that returns the renumbering in a vector + * instead. */ template - std::vector - hierarchic_to_lexicographic_numbering(unsigned int degree); + DEAL_II_DEPRECATED void + hierarchic_to_lexicographic_numbering(unsigned int degree, + std::vector &h2l); /** * Like the previous functions but using a FiniteElementData instead of the * polynomial degree. + * + * @deprecated Use the function that returns the renumbering in a vector and + * uses the degree of the basis as an argument instead. */ template - void + DEAL_II_DEPRECATED void hierarchic_to_lexicographic_numbering(const FiniteElementData &fe_data, std::vector & h2l); /** - * Like the previous function but instead of returning its result through - * the last argument return it as a value. + * @deprecated Use the function that uses the degree of the basis as an + * argument instead. */ template - std::vector - hierarchic_to_lexicographic_numbering(const FiniteElementData &fe_data); + DEAL_II_DEPRECATED std::vector + hierarchic_to_lexicographic_numbering(const FiniteElementData &fe_data); /** * This is the reverse function to the above one, generating the map from - * the lexicographic to the hierarchical numbering. All the remarks made - * about the above function are also valid here. + * the lexicographic to the hierarchical numbering for a given polynomial + * degree of a continuous finite element. All the remarks made about the + * above function are also valid here. */ template std::vector lexicographic_to_hierarchic_numbering(unsigned int degree); /** - * Like the previous functions but using a FiniteElementData instead of the - * polynomial degree. + * @deprecated Use the function that returns the renumbering in a vector and + * uses the degree of the basis as an argument instead. */ template - void + DEAL_II_DEPRECATED void lexicographic_to_hierarchic_numbering(const FiniteElementData &fe_data, std::vector & l2h); /** - * Like the previous function but instead of returning its result through - * the last argument return it as a value. + * @deprecated Use the function that uses the degree of the basis as an + * argument instead. */ template - std::vector - lexicographic_to_hierarchic_numbering(const FiniteElementData &fe_data); + DEAL_II_DEPRECATED std::vector + lexicographic_to_hierarchic_numbering(const FiniteElementData &fe_data); /** diff --git a/include/deal.II/fe/fe_tools.templates.h b/include/deal.II/fe/fe_tools.templates.h index b4a594fa3a..a20409b117 100644 --- a/include/deal.II/fe/fe_tools.templates.h +++ b/include/deal.II/fe/fe_tools.templates.h @@ -3123,17 +3123,15 @@ namespace FETools template - void - hierarchic_to_lexicographic_numbering(const unsigned int degree, - std::vector &h2l) + std::vector + hierarchic_to_lexicographic_numbering(const unsigned int degree) { // number of support points in each direction const unsigned int n = degree + 1; const unsigned int dofs_per_cell = Utilities::fixed_power(n); - // Assert size matches degree - AssertDimension(h2l.size(), dofs_per_cell); + std::vector h2l(dofs_per_cell); // polynomial degree const unsigned int dofs_per_line = degree - 1; @@ -3289,17 +3287,19 @@ namespace FETools default: Assert(false, ExcNotImplemented()); } + + return h2l; } template - std::vector - hierarchic_to_lexicographic_numbering(const unsigned int degree) + void + hierarchic_to_lexicographic_numbering(const unsigned int degree, + std::vector &h2l) { - std::vector h2l(Utilities::pow(degree + 1, dim)); - hierarchic_to_lexicographic_numbering(degree, h2l); - return h2l; + AssertDimension(h2l.size(), Utilities::fixed_power(degree + 1)); + h2l = hierarchic_to_lexicographic_numbering(degree); } @@ -3321,9 +3321,7 @@ namespace FETools hierarchic_to_lexicographic_numbering(const FiniteElementData &fe) { Assert(fe.n_components() == 1, ExcInvalidFE()); - std::vector h2l(fe.dofs_per_cell); - hierarchic_to_lexicographic_numbering(fe.dofs_per_line + 1, h2l); - return h2l; + return hierarchic_to_lexicographic_numbering(fe.dofs_per_line + 1); } diff --git a/include/deal.II/matrix_free/cuda_hanging_nodes_internal.h b/include/deal.II/matrix_free/cuda_hanging_nodes_internal.h index a4c87d51d1..605313d404 100644 --- a/include/deal.II/matrix_free/cuda_hanging_nodes_internal.h +++ b/include/deal.II/matrix_free/cuda_hanging_nodes_internal.h @@ -147,7 +147,7 @@ namespace CUDAWrappers fe_q.get_subface_interpolation_matrix(fe_q, 0, interpolation_matrix); std::vector mapping = - FETools::lexicographic_to_hierarchic_numbering<1>(FE_Q<1>(fe_degree)); + FETools::lexicographic_to_hierarchic_numbering<1>(fe_degree); FullMatrix mapped_matrix(fe_q.dofs_per_face, fe_q.dofs_per_face); @@ -282,8 +282,7 @@ namespace CUDAWrappers std::vector neighbor_dofs(dofs_per_face); const auto lex_face_mapping = - FETools::lexicographic_to_hierarchic_numbering( - FE_Q(fe_degree)); + FETools::lexicographic_to_hierarchic_numbering(fe_degree); for (const unsigned int face : GeometryInfo::face_indices()) { diff --git a/source/fe/fe_bernstein.cc b/source/fe/fe_bernstein.cc index 99294631a2..db2860bf1d 100644 --- a/source/fe/fe_bernstein.cc +++ b/source/fe/fe_bernstein.cc @@ -367,10 +367,7 @@ FE_Bernstein::renumber_bases(const unsigned int deg) { TensorProductPolynomials tpp( dealii::generate_complete_bernstein_basis(deg)); - std::vector renumber(Utilities::fixed_power(deg + 1)); - const FiniteElementData fe(this->get_dpo_vector(deg), 1, deg); - FETools::hierarchic_to_lexicographic_numbering(fe, renumber); - tpp.set_numbering(renumber); + tpp.set_numbering(FETools::hierarchic_to_lexicographic_numbering(deg)); return tpp; } diff --git a/source/fe/fe_trace.cc b/source/fe/fe_trace.cc index f844990475..d40fef3a92 100644 --- a/source/fe/fe_trace.cc +++ b/source/fe/fe_trace.cc @@ -50,9 +50,8 @@ FE_TraceQ::FE_TraceQ(const unsigned int degree) Assert(degree > 0, ExcMessage("FE_Trace can only be used for polynomial degrees " "greater than zero")); - std::vector renumber(this->dofs_per_face); - FETools::hierarchic_to_lexicographic_numbering(degree, renumber); - this->poly_space.set_numbering(renumber); + this->poly_space.set_numbering( + FETools::hierarchic_to_lexicographic_numbering(degree)); // Initialize face support points this->unit_face_support_points = fe_q.get_unit_face_support_points();