From 4e76305d110a3c887952d05ebf5ea3f9a42c0b17 Mon Sep 17 00:00:00 2001 From: bangerth Date: Thu, 27 Jul 2006 00:40:02 +0000 Subject: [PATCH] First part of the implementation of hp vertex dof equivalences git-svn-id: https://svn.dealii.org/trunk@13440 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/source/dofs/hp_dof_handler.cc | 116 ++++++++++++++++-- 1 file changed, 109 insertions(+), 7 deletions(-) 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 07ad1952e8..df5bcc7b82 100644 --- a/deal.II/deal.II/source/dofs/hp_dof_handler.cc +++ b/deal.II/deal.II/source/dofs/hp_dof_handler.cc @@ -2012,16 +2012,118 @@ namespace hp std::vector user_flags; tria->save_user_flags(user_flags); const_cast &>(*tria).clear_user_flags (); - - unsigned int next_free_dof = 0; - active_cell_iterator cell = begin_active(), - endc = end(); - for (; cell != endc; ++cell) - next_free_dof = distribute_dofs_on_cell (cell, next_free_dof); - used_dofs = next_free_dof; + ///////////////////////////////// + + // Step 1: distribute DoFs on all + // active entities + { + unsigned int next_free_dof = 0; + active_cell_iterator cell = begin_active(), + endc = end(); + + for (; cell != endc; ++cell) + next_free_dof = distribute_dofs_on_cell (cell, next_free_dof); + + used_dofs = next_free_dof; + } + + ///////////////////////////////// + + // Step 2: identify certain dofs + // if the finite element tells us + // that they should have the same + // value + std::vector new_dof_indices (used_dofs); + + // Step 2a: do so for vertices + { + typedef + std::vector > Identities; + + Table<2,boost::shared_ptr > + vertex_dof_identities; + + // loop over all vertices and + // see which one we need to + // work on + for (unsigned int vertex_index=0; vertex_index 1) + { + const unsigned int + first_fe_index = nth_active_vertex_fe_index (vertex_index, 0); + + // loop over all the + // other FEs with which + // we want to identify + // the DoF indices of + // the first FE of + for (unsigned int f=1; f + (new Identities((*finite_elements)[first_fe_index] + .hp_vertex_dof_identities + ((*finite_elements)[other_fe_index]))); + + // then loop + // through the + // identities we + // have. first get + // the global + // numbers of the + // dofs we want to + // identify and + // make sure they + // are not yet + // constrained to + // anything else, + // except for to + // each other + Identities &identities + = *vertex_dof_identities[first_fe_index][other_fe_index]; + for (unsigned int i=0; i &>(*tria).load_user_flags(user_flags); } -- 2.39.5