From a1fa57cc9f5857bf4aa7b4e8346e4d3f3c08fe1f Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Wed, 20 Sep 2023 15:47:46 +0200 Subject: [PATCH] Improve performance of Cuthill-McKee reordering --- source/lac/sparsity_tools.cc | 72 +++++++++++++++++------------------- 1 file changed, 33 insertions(+), 39 deletions(-) diff --git a/source/lac/sparsity_tools.cc b/source/lac/sparsity_tools.cc index b11d573c08..65fdb2c792 100644 --- a/source/lac/sparsity_tools.cc +++ b/source/lac/sparsity_tools.cc @@ -624,31 +624,38 @@ namespace SparsityTools for (const auto &last_round_dof : last_round_dofs) new_indices[last_round_dof] = next_free_number++; + // store the indices of the dofs to be renumbered in the next round + std::vector next_round_dofs; + + // store for each coordination number the dofs with these coordination + // number + std::vector> + dofs_by_coordination; + // now do as many steps as needed to renumber all dofs while (true) { - // store the indices of the dofs to be renumbered in the next round - std::vector next_round_dofs; + next_round_dofs.clear(); // find all neighbors of the dofs numbered in the last round for (const auto dof : last_round_dofs) - for (DynamicSparsityPattern::iterator j = sparsity.begin(dof); - j < sparsity.end(dof); - ++j) - next_round_dofs.push_back(j->column()); - - // sort dof numbers - std::sort(next_round_dofs.begin(), next_round_dofs.end()); - - // delete multiple entries - next_round_dofs.erase(std::unique(next_round_dofs.begin(), - next_round_dofs.end()), - next_round_dofs.end()); + { + const unsigned int row_length = sparsity.row_length(dof); + for (unsigned int i = 0; i < row_length; ++i) + { + // skip dofs which are already numbered + const auto column = sparsity.column_number(dof, i); + if (new_indices[column] == numbers::invalid_size_type) + { + next_round_dofs.push_back(column); - // eliminate dofs which are already numbered - for (int s = next_round_dofs.size() - 1; s >= 0; --s) - if (new_indices[next_round_dofs[s]] != numbers::invalid_size_type) - next_round_dofs.erase(next_round_dofs.begin() + s); + // assign a dummy value to 'new_indices' to avoid adding + // the same index again; those will get the right number + // at the end of the outer 'while' loop + new_indices[column] = 0; + } + } + } // check whether there are any new dofs in the list. if there are // none, then we have completely numbered the current component of the @@ -680,32 +687,19 @@ namespace SparsityTools } - - // store for each coordination number the dofs with these coordination - // number - std::multimap - dofs_by_coordination; - // find coordination number for each of these dofs + dofs_by_coordination.clear(); for (const types::global_dof_index next_round_dof : next_round_dofs) - { - const DynamicSparsityPattern::size_type coordination = - sparsity.row_length(next_round_dof); - - // insert this dof at its coordination number - const std::pair - new_entry(coordination, next_round_dof); - dofs_by_coordination.insert(new_entry); - } + dofs_by_coordination.emplace_back(sparsity.row_length(next_round_dof), + next_round_dof); + std::sort(dofs_by_coordination.begin(), dofs_by_coordination.end()); // assign new DoF numbers to the elements of the present front: - std::multimap::iterator i; - for (i = dofs_by_coordination.begin(); i != dofs_by_coordination.end(); - ++i) - new_indices[i->second] = next_free_number++; + for (const auto &i : dofs_by_coordination) + new_indices[i.second] = next_free_number++; - // after that: copy this round's dofs for the next round - last_round_dofs = next_round_dofs; + // after that: use this round's dofs for the next round + last_round_dofs.swap(next_round_dofs); } // test for all indices numbered. this mostly tests whether the -- 2.39.5