]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Write wrappers for the BOOST implementation of the King reordering algorithm.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 5 Mar 2008 15:59:46 +0000 (15:59 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 5 Mar 2008 15:59:46 +0000 (15:59 +0000)
git-svn-id: https://svn.dealii.org/trunk@15861 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_renumbering.h
deal.II/deal.II/source/dofs/dof_renumbering.cc

index 14ef002a63444e40b28e708ab286163f2de3bfaa..8ba22cb130c069d16e5dcaad337ddd383f4f8b8d 100644 (file)
@@ -221,7 +221,7 @@ class DoFRenumbering
                                      * the Boost Graph Library (BGL)
                                      * by Jeremy Siek and others.
                                      *
-                                     * While often slight slower to
+                                     * While often slighty slower to
                                      * compute, the algorithms using
                                      * often lead to matrices with
                                      * smaller bandwidths and sparse
@@ -264,6 +264,40 @@ 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 King
+                                         * 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 <class DH>
+       static void 
+       king_ordering (DH&                              dof_handler,
+                      const bool                       reversed_numbering = false,
+                      const bool                       use_constraints    = false);
+
+                                        /**
+                                         * Compute the renumbering
+                                         * for the King algorithm but
+                                         * do not actually renumber
+                                         * the degrees of freedom in
+                                         * the DoF handler argument.
+                                         */
+       template <class DH>
+       static void
+       compute_king_ordering (std::vector<unsigned int>& new_dof_indices,
+                              const DH&,
+                              const bool reversed_numbering = false,
+                              const bool use_constraints    = false);
     };
     
                                     /**
index 36825e3645c2531fdc99585998c88916013c64d3..1d93fd6cc5678302b2868a4a7e2a0b73636d5410 100644 (file)
@@ -38,7 +38,7 @@
 #include <boost/config.hpp>
 #include <boost/graph/adjacency_list.hpp>
 #include <boost/graph/cuthill_mckee_ordering.hpp>
-#include "boost/graph/minimum_degree_ordering.hpp"
+#include <boost/graph/king_ordering.hpp>
 #include <boost/graph/properties.hpp>
 #include <boost/graph/bandwidth.hpp>
 
@@ -259,6 +259,61 @@ compute_Cuthill_McKee (std::vector<unsigned int>& new_dof_indices,
 }
 
 
+
+template <class DH>
+void
+DoFRenumbering::boost::
+king_ordering (DH&              dof_handler,
+              const bool       reversed_numbering,
+              const bool       use_constraints)
+{
+  std::vector<unsigned int> renumbering(dof_handler.n_dofs(),
+                                       DH::invalid_dof_index);
+  compute_king_ordering(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 <class DH>
+void
+DoFRenumbering::boost::
+compute_king_ordering (std::vector<unsigned int>& new_dof_indices,
+                      const DH        &dof_handler,
+                      const bool       reversed_numbering,
+                      const bool       use_constraints)
+{
+  types::Graph
+    graph(dof_handler.n_dofs());
+  types::property_map<types::Graph,types::vertex_degree_t>::type
+    graph_degree;
+
+  create_graph (dof_handler, use_constraints, graph, graph_degree);
+    
+  types::property_map<types::Graph, types::vertex_index_t>::type
+    index_map = get(::boost::vertex_index, graph);
+
+
+  std::vector<types::Vertex> inv_perm(num_vertices(graph));
+
+  if (reversed_numbering == false)
+    ::boost::king_ordering(graph, inv_perm.rbegin());
+  else
+    ::boost::king_ordering(graph, inv_perm.begin());
+    
+  for (types::size_type c = 0; c != inv_perm.size(); ++c)
+    new_dof_indices[index_map[inv_perm[c]]] = c;
+
+  Assert (std::find (new_dof_indices.begin(), new_dof_indices.end(),
+                    DH::invalid_dof_index) == new_dof_indices.end(),
+         ExcInternalError());
+}
+
+
   
 
 
@@ -1650,6 +1705,17 @@ void
 DoFRenumbering::boost::compute_Cuthill_McKee (std::vector<unsigned int> &,
                                              const DoFHandler<deal_II_dimension> &, bool, bool);
 
+template
+void
+DoFRenumbering::boost::king_ordering (DoFHandler<deal_II_dimension> &, bool, bool);
+
+template
+void
+DoFRenumbering::boost::compute_king_ordering (std::vector<unsigned int> &,
+                                             const DoFHandler<deal_II_dimension> &, bool, bool);
+
+
+// non-boost functions:
 
 template
 void DoFRenumbering::Cuthill_McKee<DoFHandler<deal_II_dimension> >

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.