]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Change arguments for SparsityTools::reorder_Cuthill_McKee and GridTools::get_face_con...
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Fri, 10 Apr 2015 07:49:58 +0000 (09:49 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Sat, 11 Apr 2015 06:22:05 +0000 (08:22 +0200)
include/deal.II/grid/grid_tools.h
include/deal.II/lac/sparsity_tools.h
source/distributed/tria.cc
source/dofs/dof_renumbering.cc
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
source/lac/sparsity_tools.cc
tests/bits/gerold_2.cc
tests/lac/sparsity_tools_01.cc

index e8429c9ea7f14474a2ae053493bed8c781abefd2..7483bb56f6e1c82bbd378a688cfd3d26a4ce334b 100644 (file)
@@ -573,7 +573,17 @@ namespace GridTools
   template <int dim, int spacedim>
   void
   get_face_connectivity_of_cells (const Triangulation<dim, spacedim> &triangulation,
-                                  SparsityPattern                    &connectivity);
+                                  DynamicSparsityPattern             &connectivity);
+
+  /**
+   * As above, but filling a SparsityPattern object instead.
+   *
+   * @deprecated
+   */
+  template <int dim, int spacedim>
+  void
+  get_face_connectivity_of_cells (const Triangulation<dim, spacedim> &triangulation,
+                                  SparsityPattern                    &connectivity) DEAL_II_DEPRECATED;
 
   /**
    * Use the METIS partitioner to generate a partitioning of the active cells
index 09c05a4dde086d7b2829cb607a827b53c1cd667a..1c5fb872485591ea9091f571caff4cd6c5a7c3aa 100644 (file)
@@ -20,6 +20,7 @@
 #include <deal.II/base/config.h>
 #include <deal.II/base/exceptions.h>
 #include <deal.II/lac/sparsity_pattern.h>
+#include <deal.II/lac/dynamic_sparsity_pattern.h>
 
 #include <vector>
 
@@ -129,9 +130,19 @@ namespace SparsityTools
    * part of the algorithm that chooses starting indices.
    */
   void
-  reorder_Cuthill_McKee (const SparsityPattern     &sparsity,
+  reorder_Cuthill_McKee (const DynamicSparsityPattern &sparsity,
+                         std::vector<DynamicSparsityPattern::size_type> &new_indices,
+                         const std::vector<DynamicSparsityPattern::size_type> &starting_indices = std::vector<DynamicSparsityPattern::size_type>());
+
+  /**
+   * As above, but taking a SparsityPattern object instead.
+   *
+   * @deprecated
+   */
+  void
+  reorder_Cuthill_McKee (const SparsityPattern &sparsity,
                          std::vector<SparsityPattern::size_type> &new_indices,
-                         const std::vector<SparsityPattern::size_type> &starting_indices = std::vector<SparsityPattern::size_type>());
+                         const std::vector<SparsityPattern::size_type> &starting_indices = std::vector<SparsityPattern::size_type>()) DEAL_II_DEPRECATED;
 
 
 #ifdef DEAL_II_WITH_MPI
index f56143f2f25535837a4474d432ebca0004156a24..7eb7533fe5a308bee5f0480c3252c11e728fadcb 100644 (file)
@@ -18,7 +18,7 @@
 #include <deal.II/base/memory_consumption.h>
 #include <deal.II/base/logstream.h>
 #include <deal.II/lac/sparsity_tools.h>
-#include <deal.II/lac/sparsity_pattern.h>
+#include <deal.II/lac/dynamic_sparsity_pattern.h>
 #include <deal.II/grid/tria.h>
 #include <deal.II/grid/tria_accessor.h>
 #include <deal.II/grid/tria_iterator.h>
@@ -2489,7 +2489,7 @@ namespace parallel
     void
     Triangulation<dim,spacedim>::setup_coarse_cell_to_p4est_tree_permutation ()
     {
-      SparsityPattern cell_connectivity;
+      DynamicSparsityPattern cell_connectivity;
       GridTools::get_face_connectivity_of_cells (*this, cell_connectivity);
       coarse_cell_to_p4est_tree_permutation.resize (this->n_cells(0));
       SparsityTools::
index f45484c5f11ba31544f89b19bf697c1b4cc857fa..e209c33cffaf12d2c36cfc29601b32460a4db7b6 100644 (file)
@@ -398,63 +398,48 @@ namespace DoFRenumbering
     constraints.close ();
 
     IndexSet locally_owned = dof_handler.locally_owned_dofs();
-    SparsityPattern sparsity;
 
     // otherwise compute the Cuthill-McKee permutation
-    if (DH::dimension == 1)
-      {
-        sparsity.reinit (dof_handler.n_dofs(),
-                         dof_handler.n_dofs(),
-                         dof_handler.max_couplings_between_dofs());
-        DoFTools::make_sparsity_pattern (dof_handler, sparsity, constraints);
-        sparsity.compress();
-      }
-    else
-      {
-        DynamicSparsityPattern dsp (dof_handler.n_dofs(),
-                                    dof_handler.n_dofs(),
-                                    dof_handler.locally_owned_dofs());
-        DoFTools::make_sparsity_pattern (dof_handler, dsp, constraints);
-
-        // If the index set is not complete, need to get indices in local
-        // index space.
-        if (dof_handler.locally_owned_dofs().n_elements() !=
-            dof_handler.locally_owned_dofs().size())
-          {
-            // Create sparsity pattern from dsp by transferring its indices to
-            // processor-local index space and doing Cuthill-McKee there
-            std::vector<unsigned int> row_lengths(locally_owned.n_elements());
-            for (unsigned int i=0; i<locally_owned.n_elements(); ++i)
-              row_lengths[i] = dsp.row_length(locally_owned.nth_index_in_set(i));
-            sparsity.reinit(locally_owned.n_elements(), locally_owned.n_elements(),
-                            row_lengths);
-            std::vector<types::global_dof_index> row_entries;
-            for (unsigned int i=0; i<locally_owned.n_elements(); ++i)
-              {
-                const types::global_dof_index row = locally_owned.nth_index_in_set(i);
-                row_entries.resize(0);
-                for (DynamicSparsityPattern::row_iterator it =
-                       dsp.row_begin(row); it != dsp.row_end(row); ++it)
-                  if (*it != row && locally_owned.is_element(*it))
-                    row_entries.push_back(locally_owned.index_within_set(*it));
-                sparsity.add_entries(i, row_entries.begin(), row_entries.end(),
-                                     true);
-              }
-            sparsity.compress();
-          }
-        else
-          sparsity.copy_from(dsp);
-      }
+    DynamicSparsityPattern dsp (dof_handler.n_dofs(),
+                                dof_handler.n_dofs(),
+                                dof_handler.locally_owned_dofs());
+    DoFTools::make_sparsity_pattern (dof_handler, dsp, constraints);
 
     // constraints are not needed anymore
     constraints.clear ();
 
-    Assert(new_indices.size() == sparsity.n_rows(),
-           ExcDimensionMismatch(new_indices.size(),
-                                sparsity.n_rows()));
+    // If the index set is not complete, need to get indices in local index
+    // space.
+    if (dof_handler.locally_owned_dofs().n_elements() !=
+        dof_handler.locally_owned_dofs().size())
+      {
+        // Create sparsity pattern from dsp by transferring its indices to
+        // processor-local index space and doing Cuthill-McKee there
+        DynamicSparsityPattern sparsity(locally_owned.n_elements(),
+                                        locally_owned.n_elements());
+        std::vector<types::global_dof_index> row_entries;
+        for (unsigned int i=0; i<locally_owned.n_elements(); ++i)
+          {
+            const types::global_dof_index row = locally_owned.nth_index_in_set(i);
+            row_entries.resize(0);
+            for (DynamicSparsityPattern::row_iterator it =
+                   dsp.row_begin(row); it != dsp.row_end(row); ++it)
+              if (*it != row && locally_owned.is_element(*it))
+                row_entries.push_back(locally_owned.index_within_set(*it));
+            sparsity.add_entries(i, row_entries.begin(), row_entries.end(),
+                                 true);
+          }
 
-    SparsityTools::reorder_Cuthill_McKee (sparsity, new_indices,
-                                          starting_indices);
+        AssertDimension(new_indices.size(), sparsity.n_rows());
+        SparsityTools::reorder_Cuthill_McKee (sparsity, new_indices,
+                                              starting_indices);
+      }
+    else
+      {
+        AssertDimension(new_indices.size(), dsp.n_rows());
+        SparsityTools::reorder_Cuthill_McKee (dsp, new_indices,
+                                              starting_indices);
+      }
 
     if (reversed_numbering)
       new_indices = Utilities::reverse_permutation (new_indices);
@@ -476,25 +461,12 @@ namespace DoFRenumbering
            ExcNotInitialized());
 
     // make the connection graph
-    SparsityPattern sparsity;
-    if (DH::dimension < 2)
-      {
-        sparsity.reinit (dof_handler.n_dofs(level),
-                         dof_handler.n_dofs(level),
-                         dof_handler.max_couplings_between_dofs());
-        MGTools::make_sparsity_pattern (dof_handler, sparsity, level);
-        sparsity.compress();
-      }
-    else
-      {
-        DynamicSparsityPattern dsp (dof_handler.n_dofs(level),
-                                    dof_handler.n_dofs(level));
-        MGTools::make_sparsity_pattern (dof_handler, dsp, level);
-        sparsity.copy_from (dsp);
-      }
+    DynamicSparsityPattern dsp (dof_handler.n_dofs(level),
+                                dof_handler.n_dofs(level));
+    MGTools::make_sparsity_pattern (dof_handler, dsp, level);
 
-    std::vector<types::global_dof_index> new_indices(sparsity.n_rows());
-    SparsityTools::reorder_Cuthill_McKee (sparsity, new_indices,
+    std::vector<types::global_dof_index> new_indices(dsp.n_rows());
+    SparsityTools::reorder_Cuthill_McKee (dsp, new_indices,
                                           starting_indices);
 
     if (reversed_numbering)
index edd0ddfc9e4e6e9cc2845dabb63a19e18a332f9d..adf7da1a7dd565d2c1f0d3060db03a966da74e3b 100644 (file)
@@ -1520,21 +1520,10 @@ next_cell:
   template <int dim, int spacedim>
   void
   get_face_connectivity_of_cells (const Triangulation<dim,spacedim> &triangulation,
-                                  SparsityPattern          &cell_connectivity)
+                                  DynamicSparsityPattern            &cell_connectivity)
   {
-    // as built in this function, we
-    // only consider face neighbors,
-    // which leads to a fixed number of
-    // entries per row (don't forget
-    // that each cell couples with
-    // itself, and that neighbors can
-    // be refined)
     cell_connectivity.reinit (triangulation.n_active_cells(),
-                              triangulation.n_active_cells(),
-                              GeometryInfo<dim>::faces_per_cell
-                              * GeometryInfo<dim>::max_children_per_face
-                              +
-                              1);
+                              triangulation.n_active_cells());
 
     // create a map pair<lvl,idx> -> SparsityPattern index
     // TODO: we are no longer using user_indices for this because we can get
@@ -1548,22 +1537,13 @@ next_cell:
          cell != triangulation.end(); ++cell, ++index)
       indexmap[std::pair<unsigned int,unsigned int>(cell->level(),cell->index())] = index;
 
-    // next loop over all cells and
-    // their neighbors to build the
-    // sparsity pattern. note that it's
-    // a bit hard to enter all the
-    // connections when a neighbor has
-    // children since we would need to
-    // find out which of its children
-    // is adjacent to the current
-    // cell. this problem can be
-    // omitted if we only do something
-    // if the neighbor has no children
-    // -- in that case it is either on
-    // the same or a coarser level than
-    // we are. in return, we have to
-    // add entries in both directions
-    // for both cells
+    // next loop over all cells and their neighbors to build the sparsity
+    // pattern. note that it's a bit hard to enter all the connections when a
+    // neighbor has children since we would need to find out which of its
+    // children is adjacent to the current cell. this problem can be omitted
+    // if we only do something if the neighbor has no children -- in that case
+    // it is either on the same or a coarser level than we are. in return, we
+    // have to add entries in both directions for both cells
     index = 0;
     for (typename dealii::internal::ActiveCellIterator<dim, spacedim, Triangulation<dim, spacedim> >::type
          cell = triangulation.begin_active();
@@ -1581,12 +1561,22 @@ next_cell:
               cell_connectivity.add (other_index, index);
             }
       }
+  }
+
 
-    // now compress the so-built connectivity pattern
-    cell_connectivity.compress ();
+
+  template <int dim, int spacedim>
+  void
+  get_face_connectivity_of_cells (const Triangulation<dim,spacedim> &triangulation,
+                                  SparsityPattern                   &cell_connectivity)
+  {
+    DynamicSparsityPattern dsp;
+    get_face_connectivity_of_cells(triangulation, dsp);
+    cell_connectivity.copy_from(dsp);
   }
 
 
+
   template <int dim, int spacedim>
   void
   partition_triangulation (const unsigned int           n_partitions,
index 870714d4e33b2efac0e6c4f37d018ad3d2b6055a..63f0288e80ef864a5389ccd6d5402b31ac319d7d 100644 (file)
@@ -158,6 +158,11 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS
                                     Triangulation<deal_II_dimension, deal_II_space_dimension> &,
                                     const bool);
 
+    template
+      void get_face_connectivity_of_cells
+      (const Triangulation<deal_II_dimension, deal_II_space_dimension> &triangulation,
+       DynamicSparsityPattern   &cell_connectivity);
+
     template
       void get_face_connectivity_of_cells
       (const Triangulation<deal_II_dimension, deal_II_space_dimension> &triangulation,
index d28a31345e8577fa2d442952d46a20fd6e9d06ae..1fba415d92fffe1ba465a3599c4632bc65e58a11 100644 (file)
@@ -144,29 +144,20 @@ namespace SparsityTools
      * invalid_size_type indicates that a node has not been numbered yet),
      * pick a valid starting index among the as-yet unnumbered one.
      */
-    SparsityPattern::size_type
-    find_unnumbered_starting_index (const SparsityPattern     &sparsity,
-                                    const std::vector<SparsityPattern::size_type> &new_indices)
+    DynamicSparsityPattern::size_type
+    find_unnumbered_starting_index (const DynamicSparsityPattern &sparsity,
+                                    const std::vector<DynamicSparsityPattern::size_type> &new_indices)
     {
       {
-        SparsityPattern::size_type starting_point   = numbers::invalid_size_type;
-        SparsityPattern::size_type min_coordination = sparsity.n_rows();
-        for (SparsityPattern::size_type row=0; row<sparsity.n_rows(); ++row)
+        DynamicSparsityPattern::size_type starting_point   = numbers::invalid_size_type;
+        DynamicSparsityPattern::size_type min_coordination = sparsity.n_rows();
+        for (DynamicSparsityPattern::size_type row=0; row<sparsity.n_rows(); ++row)
           // look over all as-yet unnumbered indices
           if (new_indices[row] == numbers::invalid_size_type)
             {
-              SparsityPattern::iterator j = sparsity.begin(row);
-
-              // loop until we hit the end of this row's entries
-              for ( ; j<sparsity.end(row); ++j)
-                if (j->is_valid_entry() == false)
-                  break;
-              // post-condition after loop: coordination, i.e. the number of
-              // entries in this row is now j-rowstart[row]
-              if (static_cast<SparsityPattern::size_type>(j-sparsity.begin(row)) <
-                  min_coordination)
+              if (sparsity.row_length(row) < min_coordination)
                 {
-                  min_coordination = j-sparsity.begin(row);
+                  min_coordination = sparsity.row_length(row);
                   starting_point   = row;
                 }
             }
@@ -180,7 +171,7 @@ namespace SparsityTools
         // starting point, e.g. the first unnumbered one
         if (starting_point == numbers::invalid_size_type)
           {
-            for (SparsityPattern::size_type i=0; i<new_indices.size(); ++i)
+            for (DynamicSparsityPattern::size_type i=0; i<new_indices.size(); ++i)
               if (new_indices[i] == numbers::invalid_size_type)
                 {
                   starting_point = i;
@@ -197,10 +188,11 @@ namespace SparsityTools
   }
 
 
+
   void
-  reorder_Cuthill_McKee (const SparsityPattern                         &sparsity,
-                         std::vector<SparsityPattern::size_type>       &new_indices,
-                         const std::vector<SparsityPattern::size_type> &starting_indices)
+  reorder_Cuthill_McKee (const DynamicSparsityPattern                         &sparsity,
+                         std::vector<DynamicSparsityPattern::size_type>       &new_indices,
+                         const std::vector<DynamicSparsityPattern::size_type> &starting_indices)
   {
     Assert (sparsity.n_rows() == sparsity.n_cols(),
             ExcDimensionMismatch (sparsity.n_rows(), sparsity.n_cols()));
@@ -214,20 +206,20 @@ namespace SparsityTools
 
     // store the indices of the dofs renumbered in the last round. Default to
     // starting points
-    std::vector<SparsityPattern::size_type> last_round_dofs (starting_indices);
+    std::vector<DynamicSparsityPattern::size_type> last_round_dofs (starting_indices);
 
     // initialize the new_indices array with invalid values
     std::fill (new_indices.begin(), new_indices.end(),
                numbers::invalid_size_type);
 
     // delete disallowed elements
-    for (SparsityPattern::size_type i=0; i<last_round_dofs.size(); ++i)
+    for (DynamicSparsityPattern::size_type i=0; i<last_round_dofs.size(); ++i)
       if ((last_round_dofs[i]==numbers::invalid_size_type) ||
           (last_round_dofs[i]>=sparsity.n_rows()))
         last_round_dofs[i] = numbers::invalid_size_type;
 
     std::remove_if (last_round_dofs.begin(), last_round_dofs.end(),
-                    std::bind2nd(std::equal_to<SparsityPattern::size_type>(),
+                    std::bind2nd(std::equal_to<DynamicSparsityPattern::size_type>(),
                                  numbers::invalid_size_type));
 
     // now if no valid points remain: find dof with lowest coordination number
@@ -237,81 +229,55 @@ namespace SparsityTools
                                                             new_indices));
 
     // store next free dof index
-    SparsityPattern::size_type next_free_number = 0;
+    DynamicSparsityPattern::size_type next_free_number = 0;
 
     // enumerate the first round dofs
-    for (SparsityPattern::size_type i=0; i!=last_round_dofs.size(); ++i)
+    for (DynamicSparsityPattern::size_type i=0; i!=last_round_dofs.size(); ++i)
       new_indices[last_round_dofs[i]] = next_free_number++;
 
-    // now do as many steps as needed to
-    // renumber all dofs
+    // 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<SparsityPattern::size_type> next_round_dofs;
-
-        // find all neighbors of the
-        // dofs numbered in the last
-        // round
-        for (SparsityPattern::size_type i=0; i<last_round_dofs.size(); ++i)
-          for (SparsityPattern::iterator j=sparsity.begin(last_round_dofs[i]);
-               j<sparsity.end(last_round_dofs[i]); ++j)
-            if (j->is_valid_entry() == false)
-              break;
-            else
-              next_round_dofs.push_back (j->column());
+        // store the indices of the dofs to be renumbered in the next round
+        std::vector<DynamicSparsityPattern::size_type> next_round_dofs;
+
+        // find all neighbors of the dofs numbered in the last round
+        for (DynamicSparsityPattern::size_type i=0; i<last_round_dofs.size(); ++i)
+          for (DynamicSparsityPattern::row_iterator j=sparsity.row_begin(last_round_dofs[i]);
+               j<sparsity.row_end(last_round_dofs[i]); ++j)
+            next_round_dofs.push_back (*j);
 
         // sort dof numbers
         std::sort (next_round_dofs.begin(), next_round_dofs.end());
 
         // delete multiple entries
-        std::vector<SparsityPattern::size_type>::iterator end_sorted;
+        std::vector<DynamicSparsityPattern::size_type>::iterator end_sorted;
         end_sorted = std::unique (next_round_dofs.begin(), next_round_dofs.end());
         next_round_dofs.erase (end_sorted, next_round_dofs.end());
 
-        // eliminate dofs which are
-        // already numbered
+        // 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);
 
-        // check whether there are
-        // any new dofs in the
-        // list. if there are none,
-        // then we have completely
-        // numbered the current
-        // component of the
-        // graph. check if there are
-        // as yet unnumbered
-        // components of the graph
-        // that we would then have to
-        // do next
+        // check whether there are any new dofs in the list. if there are
+        // none, then we have completely numbered the current component of the
+        // graph. check if there are as yet unnumbered components of the graph
+        // that we would then have to do next
         if (next_round_dofs.empty())
           {
             if (std::find (new_indices.begin(), new_indices.end(),
                            numbers::invalid_size_type)
                 ==
                 new_indices.end())
-              // no unnumbered
-              // indices, so we can
-              // leave now
+              // no unnumbered indices, so we can leave now
               break;
 
-            // otherwise find a valid
-            // starting point for the
-            // next component of the
-            // graph and continue
-            // with numbering that
-            // one. we only do so if
-            // no starting indices
-            // were provided by the
-            // user (see the
-            // documentation of this
-            // function) so produce
-            // an error if we got
-            // here and starting
-            // indices were given
+            // otherwise find a valid starting point for the next component of
+            // the graph and continue with numbering that one. we only do so
+            // if no starting indices were provided by the user (see the
+            // documentation of this function) so produce an error if we got
+            // here and starting indices were given
             Assert (starting_indices.empty(),
                     ExcMessage ("The input graph appears to have more than one "
                                 "component, but as stated in the documentation "
@@ -326,48 +292,36 @@ namespace SparsityTools
 
 
 
-        // store for each coordination
-        // number the dofs with these
-        // coordination number
-        std::multimap<SparsityPattern::size_type, int> dofs_by_coordination;
+        // store for each coordination number the dofs with these coordination
+        // number
+        std::multimap<DynamicSparsityPattern::size_type, int> dofs_by_coordination;
 
-        // find coordination number for
-        // each of these dofs
-        for (std::vector<SparsityPattern::size_type>::iterator s=next_round_dofs.begin();
+        // find coordination number for each of these dofs
+        for (std::vector<DynamicSparsityPattern::size_type>::iterator s=next_round_dofs.begin();
              s!=next_round_dofs.end(); ++s)
           {
-            SparsityPattern::size_type coordination = 0;
-            for (SparsityPattern::iterator j=sparsity.begin(*s);
-                 j<sparsity.end(*s); ++j)
-              if (j->is_valid_entry() == false)
-                break;
-              else
-                ++coordination;
-
-            // insert this dof at its
-            // coordination number
-            const std::pair<const SparsityPattern::size_type, int> new_entry (coordination, *s);
+            DynamicSparsityPattern::size_type coordination = 0;
+            for (DynamicSparsityPattern::row_iterator j=sparsity.row_begin(*s);
+                 j<sparsity.row_end(*s); ++j)
+              ++coordination;
+
+            // insert this dof at its coordination number
+            const std::pair<const DynamicSparsityPattern::size_type, int> new_entry (coordination, *s);
             dofs_by_coordination.insert (new_entry);
           }
 
-        // assign new DoF numbers to
-        // the elements of the present
-        // front:
-        std::multimap<SparsityPattern::size_type, int>::iterator i;
+        // assign new DoF numbers to the elements of the present front:
+        std::multimap<DynamicSparsityPattern::size_type, int>::iterator i;
         for (i = dofs_by_coordination.begin(); i!=dofs_by_coordination.end(); ++i)
           new_indices[i->second] = next_free_number++;
 
-        // after that: copy this round's
-        // dofs for the next round
+        // after that: copy this round's dofs for the next round
         last_round_dofs = next_round_dofs;
       }
 
-    // test for all indices
-    // numbered. this mostly tests
-    // whether the
-    // front-marching-algorithm (which
-    // Cuthill-McKee actually is) has
-    // reached all points.
+    // test for all indices numbered. this mostly tests whether the
+    // front-marching-algorithm (which Cuthill-McKee actually is) has reached
+    // all points.
     Assert ((std::find (new_indices.begin(), new_indices.end(), numbers::invalid_size_type)
              ==
              new_indices.end())
@@ -376,6 +330,25 @@ namespace SparsityTools
             ExcInternalError());
   }
 
+
+
+  void
+  reorder_Cuthill_McKee (const SparsityPattern                   &sparsity,
+                         std::vector<SparsityPattern::size_type> &new_indices,
+                         const std::vector<SparsityPattern::size_type> &starting_indices)
+  {
+    DynamicSparsityPattern dsp(sparsity.n_rows(), sparsity.n_cols());
+    for (unsigned int row=0; row<sparsity.n_rows(); ++row)
+      {
+        for (SparsityPattern::iterator it=sparsity.begin(row); it!=sparsity.end(row)
+               && it->is_valid_entry() ; ++it)
+          dsp.add(row, it->column());
+      }
+    reorder_Cuthill_McKee(dsp, new_indices, starting_indices);
+  }
+
+
+
 #ifdef DEAL_II_WITH_MPI
   template <class DSP_t>
   void distribute_sparsity_pattern(DSP_t &dsp,
index b07676abc950ef52bb9b812e624a7d37f57fbda3..55cf7110cc5a70c6dc972e90fd845bf8b68dee60 100644 (file)
@@ -25,7 +25,7 @@
 #include <deal.II/grid/grid_generator.h>
 #include <deal.II/grid/grid_tools.h>
 
-#include <deal.II/lac/sparsity_pattern.h>
+#include <deal.II/lac/dynamic_sparsity_pattern.h>
 #include <deal.II/lac/sparsity_tools.h>
 
 #include <fstream>
@@ -56,7 +56,7 @@ void LaplaceProblem<dim>::run ()
   std::ifstream input_file(SOURCE_DIR "/gerold_1.inp");
   grid_in.read_ucd(input_file);
 
-  SparsityPattern cell_connectivity;
+  DynamicSparsityPattern cell_connectivity;
   GridTools::get_face_connectivity_of_cells (triangulation,
                                              cell_connectivity);
   std::vector<types::global_dof_index> permutation(triangulation.n_active_cells());
index d99c32280b628f64fe1e4ce518482f1bc676337f..a6f232ea195854e356da90cc93c44aea62d7264e 100644 (file)
@@ -21,7 +21,7 @@
 
 #include "../tests.h"
 #include <deal.II/lac/sparsity_tools.h>
-#include <deal.II/lac/sparsity_pattern.h>
+#include <deal.II/lac/dynamic_sparsity_pattern.h>
 
 #include <iomanip>
 #include <fstream>
@@ -35,24 +35,20 @@ int main ()
   deallog.depth_console(0);
   deallog.threshold_double(1.e-10);
 
-  SparsityPattern sp (4,4,4);
+  DynamicSparsityPattern dsp (4,4);
   for (unsigned int i=0; i<4; ++i)
-    sp.add (i,i);
+    dsp.add (i,i);
 
-  // create a graph with components
-  // 0,2 and 1,3 that are
-  // disconnected
-  sp.add (0,2);
-  sp.add (2,0);
+  // create a graph with components 0,2 and 1,3 that are disconnected
+  dsp.add (0,2);
+  dsp.add (2,0);
 
-  sp.add (1,3);
-  sp.add (3,1);
-
-  sp.compress ();
+  dsp.add (1,3);
+  dsp.add (3,1);
 
   // now find permutation
   std::vector<types::global_dof_index> permutation(4);
-  SparsityTools::reorder_Cuthill_McKee (sp, permutation);
+  SparsityTools::reorder_Cuthill_McKee (dsp, permutation);
 
   for (unsigned int i=0; i<permutation.size(); ++i)
     deallog << permutation[i] << std::endl;

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.