]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Split the SparsityPattern::partition function into a separate function in a new names...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 31 Oct 2008 16:13:19 +0000 (16:13 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 31 Oct 2008 16:13:19 +0000 (16:13 +0000)
git-svn-id: https://svn.dealii.org/trunk@17425 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/source/grid/grid_tools.cc
deal.II/doc/news/changes.h
deal.II/lac/include/lac/sparsity_pattern.h
deal.II/lac/include/lac/sparsity_tools.h [new file with mode: 0644]
deal.II/lac/source/sparsity_pattern.cc
deal.II/lac/source/sparsity_tools.cc [new file with mode: 0644]

index cfa7f16573b3a314aeda925bf8fe1452849aa7ca..7a3861dbec063f8544777d457c4e91b9238917c4 100644 (file)
@@ -19,6 +19,7 @@
 #include <grid/tria_iterator.h>
 #include <grid/intergrid_map.h>
 #include <lac/sparsity_pattern.h>
+#include <lac/sparsity_tools.h>
 #include <lac/compressed_sparsity_pattern.h>
 #include <dofs/dof_handler.h>
 #include <dofs/dof_accessor.h>
@@ -942,7 +943,7 @@ partition_triangulation (const unsigned int     n_partitions,
                                    // of freedom (which is associated with a
                                    // cell)
   std::vector<unsigned int> partition_indices (triangulation.n_active_cells());
-  cell_connection_graph.partition (n_partitions,  partition_indices);
+  SparsityTools::partition (cell_connection_graph, n_partitions,  partition_indices);
 
                                    // finally loop over all cells and set the
                                    // subdomain ids
index cf68f7ea23ba6d008ac9e39b56658a20ae7f2079..bf6e548013dc02b1396dc5a8fca5ad9012659978 100644 (file)
@@ -299,6 +299,15 @@ inconvenience this causes.
 <h3>lac</h3>
 
 <ol>
+  <li>
+  <p>
+  Changed: The function SparsityPattern::partition has been deprecated. It
+  is now available in a new namespace SparsityTools that collects algorithms
+  that work on sparsity patterns or connectivity graphs.
+  <br>
+  (WB 2008/10/31)
+  </p>
+
   <li>
   <p>
   Fixed: Whereas the Vector class copy operator resized the left hand side
index 6b29135671b659df5bc87acf95a55bc5d068f955..fdb6b5f534bec587fdcec646b67dbf4999c2e3da 100644 (file)
@@ -1205,6 +1205,11 @@ class SparsityPattern : public Subscriptor
     bool stores_only_added_elements () const;
     
                                      /**
+                                     * @deprecated
+                                     *
+                                     * This function is deprecated. Use
+                                     * SparsityTools::partition instead.
+                                     *
                                       * Use the METIS partitioner to generate
                                       * a partitioning of the degrees of
                                       * freedom represented by this sparsity
@@ -1472,10 +1477,6 @@ class SparsityPattern : public Subscriptor
                                      /**
                                       * Exception
                                       */
-    DeclException0 (ExcMETISNotInstalled);
-                                     /**
-                                      * Exception
-                                      */
     DeclException1 (ExcInvalidNumberOfPartitions,
                     int,
                     << "The number of partitions you gave is " << arg1
diff --git a/deal.II/lac/include/lac/sparsity_tools.h b/deal.II/lac/include/lac/sparsity_tools.h
new file mode 100644 (file)
index 0000000..b448182
--- /dev/null
@@ -0,0 +1,124 @@
+//---------------------------------------------------------------------------
+//    $Id: sparsity_tools.h 16733 2008-09-03 14:22:26Z heister $
+//    Version: $Name$
+//
+//    Copyright (C) 2008 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//---------------------------------------------------------------------------
+#ifndef __deal2__sparsity_tools_h
+#define __deal2__sparsity_tools_h
+
+
+#include <base/config.h>
+#include <base/exceptions.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+class SparsityPattern;
+
+
+
+/*! @addtogroup Sparsity
+ *@{
+ */
+
+/**
+ * A namespace for functions that deal with things that one can do on sparsity
+ * patterns, such as renumbering rows and columns (or degrees of freedom if
+ * you want) according to the connectivity, or partitioning degrees of
+ * freedom.
+*/
+namespace SparsityTools
+{
+                                  /**
+                                   * Use the METIS partitioner to generate
+                                   * a partitioning of the degrees of
+                                   * freedom represented by this sparsity
+                                   * pattern. In effect, we view this
+                                   * sparsity pattern as a graph of
+                                   * connections between various degrees of
+                                   * freedom, where each nonzero entry in
+                                   * the sparsity pattern corresponds to an
+                                   * edge between two nodes in the
+                                   * connection graph. The goal is then to
+                                   * decompose this graph into groups of
+                                   * nodes so that a minimal number of
+                                   * edges are cut by the boundaries
+                                   * between node groups. This partitioning
+                                   * is done by METIS. Note that METIS can
+                                   * only partition symmetric sparsity
+                                   * patterns, and that of course the
+                                   * sparsity pattern has to be square. We
+                                   * do not check for symmetry of the
+                                   * sparsity pattern, since this is an
+                                   * expensive operation, but rather leave
+                                   * this as the responsibility of caller
+                                   * of this function.
+                                   *
+                                   * After calling this function, the
+                                   * output array will have values between
+                                   * zero and @p n_partitions-1 for each
+                                   * node (i.e. row or column of the
+                                   * matrix).
+                                   *
+                                   * This function will generate an error
+                                   * if METIS is not installed unless
+                                   * @p n_partitions is one. I.e., you can
+                                   * write a program so that it runs in the
+                                   * single-processor single-partition case
+                                   * without METIS installed, and only
+                                   * requires METIS when multiple
+                                   * partitions are required.
+                                   *
+                                   * Note that the sparsity pattern itself
+                                   * is not changed by calling this
+                                   * function. However, you will likely use
+                                   * the information generated by calling
+                                   * this function to renumber degrees of
+                                   * freedom, after which you will of
+                                   * course have to regenerate the sparsity
+                                   * pattern.
+                                   *
+                                   * This function will rarely be called
+                                   * separately, since in finite element
+                                   * methods you will want to partition the
+                                   * mesh, not the matrix. This can be done
+                                   * by calling
+                                   * @p GridTools::partition_triangulation.
+                                   */
+  void partition (const SparsityPattern     &sparsity_pattern,
+                 const unsigned int         n_partitions,
+                 std::vector<unsigned int> &partition_indices);
+
+                                  /**
+                                   * Exception
+                                   */
+  DeclException0 (ExcMETISNotInstalled);
+                                  /**
+                                   * Exception
+                                   */
+  DeclException1 (ExcInvalidNumberOfPartitions,
+                 int,
+                 << "The number of partitions you gave is " << arg1
+                 << ", but must be greater than zero.");
+                                  /**
+                                   * Exception
+                                   */
+  DeclException2 (ExcInvalidArraySize,
+                 int, int,
+                 << "The array has size " << arg1 << " but should have size "
+                 << arg2);
+}
+
+/**
+ *@}
+ */
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index 27615135bb90cdfd7311ace0b2fbd2534f9732f4..4165e97f32f05002ec26e6f4e14bae1640e89252 100644 (file)
@@ -13,6 +13,7 @@
 
 #include <base/vector_slice.h>
 #include <lac/sparsity_pattern.h>
+#include <lac/sparsity_tools.h>
 #include <lac/full_matrix.h>
 #include <lac/compressed_sparsity_pattern.h>
 #include <lac/compressed_set_sparsity_pattern.h>
 
 DEAL_II_NAMESPACE_OPEN
 
-#ifdef DEAL_II_USE_METIS
-extern "C" {
-#  include <metis.h>
-}
-#endif
 
 const unsigned int SparsityPattern::invalid_entry;
 
@@ -1118,73 +1114,7 @@ SparsityPattern::
 partition (const unsigned int         n_partitions,
            std::vector<unsigned int> &partition_indices) const
 {
-  Assert (rows==cols, ExcNotQuadratic());
-  Assert (is_compressed(), ExcNotCompressed());
-
-  Assert (n_partitions > 0, ExcInvalidNumberOfPartitions(n_partitions));
-  Assert (partition_indices.size() == rows,
-          ExcInvalidArraySize (partition_indices.size(),rows));
-
-                                   // check for an easy return
-  if (n_partitions == 1)
-    {
-      std::fill_n (partition_indices.begin(), partition_indices.size(), 0U);
-      return;
-    }
-
-                                   // Make sure that METIS is actually
-                                   // installed and detected
-#ifndef DEAL_II_USE_METIS
-  AssertThrow (false, ExcMETISNotInstalled());
-#else
-
-                                   // generate the data structures for
-                                   // METIS. Note that this is particularly
-                                   // simple, since METIS wants exactly our
-                                   // compressed row storage format. we only
-                                   // have to set up a few auxiliary arrays
-  int
-    n = static_cast<signed int>(rows),
-    wgtflag = 0,                          // no weights on nodes or edges
-    numflag = 0,                          // C-style 0-based numbering
-    nparts  = static_cast<int>(n_partitions), // number of subdomains to create
-    dummy;                                // the numbers of edges cut by the
-                                          //   resulting partition
-
-                                   // use default options for METIS
-  int options[5] = { 0,0,0,0,0 };
-
-                                   // one more nuisance: we have to copy our
-                                   // own data to arrays that store signed
-                                   // integers :-(
-  std::vector<idxtype> int_rowstart (rowstart, rowstart+cols+1);
-  std::vector<idxtype> int_colnums (colnums, colnums+max_vec_len+1);
-
-  std::vector<idxtype> int_partition_indices (rows);
-
-                                   // Select which type of partitioning to create
-
-                                   // Use recursive if the number of
-                                   // partitions is less than or equal to 8
-  if (n_partitions <= 8)
-    METIS_PartGraphRecursive(&n, &int_rowstart[0], &int_colnums[0],
-                             NULL, NULL,
-                             &wgtflag, &numflag, &nparts, &options[0],
-                             &dummy, &int_partition_indices[0]);
-  
-                                   // Otherwise  use kway
-  else
-    METIS_PartGraphKway(&n, &int_rowstart[0], &int_colnums[0],
-                        NULL, NULL,
-                        &wgtflag, &numflag, &nparts, &options[0],
-                        &dummy, &int_partition_indices[0]);
-
-                                   // now copy back generated indices into the
-                                   // output array
-  std::copy (int_partition_indices.begin(),
-             int_partition_indices.end(),
-             partition_indices.begin());
-#endif
+  SparsityTools::partition (*this, n_partitions, partition_indices);
 }
 
 
diff --git a/deal.II/lac/source/sparsity_tools.cc b/deal.II/lac/source/sparsity_tools.cc
new file mode 100644 (file)
index 0000000..9d844aa
--- /dev/null
@@ -0,0 +1,112 @@
+//---------------------------------------------------------------------------
+//    $Id: sparsity_pattern.cc 16733 2008-09-03 14:22:26Z heister $
+//    Version: $Name$
+//
+//    Copyright (C) 2008 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//---------------------------------------------------------------------------
+
+#include <base/exceptions.h>
+#include <lac/exceptions.h>
+#include <lac/sparsity_pattern.h>
+#include <lac/sparsity_tools.h>
+
+#ifdef DEAL_II_USE_METIS
+extern "C" {
+#  include <metis.h>
+}
+#endif
+
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace SparsityTools
+{
+  
+  void partition (const SparsityPattern     &sparsity_pattern,
+                 const unsigned int         n_partitions,
+                 std::vector<unsigned int> &partition_indices)
+  {
+    Assert (sparsity_pattern.n_rows()==sparsity_pattern.n_cols(),
+           ExcNotQuadratic());
+    Assert (sparsity_pattern.is_compressed(),
+           SparsityPattern::ExcNotCompressed());
+
+    Assert (n_partitions > 0, ExcInvalidNumberOfPartitions(n_partitions));
+    Assert (partition_indices.size() == sparsity_pattern.n_rows(),
+           ExcInvalidArraySize (partition_indices.size(),
+                                sparsity_pattern.n_rows()));
+
+                                    // check for an easy return
+    if (n_partitions == 1)
+      {
+       std::fill_n (partition_indices.begin(), partition_indices.size(), 0U);
+       return;
+      }
+
+                                    // Make sure that METIS is actually
+                                    // installed and detected
+#ifndef DEAL_II_USE_METIS
+    AssertThrow (false, ExcMETISNotInstalled());
+#else
+
+                                    // generate the data structures for
+                                    // METIS. Note that this is particularly
+                                    // simple, since METIS wants exactly our
+                                    // compressed row storage format. we only
+                                    // have to set up a few auxiliary arrays
+    int
+      n = static_cast<signed int>(sparsity_pattern.n_rows()),
+      wgtflag = 0,                          // no weights on nodes or edges
+      numflag = 0,                          // C-style 0-based numbering
+      nparts  = static_cast<int>(n_partitions), // number of subdomains to create
+      dummy;                                // the numbers of edges cut by the
+                                    //   resulting partition
+
+                                    // use default options for METIS
+    int options[5] = { 0,0,0,0,0 };
+
+                                    // one more nuisance: we have to copy our
+                                    // own data to arrays that store signed
+                                    // integers :-(
+    std::vector<idxtype> int_rowstart (sparsity_pattern.get_rowstart_indices(),
+                                      sparsity_pattern.get_rowstart_indices() +
+                                      sparsity_pattern.n_cols()+1);
+    std::vector<idxtype> int_colnums (sparsity_pattern.get_colnums(),
+                                     sparsity_pattern.get_colnums()+max_vec_len+1);
+
+    std::vector<idxtype> int_partition_indices (sparsity_pattern.n_rows());
+
+                                    // Select which type of partitioning to create
+
+                                    // Use recursive if the number of
+                                    // partitions is less than or equal to 8
+    if (n_partitions <= 8)
+      METIS_PartGraphRecursive(&n, &int_rowstart[0], &int_colnums[0],
+                              NULL, NULL,
+                              &wgtflag, &numflag, &nparts, &options[0],
+                              &dummy, &int_partition_indices[0]);
+  
+                                    // Otherwise  use kway
+    else
+      METIS_PartGraphKway(&n, &int_rowstart[0], &int_colnums[0],
+                         NULL, NULL,
+                         &wgtflag, &numflag, &nparts, &options[0],
+                         &dummy, &int_partition_indices[0]);
+
+                                    // now copy back generated indices into the
+                                    // output array
+    std::copy (int_partition_indices.begin(),
+              int_partition_indices.end(),
+              partition_indices.begin());
+#endif
+  }
+
+}
+
+DEAL_II_NAMESPACE_CLOSE

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.