From a7f474034a912fc389abfe96deb4b1c606ad0025 Mon Sep 17 00:00:00 2001 From: Guido Kanschat Date: Fri, 18 Apr 2003 18:21:41 +0000 Subject: [PATCH] change of return parameter git-svn-id: https://svn.dealii.org/trunk@7408 0785d39b-7218-0410-832d-ea1e28bc413d --- .../include/numerics/dof_renumbering.h | 21 ++- .../deal.II/source/dofs/dof_renumbering.cc | 162 ++++++++++-------- 2 files changed, 102 insertions(+), 81 deletions(-) diff --git a/deal.II/deal.II/include/numerics/dof_renumbering.h b/deal.II/deal.II/include/numerics/dof_renumbering.h index 344bdaae36..f7ae3a3ef6 100644 --- a/deal.II/deal.II/include/numerics/dof_renumbering.h +++ b/deal.II/deal.II/include/numerics/dof_renumbering.h @@ -227,10 +227,11 @@ class DoFRenumbering * vector. */ template - static std::vector - compute_Cuthill_McKee (DoFHandler &dof_handler, - const bool reversed_numbering = false, - const bool use_constraints = false, + static void + compute_Cuthill_McKee (std::vector&, + const DoFHandler&, + const bool reversed_numbering = false, + const bool use_constraints = false, const std::vector &starting_indices = std::vector()); /** @@ -487,9 +488,10 @@ class DoFRenumbering * vector. */ template - static std::vector - compute_sort_selected_dofs_back (DoFHandler &dof_handler, - const std::vector &selected_dofs); + static void + compute_sort_selected_dofs_back (std::vector&, + const DoFHandler&, + const std::vector&selected_dofs); /** * Renumber the degrees of @@ -509,8 +511,9 @@ class DoFRenumbering * vector. */ template - static std::vector - compute_random (DoFHandler &dof_handler); + static void + compute_random (std::vector&, + const DoFHandler &dof_handler); /** * Exception diff --git a/deal.II/deal.II/source/dofs/dof_renumbering.cc b/deal.II/deal.II/source/dofs/dof_renumbering.cc index cd2393970e..40c373b3e7 100644 --- a/deal.II/deal.II/source/dofs/dof_renumbering.cc +++ b/deal.II/deal.II/source/dofs/dof_renumbering.cc @@ -53,14 +53,16 @@ extern "C" long int lrand48 (void); template void -DoFRenumbering::Cuthill_McKee (DoFHandler &dof_handler, - const bool reversed_numbering, - const bool use_constraints, - const std::vector &starting_indices) +DoFRenumbering::Cuthill_McKee ( + DoFHandler& dof_handler, + const bool reversed_numbering, + const bool use_constraints, + const std::vector& starting_indices) { - std::vector renumbering - =compute_Cuthill_McKee(dof_handler, reversed_numbering, - use_constraints, starting_indices); + std::vector renumbering(dof_handler.n_dofs(), + DoFHandler::invalid_dof_index); + compute_Cuthill_McKee(renumbering, dof_handler, reversed_numbering, + use_constraints, starting_indices); // actually perform renumbering; // this is dimension specific and @@ -71,11 +73,13 @@ DoFRenumbering::Cuthill_McKee (DoFHandler &dof_handler, template -std::vector -DoFRenumbering::compute_Cuthill_McKee (DoFHandler &dof_handler, - const bool reversed_numbering, - const bool use_constraints, - const std::vector &starting_indices) +void +DoFRenumbering::compute_Cuthill_McKee ( + std::vector& new_indices, + const DoFHandler& dof_handler, + const bool reversed_numbering, + const bool use_constraints, + const std::vector& starting_indices) { // make the connection graph SparsityPattern sparsity (dof_handler.n_dofs(), @@ -93,8 +97,8 @@ DoFRenumbering::compute_Cuthill_McKee (DoFHandler &dof_hand const unsigned int n_dofs = sparsity.n_rows(); // store the new dof numbers; invalid_dof_index means // that no new number was chosen yet - std::vector new_number(sparsity.n_rows(), - DoFHandler::invalid_dof_index); + Assert(new_indices.size() == n_dofs, + ExcDimensionMismatch(new_indices.size(), n_dofs)); // store the indices of the dofs renumbered // in the last round. Default to starting @@ -163,7 +167,7 @@ DoFRenumbering::compute_Cuthill_McKee (DoFHandler &dof_hand // enumerate the first round dofs for (unsigned int i=0; i!=last_round_dofs.size(); ++i) - new_number[last_round_dofs[i]] = next_free_number++; + new_indices[last_round_dofs[i]] = next_free_number++; bool all_dofs_renumbered = false; @@ -197,7 +201,7 @@ DoFRenumbering::compute_Cuthill_McKee (DoFHandler &dof_hand // eliminate dofs which are // already numbered for (int s=next_round_dofs.size()-1; s>=0; --s) - if (new_number[next_round_dofs[s]] != DoFHandler::invalid_dof_index) + if (new_indices[next_round_dofs[s]] != DoFHandler::invalid_dof_index) next_round_dofs.erase (next_round_dofs.begin() + s); // check whether there are any new @@ -237,7 +241,7 @@ DoFRenumbering::compute_Cuthill_McKee (DoFHandler &dof_hand // front: std::multimap::iterator i; for (i = dofs_by_coordination.begin(); i!=dofs_by_coordination.end(); ++i) - new_number[i->second] = next_free_number++; + new_indices[i->second] = next_free_number++; // after that: copy this round's // dofs for the next round @@ -261,29 +265,28 @@ DoFRenumbering::compute_Cuthill_McKee (DoFHandler &dof_hand // In any case, if not all DoFs // have been reached, renumbering // will not be possible - if (std::find (new_number.begin(), new_number.end(), DoFHandler::invalid_dof_index) + if (std::find (new_indices.begin(), new_indices.end(), DoFHandler::invalid_dof_index) != - new_number.end()) + new_indices.end()) Assert (false, ExcRenumberingIncomplete()); Assert (next_free_number == n_dofs, ExcRenumberingIncomplete()); #endif if (reversed_numbering) - for (std::vector::iterator i=new_number.begin(); - i!=new_number.end(); ++i) + for (std::vector::iterator i=new_indices.begin(); + i!=new_indices.end(); ++i) *i = n_dofs-*i-1; - - return new_number; } #ifdef ENABLE_MULTIGRID template -void DoFRenumbering::Cuthill_McKee (MGDoFHandler &dof_handler, - const unsigned int level, - const bool reversed_numbering, - const std::vector &starting_indices) +void DoFRenumbering::Cuthill_McKee ( + MGDoFHandler &dof_handler, + const unsigned int level, + const bool reversed_numbering, + const std::vector &starting_indices) { // make the connection graph SparsityPattern sparsity (dof_handler.n_dofs(level), @@ -293,7 +296,7 @@ void DoFRenumbering::Cuthill_McKee (MGDoFHandler &dof_handler const unsigned int n_dofs = sparsity.n_rows(); // store the new dof numbers; invalid_dof_index means // that no new number was chosen yet - std::vector new_number(n_dofs, DoFHandler::invalid_dof_index); + std::vector new_indices(n_dofs, DoFHandler::invalid_dof_index); // store the indices of the dofs renumbered // in the last round. Default to starting @@ -344,7 +347,7 @@ void DoFRenumbering::Cuthill_McKee (MGDoFHandler &dof_handler // enumerate the first round dofs for (unsigned int i=0; i!=last_round_dofs.size(); ++i) - new_number[last_round_dofs[i]] = next_free_number++; + new_indices[last_round_dofs[i]] = next_free_number++; bool all_dofs_renumbered = false; @@ -378,7 +381,7 @@ void DoFRenumbering::Cuthill_McKee (MGDoFHandler &dof_handler // eliminate dofs which are // already numbered for (int s=next_round_dofs.size()-1; s>=0; --s) - if (new_number[next_round_dofs[s]] != DoFHandler::invalid_dof_index) + if (new_indices[next_round_dofs[s]] != DoFHandler::invalid_dof_index) next_round_dofs.erase (next_round_dofs.begin() + s); // check whether there are any new @@ -416,7 +419,7 @@ void DoFRenumbering::Cuthill_McKee (MGDoFHandler &dof_handler //// std::multimap::iterator i; for (i = dofs_by_coordination.begin(); i!=dofs_by_coordination.end(); ++i) - new_number[i->second] = next_free_number++; + new_indices[i->second] = next_free_number++; // after that: copy this round's // dofs for the next round @@ -425,31 +428,32 @@ void DoFRenumbering::Cuthill_McKee (MGDoFHandler &dof_handler #ifdef DEBUG // test for all indices numbered - if (std::find (new_number.begin(), new_number.end(), + if (std::find (new_indices.begin(), new_indices.end(), DoFHandler::invalid_dof_index) != - new_number.end()) + new_indices.end()) Assert (false, ExcRenumberingIncomplete()); Assert (next_free_number == n_dofs, ExcRenumberingIncomplete()); #endif if (reversed_numbering) - for (std::vector::iterator i=new_number.begin(); i!=new_number.end(); ++i) + for (std::vector::iterator i=new_indices.begin(); i!=new_indices.end(); ++i) *i = n_dofs-*i; // actually perform renumbering; // this is dimension specific and // thus needs an own function - dof_handler.renumber_dofs (level, new_number); + dof_handler.renumber_dofs (level, new_indices); } #endif template void -DoFRenumbering::component_wise (DoFHandler &dof_handler, - const std::vector &component_order_arg) +DoFRenumbering::component_wise ( + DoFHandler &dof_handler, + const std::vector &component_order_arg) { std::vector renumbering (dof_handler.n_dofs(), DoFHandler::invalid_dof_index); @@ -464,9 +468,10 @@ DoFRenumbering::component_wise (DoFHandler &dof_handler, template void -DoFRenumbering::compute_component_wise (std::vector& new_indices, - const DoFHandler& dof_handler, - const std::vector &component_order_arg) +DoFRenumbering::compute_component_wise ( + std::vector& new_indices, + const DoFHandler& dof_handler, + 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(); @@ -618,9 +623,10 @@ DoFRenumbering::compute_component_wise (std::vector& new_indices, template -void DoFRenumbering::component_wise (MGDoFHandler& dof_handler, - unsigned int level, - const std::vector &component_order_arg) +void DoFRenumbering::component_wise ( + MGDoFHandler& dof_handler, + 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(); @@ -770,11 +776,13 @@ void DoFRenumbering::component_wise (MGDoFHandler& dof_handler, template void -DoFRenumbering::sort_selected_dofs_back (DoFHandler &dof_handler, - const std::vector &selected_dofs) +DoFRenumbering::sort_selected_dofs_back ( + DoFHandler &dof_handler, + const std::vector &selected_dofs) { - std::vector renumbering - =compute_sort_selected_dofs_back(dof_handler, selected_dofs); + std::vector renumbering(dof_handler.n_dofs(), + DoFHandler::invalid_dof_index); + compute_sort_selected_dofs_back(renumbering, dof_handler, selected_dofs); dof_handler.renumber_dofs(renumbering); } @@ -782,9 +790,11 @@ DoFRenumbering::sort_selected_dofs_back (DoFHandler &dof_handler, template -std::vector -DoFRenumbering::compute_sort_selected_dofs_back (DoFHandler &dof_handler, - const std::vector &selected_dofs) +void +DoFRenumbering::compute_sort_selected_dofs_back ( + std::vector& new_indices, + const DoFHandler& dof_handler, + const std::vector& selected_dofs) { const unsigned int n_dofs = dof_handler.n_dofs(); Assert (selected_dofs.size() == n_dofs, @@ -792,7 +802,9 @@ DoFRenumbering::compute_sort_selected_dofs_back (DoFHandler &dof_ha // re-sort the dofs according to // their selection state - std::vector new_dof_indices (n_dofs); + Assert (new_indices.size() == n_dofs, + ExcDimensionMismatch(new_indices.size(), n_dofs)); + const unsigned int n_selected_dofs = count (selected_dofs.begin(), selected_dofs.end(), false); @@ -802,18 +814,18 @@ DoFRenumbering::compute_sort_selected_dofs_back (DoFHandler &dof_ha for (unsigned int i=0; i -void DoFRenumbering::cell_wise_dg (MGDoFHandler& dof, - const unsigned int level, - const typename std::vector::cell_iterator>& cells) +void DoFRenumbering::cell_wise_dg ( + MGDoFHandler& dof, + const unsigned int level, + const typename std::vector::cell_iterator>& cells) { Assert(cells.size() == dof.get_tria().n_cells(level), ExcDimensionMismatch(cells.size(), @@ -977,7 +990,6 @@ struct CompCells } }; -//TODO:[RH or other] Remove necessity of copying a BIG vector renumbering template void @@ -995,7 +1007,7 @@ DoFRenumbering::downstream_dg (DoFHandler& dof, template void DoFRenumbering::compute_downstream_dg ( - std::vector& renumbering, + std::vector& new_indices, const DoFHandler& dof, const Point& direction) { @@ -1009,7 +1021,7 @@ DoFRenumbering::compute_downstream_dg ( copy (begin, end, ordered_cells.begin()); sort (ordered_cells.begin(), ordered_cells.end(), comparator); - compute_cell_wise_dg(renumbering, dof, ordered_cells); + compute_cell_wise_dg(new_indices, dof, ordered_cells); } @@ -1039,25 +1051,28 @@ template void DoFRenumbering::random (DoFHandler &dof_handler) { - std::vector renumbering - =compute_random(dof_handler); + std::vector renumbering(dof_handler.n_dofs(), + DoFHandler::invalid_dof_index); + compute_random(renumbering, dof_handler); dof_handler.renumber_dofs(renumbering); } template -std::vector -DoFRenumbering::compute_random (DoFHandler &dof_handler) +void +DoFRenumbering::compute_random ( + std::vector& new_indices, + const DoFHandler &dof_handler) { const unsigned int n_dofs = dof_handler.n_dofs(); + Assert(new_indices.size() == n_dofs, + ExcDimensionMismatch(new_indices.size(), n_dofs)); - std::vector new_indices (n_dofs); for (unsigned i=0; i const std::vector&); template -std::vector +void DoFRenumbering::compute_Cuthill_McKee -(DoFHandler&, +(std::vector&, + const DoFHandler&, const bool, const bool, const std::vector&); @@ -1130,9 +1146,10 @@ void DoFRenumbering::sort_selected_dofs_back const std::vector &); template -std::vector +void DoFRenumbering::compute_sort_selected_dofs_back -(DoFHandler &, +(std::vector&, + const DoFHandler &, const std::vector &); template @@ -1140,9 +1157,10 @@ void DoFRenumbering::random (DoFHandler &); template -std::vector +void DoFRenumbering::compute_random -(DoFHandler &); +(std::vector&, + const DoFHandler &); #ifdef ENABLE_MULTIGRID template -- 2.39.5