From: bangerth Date: Fri, 28 Jul 2006 17:17:39 +0000 (+0000) Subject: Implement unification of DoFs X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b045198a051140f8b712eb8e40cc93002d8ae7dd;p=dealii-svn.git Implement unification of DoFs git-svn-id: https://svn.dealii.org/trunk@13475 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/hp_dof_handler.h b/deal.II/deal.II/include/dofs/hp_dof_handler.h index 38f8d4aed0..d4b3315b62 100644 --- a/deal.II/deal.II/include/dofs/hp_dof_handler.h +++ b/deal.II/deal.II/include/dofs/hp_dof_handler.h @@ -17,6 +17,7 @@ #include #include +#include #include #include #include @@ -175,7 +176,7 @@ namespace hp virtual void clear (); /** - * Actually do the renumbering based on + * Renumber degrees of freedom based on * a list of new dof numbers for all the * dofs. * @@ -185,13 +186,25 @@ namespace hp * indices after renumbering in the * order of the old indices. * - * This function is called by the - * @p renumber_dofs function after computing - * the ordering of the degrees of freedom. - * However, you can call this function - * yourself, which is necessary if a user - * wants to implement an ordering scheme - * herself, for example downwind numbering. + * This function is called by + * the functions in + * DoFRenumbering function + * after computing the ordering + * of the degrees of freedom. + * However, you can call this + * function yourself, which is + * necessary if a user wants to + * implement an ordering scheme + * herself, for example + * downwind numbering. + * + * The @p new_number array must + * have a size equal to the + * number of degrees of + * freedom. Each entry must + * state the new global DoF + * number of the degree of + * freedom referenced. */ void renumber_dofs (const std::vector &new_numbers); @@ -1132,7 +1145,44 @@ namespace hp */ void compute_quad_dof_identities (std::vector &new_dof_indices) const; - + + /** + * Renumber the objects with + * the given and all lower + * structural dimensions, + * i.e. renumber vertices by + * giving a template argument + * of zero to the int2type + * argument, lines and vertices + * with one, etc. + * + * Note that in contrast to the + * public renumber_dofs() + * function, these internal + * functions do not ensure that + * the new DoFs are + * contiguously numbered. The + * function may therefore also + * be used to assign different + * DoFs the same number, for + * example to unify hp DoFs + * corresponding to different + * finite elements but + * co-located on the same + * entity. + */ + void renumber_dofs_internal (const std::vector &new_numbers, + internal::int2type<0>); + + void renumber_dofs_internal (const std::vector &new_numbers, + internal::int2type<1>); + + void renumber_dofs_internal (const std::vector &new_numbers, + internal::int2type<2>); + + void renumber_dofs_internal (const std::vector &new_numbers, + internal::int2type<3>); + /** * Space to store the DoF * numbers for the different @@ -1351,9 +1401,6 @@ namespace hp template <> unsigned int DoFHandler<3>::distribute_dofs_on_cell (active_cell_iterator &cell, unsigned int next_free_dof); - template <> void DoFHandler<1>::renumber_dofs (const std::vector &new_numbers); - template <> void DoFHandler<2>::renumber_dofs (const std::vector &new_numbers); - template <> void DoFHandler<3>::renumber_dofs (const std::vector &new_numbers); template <> void DoFHandler<1>::reserve_space (); template <> void DoFHandler<2>::reserve_space (); template <> void DoFHandler<3>::reserve_space (); diff --git a/deal.II/deal.II/source/dofs/hp_dof_handler.cc b/deal.II/deal.II/source/dofs/hp_dof_handler.cc index 51005fb7e5..61c71f5b92 100644 --- a/deal.II/deal.II/source/dofs/hp_dof_handler.cc +++ b/deal.II/deal.II/source/dofs/hp_dof_handler.cc @@ -2444,22 +2444,54 @@ namespace hp // faces and other // lower-dimensional objects // where elements come together - std::vector new_dof_indices (used_dofs, - deal_II_numbers::invalid_unsigned_int); - compute_vertex_dof_identities (new_dof_indices); - compute_line_dof_identities (new_dof_indices); - compute_quad_dof_identities (new_dof_indices); + std::vector + constrained_indices (used_dofs, deal_II_numbers::invalid_unsigned_int); + compute_vertex_dof_identities (constrained_indices); + compute_line_dof_identities (constrained_indices); + compute_quad_dof_identities (constrained_indices); + + const unsigned int used_dofs_before = used_dofs; + + // loop over all dofs and assign + // new numbers to those which are + // not constrained + std::vector + new_dof_indices (used_dofs, deal_II_numbers::invalid_unsigned_int); + unsigned int next_free_dof = 0; + for (unsigned int i=0; i()); + + used_dofs = next_free_dof; // finally restore the user flags @@ -2654,10 +2686,8 @@ namespace hp #endif -#if deal_II_dimension == 1 - - template <> - void DoFHandler<1>::renumber_dofs (const std::vector &new_numbers) + template + void DoFHandler::renumber_dofs (const std::vector &new_numbers) { Assert (new_numbers.size() == n_dofs(), ExcRenumberingIncomplete()); #ifdef DEBUG @@ -2671,128 +2701,135 @@ namespace hp unsigned int i = 0; for (; p!=tmp.end(); ++p, ++i) Assert (*p == i, ExcNewNumbersNotConsecutive(i)); - }; + } #endif - // note that we can not use cell iterators - // in this function since then we would - // renumber the dofs on the interface of - // two cells more than once. Anyway, this - // way it's not only more correct but also - // faster; note, however, that dof numbers - // may be invalid_dof_index, namely when the appropriate - // vertex/line/etc is unused - for (std::vector::iterator i=vertex_dofs.begin(); i!=vertex_dofs.end(); ++i) - if (*i != invalid_dof_index) - *i = new_numbers[*i]; - - for (unsigned int level=0; level::iterator i=levels[level]->lines.dofs.begin(); - i!=levels[level]->lines.dofs.end(); ++i) - if (*i != invalid_dof_index) - *i = new_numbers[*i]; + renumber_dofs_internal (new_numbers, internal::int2type()); } -#endif - -#if deal_II_dimension == 2 - template <> - void DoFHandler<2>::renumber_dofs (const std::vector &new_numbers) + template + void + DoFHandler:: + renumber_dofs_internal (const std::vector &new_numbers, + internal::int2type<0>) { Assert (new_numbers.size() == n_dofs(), ExcRenumberingIncomplete()); -#ifdef DEBUG - // assert that the new indices are - // consecutively numbered - if (true) - { - std::vector tmp(new_numbers); - std::sort (tmp.begin(), tmp.end()); - std::vector::const_iterator p = tmp.begin(); - unsigned int i = 0; - for (; p!=tmp.end(); ++p, ++i) - Assert (*p == i, ExcNewNumbersNotConsecutive(i)); - }; -#endif - // note that we can not use cell iterators - // in this function since then we would - // renumber the dofs on the interface of - // two cells more than once. Anyway, this - // way it's not only more correct but also - // faster; note, however, that dof numbers - // may be invalid_dof_index, namely when the appropriate - // vertex/line/etc is unused - for (std::vector::iterator i=vertex_dofs.begin(); i!=vertex_dofs.end(); ++i) - if (*i != invalid_dof_index) - *i = new_numbers[*i]; - - for (unsigned int level=0; level::iterator i=levels[level]->quads.dofs.begin(); - i!=levels[level]->quads.dofs.end(); ++i) - if (*i != invalid_dof_index) - *i = new_numbers[*i]; - } - for (std::vector::iterator i=faces->lines.dofs.begin(); - i!=faces->lines.dofs.end(); ++i) - if (*i != invalid_dof_index) - *i = new_numbers[*i]; + const unsigned int n_active_fe_indices + = n_active_vertex_fe_indices (vertex_index); + for (unsigned int f=0; f + void + DoFHandler:: + renumber_dofs_internal (const std::vector &new_numbers, + internal::int2type<1>) + { + Assert (new_numbers.size() == n_dofs(), ExcRenumberingIncomplete()); - template <> - void DoFHandler<3>::renumber_dofs (const std::vector &new_numbers) + renumber_dofs_internal (new_numbers, internal::int2type<0>()); + + for (line_iterator line=begin_line(); line!=end_line(); ++line) + { + const unsigned int n_active_fe_indices + = line->n_active_fe_indices (); + + for (unsigned int f=0; fnth_active_fe_index (f); + + for (unsigned int d=0; d<(*finite_elements)[fe_index].dofs_per_line; ++d) + line->set_dof_index (d, + new_numbers[line->dof_index(d,fe_index)], + fe_index); + } + } + } + + +#if deal_II_dimension >= 2 + + template + void + DoFHandler:: + renumber_dofs_internal (const std::vector &new_numbers, + internal::int2type<2>) { Assert (new_numbers.size() == n_dofs(), ExcRenumberingIncomplete()); -#ifdef DEBUG - // assert that the new indices are - // consecutively numbered - if (true) + + renumber_dofs_internal (new_numbers, internal::int2type<1>()); + + for (quad_iterator quad=begin_quad(); quad!=end_quad(); ++quad) { - std::vector tmp(new_numbers); - std::sort (tmp.begin(), tmp.end()); - std::vector::const_iterator p = tmp.begin(); - unsigned int i = 0; - for (; p!=tmp.end(); ++p, ++i) - Assert (*p == i, ExcNewNumbersNotConsecutive(i)); - }; + const unsigned int n_active_fe_indices + = quad->n_active_fe_indices (); + + for (unsigned int f=0; fnth_active_fe_index (f); + + for (unsigned int d=0; d<(*finite_elements)[fe_index].dofs_per_quad; ++d) + quad->set_dof_index (d, + new_numbers[quad->dof_index(d,fe_index)], + fe_index); + } + } + } + #endif - // note that we can not use cell iterators - // in this function since then we would - // renumber the dofs on the interface of - // two cells more than once. Anyway, this - // way it's not only more correct but also - // faster; note, however, that dof numbers - // may be invalid_dof_index, namely when the appropriate - // vertex/line/etc is unused - for (std::vector::iterator i=vertex_dofs.begin(); i!=vertex_dofs.end(); ++i) - if (*i != invalid_dof_index) - *i = new_numbers[*i]; - - for (unsigned int level=0; level= 3 + + template + void + DoFHandler:: + renumber_dofs_internal (const std::vector &new_numbers, + internal::int2type<3>) + { + Assert (new_numbers.size() == n_dofs(), ExcRenumberingIncomplete()); + + renumber_dofs_internal (new_numbers, internal::int2type<2>()); + + for (hex_iterator hex=begin_hex(); hex!=end_hex(); ++hex) { - for (std::vector::iterator i=levels[level]->hexes.dofs.begin(); - i!=levels[level]->hexes.dofs.end(); ++i) - if (*i != invalid_dof_index) - *i = new_numbers[*i]; + const unsigned int n_active_fe_indices + = hex->n_active_fe_indices (); + + for (unsigned int f=0; fnth_active_fe_index (f); + + for (unsigned int d=0; d<(*finite_elements)[fe_index].dofs_per_hex; ++d) + hex->set_dof_index (d, + new_numbers[hex->dof_index(d,fe_index)], + fe_index); + } } - - for (std::vector::iterator i=faces->lines.dofs.begin(); - i!=faces->lines.dofs.end(); ++i) - if (*i != invalid_dof_index) - *i = new_numbers[*i]; - for (std::vector::iterator i=faces->quads.dofs.begin(); - i!=faces->quads.dofs.end(); ++i) - if (*i != invalid_dof_index) - *i = new_numbers[*i]; } #endif