From 5b22356f8dd2e58582a58b24cae5610c223a3a5b Mon Sep 17 00:00:00 2001 From: Guido Kanschat Date: Thu, 3 Jul 2003 17:02:08 +0000 Subject: [PATCH] component_wise unified git-svn-id: https://svn.dealii.org/trunk@7839 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/include/dofs/dof_renumbering.h | 2 +- .../deal.II/source/dofs/dof_renumbering.cc | 177 +++--------------- 2 files changed, 32 insertions(+), 147 deletions(-) diff --git a/deal.II/deal.II/include/dofs/dof_renumbering.h b/deal.II/deal.II/include/dofs/dof_renumbering.h index ea31a9d14c..da9735c5cb 100644 --- a/deal.II/deal.II/include/dofs/dof_renumbering.h +++ b/deal.II/deal.II/include/dofs/dof_renumbering.h @@ -330,7 +330,7 @@ class DoFRenumbering * vector. */ template - static void + static unsigned int compute_component_wise (std::vector& new_index, ITERATOR& start, const ENDITERATOR& end, diff --git a/deal.II/deal.II/source/dofs/dof_renumbering.cc b/deal.II/deal.II/source/dofs/dof_renumbering.cc index 3c72818669..99c3de6323 100644 --- a/deal.II/deal.II/source/dofs/dof_renumbering.cc +++ b/deal.II/deal.II/source/dofs/dof_renumbering.cc @@ -59,6 +59,7 @@ template class WrapMGDoFIterator : public T { public: + WrapMGDoFIterator (const T& t) : T(t) {} void get_dof_indices (std::vector& v) const { T::get_mg_dof_indices(v); @@ -479,21 +480,26 @@ DoFRenumbering::component_wise ( start = dof_handler.begin_active(); const typename DoFHandler::cell_iterator end = dof_handler.end(); - + + unsigned int result = compute_component_wise::active_cell_iterator, typename DoFHandler::cell_iterator>(renumbering, start, end, component_order_arg); - if (renumbering.size()!=0) - dof_handler.renumber_dofs (renumbering); + if (result == 0) return; + + Assert (result == dof_handler.n_dofs(), + ExcRenumberingIncomplete()); + + dof_handler.renumber_dofs (renumbering); } template -void +unsigned int DoFRenumbering::compute_component_wise ( std::vector& new_indices, ITERATOR& start, @@ -513,7 +519,7 @@ DoFRenumbering::compute_component_wise ( if (fe.n_components() == 1) { new_indices.resize(0); - return; + return 0; } // Copy last argument into a @@ -648,8 +654,7 @@ DoFRenumbering::compute_component_wise ( new_indices[*dof_index] = next_free_index++; }; -// Assert (next_free_index == dof_handler.n_dofs(), -// ExcInternalError()); + return next_free_index; } @@ -660,148 +665,28 @@ void DoFRenumbering::component_wise ( unsigned int level, const std::vector &component_order_arg) { - const unsigned int dofs_per_cell = dof_handler.get_fe().dofs_per_cell; - const FiniteElement &fe = dof_handler.get_fe(); + std::vector renumbering (dof_handler.n_dofs(level), + DoFHandler::invalid_dof_index); - // do nothing if the FE has only - // one component - if (fe.n_components() == 1) - return; + typedef + WrapMGDoFIterator::cell_iterator> ITERATOR; - std::vector component_order (component_order_arg); - if (component_order.size() == 0) - for (unsigned int i=0; i local_dof_indices(dofs_per_cell); - - // prebuilt list to which component - // a given dof on a cell - // belongs. note that we get into - // trouble here if the shape - // function is not primitive, since - // then there is no single vector - // component to which it - // belongs. in this case, assign it - // to the first vector component to - // which it belongs - std::vector component_list (dofs_per_cell); - for (unsigned int i=0; i > component_to_dof_map (fe.n_components()); - for (typename MGDoFHandler::cell_iterator - cell=dof_handler.begin(level); - cell!=dof_handler.end(level); ++cell) - { - // on each cell: get dof indices - // and insert them into the global - // list using their component - cell->get_mg_dof_indices (local_dof_indices); - for (unsigned int i=0; i::cell_iterator + istart =dof_handler.begin(level); + ITERATOR start = istart; + typename MGDoFHandler::cell_iterator + iend = dof_handler.end(level); + const ITERATOR end = iend; + + unsigned int result = + compute_component_wise( + renumbering, start, end, component_order_arg); + + Assert (result == dof_handler.n_dofs(level), + ExcRenumberingIncomplete()); - // now we've got all indices sorted - // into buckets labelled with their - // component number. we've only got - // to traverse this list and assign - // the new indices - // - // however, we first want to sort - // the indices entered into the - // buckets to preserve the order - // within each component and during - // this also remove duplicate - // entries - // - // note that we no - // longer have to care about - // non-primitive shape functions - // since the buckets corresponding - // to the second and following - // vector components of a - // non-primitive FE will simply be - // empty, everything being shoved - // into the first one - for (unsigned int component=0; component new_indices (dof_handler.n_dofs(level), - DoFHandler::invalid_dof_index); - for (unsigned int c=0; c::const_iterator - begin_of_component = component_to_dof_map[component].begin(), - end_of_component = component_to_dof_map[component].end(); - - for (typename std::vector::const_iterator - dof_index = begin_of_component; - dof_index != end_of_component; ++dof_index) - new_indices[*dof_index] = next_free_index++; - }; - - Assert (next_free_index == dof_handler.n_dofs(level), - ExcInternalError()); - - dof_handler.renumber_dofs (level, new_indices); + if (renumbering.size()!=0) + dof_handler.renumber_dofs (renumbering); } -- 2.39.5