From 5e1d4effcf01ae5b4ac9e36946febd7d408787eb Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Thu, 6 Mar 2008 21:38:12 +0000 Subject: [PATCH] Interface to the minimum degree algorithm in boost to reorder matrices. git-svn-id: https://svn.dealii.org/trunk@15865 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/include/dofs/dof_renumbering.h | 36 +++++ .../deal.II/source/dofs/dof_renumbering.cc | 130 ++++++++++++++++++ 2 files changed, 166 insertions(+) diff --git a/deal.II/deal.II/include/dofs/dof_renumbering.h b/deal.II/deal.II/include/dofs/dof_renumbering.h index 8ba22cb130..b72a0344b8 100644 --- a/deal.II/deal.II/include/dofs/dof_renumbering.h +++ b/deal.II/deal.II/include/dofs/dof_renumbering.h @@ -298,6 +298,42 @@ class DoFRenumbering const DH&, const bool reversed_numbering = false, const bool use_constraints = false); + + /** + * Renumber the degrees of + * freedom based on the BOOST + * implementation of the + * minimum degree + * algorithm. This often + * results in slightly larger + * (by a few percent) + * bandwidths than the + * Cuthill-McKee algorithm, + * but sparse ILUs are often + * slightly (also by a few + * percent) better + * preconditioners. + */ + template + static void + minimum_degree (DH& dof_handler, + const bool reversed_numbering = false, + const bool use_constraints = false); + + /** + * Compute the renumbering + * for the minimum degree + * algorithm but do not + * actually renumber the + * degrees of freedom in the + * DoF handler argument. + */ + template + static void + compute_minimum_degree (std::vector& new_dof_indices, + const DH&, + const bool reversed_numbering = false, + const bool use_constraints = false); }; /** diff --git a/deal.II/deal.II/source/dofs/dof_renumbering.cc b/deal.II/deal.II/source/dofs/dof_renumbering.cc index 1d93fd6cc5..b0e1eebd95 100644 --- a/deal.II/deal.II/source/dofs/dof_renumbering.cc +++ b/deal.II/deal.II/source/dofs/dof_renumbering.cc @@ -39,6 +39,7 @@ #include #include #include +#include #include #include @@ -314,6 +315,126 @@ compute_king_ordering (std::vector& new_dof_indices, } + +template +void +DoFRenumbering::boost:: +minimum_degree (DH& dof_handler, + const bool reversed_numbering, + const bool use_constraints) +{ + std::vector renumbering(dof_handler.n_dofs(), + DH::invalid_dof_index); + compute_minimum_degree(renumbering, dof_handler, reversed_numbering, + use_constraints); + + // actually perform renumbering; + // this is dimension specific and + // thus needs an own function + dof_handler.renumber_dofs (renumbering); +} + + +template +void +DoFRenumbering::boost:: +compute_minimum_degree (std::vector& new_dof_indices, + const DH &dof_handler, + const bool reversed_numbering, + const bool use_constraints) +{ + Assert (use_constraints == false, ExcNotImplemented()); + + // the following code is pretty + // much a verbatim copy of the + // sample code for the + // minimum_degree_ordering manual + // page from the BOOST Graph + // Library + using namespace boost; + + int delta = 0; + + typedef double Type; + + // must be BGL directed graph now + typedef adjacency_list Graph; + typedef graph_traits::vertex_descriptor Vertex; + + int n = dof_handler.n_dofs(); + + Graph G(n); + + std::vector dofs_on_this_cell; + + typename DH::active_cell_iterator cell = dof_handler.begin_active(), + endc = dof_handler.end(); + + for (; cell!=endc; ++cell) + { + + const unsigned int dofs_per_cell = cell->get_fe().dofs_per_cell; + + dofs_on_this_cell.resize (dofs_per_cell); + + cell->get_dof_indices (dofs_on_this_cell); + for (unsigned int i=0; i dofs_on_this_cell[j]) + { + add_edge (dofs_on_this_cell[i], dofs_on_this_cell[j], G); + add_edge (dofs_on_this_cell[j], dofs_on_this_cell[i], G); + } + } + + + typedef std::vector Vector; + + + Vector inverse_perm(n, 0); + + Vector perm(n, 0); + + + Vector supernode_sizes(n, 1); + // init has to be 1 + + ::boost::property_map::type + id = get(vertex_index, G); + + + Vector degree(n, 0); + + + minimum_degree_ordering + (G, + make_iterator_property_map(°ree[0], id, degree[0]), + &inverse_perm[0], + &perm[0], + make_iterator_property_map(&supernode_sizes[0], id, supernode_sizes[0]), + delta, id); + + + for (int i=0; i &, const DoFHandler &, bool, bool); +template +void +DoFRenumbering::boost::minimum_degree (DoFHandler &, bool, bool); + +template +void +DoFRenumbering::boost::compute_minimum_degree (std::vector &, + const DoFHandler &, bool, bool); + // non-boost functions: -- 2.39.5