]> https://gitweb.dealii.org/ - dealii.git/commitdiff
color_sparsity_function added in sparsity tools
authornivesh <nidhinivesh.d@gmail.com>
Tue, 5 Dec 2017 13:11:14 +0000 (14:11 +0100)
committernivesh <nidhinivesh.d@gmail.com>
Tue, 5 Dec 2017 13:11:14 +0000 (14:11 +0100)
doc/doxygen/images/color_undirected_graph.png [new file with mode: 0755]
doc/news/changes/minor/20171204NiveshD [new file with mode: 0644]
include/deal.II/lac/sparsity_tools.h
source/lac/sparsity_tools.cc

diff --git a/doc/doxygen/images/color_undirected_graph.png b/doc/doxygen/images/color_undirected_graph.png
new file mode 100755 (executable)
index 0000000..de5e098
Binary files /dev/null and b/doc/doxygen/images/color_undirected_graph.png differ
diff --git a/doc/news/changes/minor/20171204NiveshD b/doc/news/changes/minor/20171204NiveshD
new file mode 100644 (file)
index 0000000..8f50f7e
--- /dev/null
@@ -0,0 +1,4 @@
+New: Add color_sparsity_pattern function to color a graph represented by SparsityPattern Object.
+The function uses coloring function from ZOLTAN library.
+<br>
+(Nivesh Dommaraju, 2017/12/04)
index e76d6d4b8c1629eefc9d2c160ba7f9d07d263dd9..0a2ee7cc23240ac6fa6209d64b32367d0c4c6294 100644 (file)
@@ -71,21 +71,22 @@ namespace SparsityTools
    * 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 or ZOLTAN, depending upon which partitioner is chosen in the fourth argument.
-   * The default is METIS. Note that METIS and ZOLTAN 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.
+   * done by METIS or ZOLTAN, depending upon which partitioner is chosen
+   * in the fourth argument. The default is METIS. Note that METIS and
+   * ZOLTAN 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).
    *
-   * If deal.II was not installed with packages ZOLTAN or METIS, this function will generate an error
-   * when corresponding partition method is chosen, 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 the packages installed, and only requires them installed when
+   * If deal.II was not installed with packages ZOLTAN or METIS, this
+   * function will generate an error when corresponding partition method
+   * is chosen, 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
+   * the packages installed, and only requires them installed when
    * multiple partitions are required.
    *
    * Note that the sparsity pattern itself is not changed by calling this
@@ -103,6 +104,38 @@ namespace SparsityTools
                   const Partitioner partitioner = Partitioner::metis
                  );
 
+  /**
+   * Using a coloring algorithm provided by ZOLTAN to color nodes whose
+   * connections are represented using a SparsityPattern object. In effect,
+   * we view this sparsity pattern as a graph of connections between nodes,
+   * where each nonzero entry in the @p sparsity_pattern corresponds to
+   * an edge between two nodes. The goal is to assign each node a color index
+   * such that no two directly connected nodes have the same color.
+   * The assigned colors are listed in @p color_indices indexed from one and
+   * the function returns total number of colors used. Note that this function
+   * requires that SparsityPattern be symmetric (and hence square) i.e an
+   * undirected graph representation. 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.
+   *
+   * ZOLTAN coloring algorithm is run in serial by each processor. Hence all
+   * processors have coloring information of all the nodes.
+   *
+   * @image html color_undirected_graph.png
+   * For example, the above image shows 5 nodes numbered from 0 to 4.
+   * The connections shown are bidirectional. After coloring, no two connected
+   * nodes have the same color.
+   *
+   * If deal.II was not installed with package ZOLTAN, this function will
+   * generate an error.
+   *
+   * @note The current function is an alternative to
+   * GraphColoring::make_graph_coloring() which is tailored to graph
+   * coloring arising in shared-memory parallel assembly of matrices.
+   */
+  unsigned int color_sparsity_pattern (const SparsityPattern     &sparsity_pattern,
+                                       std::vector<unsigned int> &color_indices);
+
   /**
    * For a given sparsity pattern, compute a re-enumeration of row/column
    * indices based on the algorithm by Cuthill-McKee.
index 995cdf7a0bc1b615147682ee05b73d9d4a64983a..ea3d33bd0ef408aba3e7bb7aabca84b50a695ca8 100644 (file)
@@ -165,8 +165,8 @@ namespace SparsityTools
 
       for (unsigned int i=0; i < n_dofs; i++)
         {
-          globalID[i] = i;    //global ids run from 0 to n_dofs-1
-          localID[i] = i;     //same as global ids
+          globalID[i] = i;
+          localID[i] = i;     //Same as global ids.
         }
     }
 
@@ -183,7 +183,12 @@ namespace SparsityTools
       Assert ( numEdges != nullptr , ExcInternalError() );
 
       for (int i=0; i<num_obj; ++i)
-        numEdges[i] = graph->row_length(globalID[i])-1;
+        {
+          if ( graph->exists(i,i) )     //Check if diagonal element is present
+            numEdges[i] = graph->row_length(globalID[i])-1;
+          else
+            numEdges[i] = graph->row_length(globalID[i]);
+        }
     }
 
 
@@ -326,6 +331,68 @@ namespace SparsityTools
   }
 
 
+  unsigned int color_sparsity_pattern (const SparsityPattern     &sparsity_pattern,
+                                       std::vector<unsigned int> &color_indices)
+  {
+    // Make sure that ZOLTAN is actually
+    // installed and detected
+#ifndef DEAL_II_TRILINOS_WITH_ZOLTAN
+    (void)sparsity_pattern;
+    AssertThrow (false, ExcZOLTANNotInstalled());
+    return 0;
+#else
+    //coloring algorithm is run in serial by each processor.
+    std::unique_ptr<Zoltan> zz = std_cxx14::make_unique<Zoltan> (MPI_COMM_SELF);
+
+    //Coloring parameters
+    zz->Set_Param ("COLORING_PROBLEM", "DISTANCE-1");  //Standard coloring
+    zz->Set_Param ("NUM_GID_ENTRIES", "1");  // 1 entry represents global ID
+    zz->Set_Param ("NUM_LID_ENTRIES", "1");  // 1 entry represents local ID
+    zz->Set_Param ("OBJ_WEIGHT_DIM", "0");   // object weights not used
+    zz->Set_Param ("RECOLORING_NUM_OF_ITERATIONS", "0");
+    zz->Set_Param ("DEBUG_LEVEL", "0" );     //level of debug info
+
+    //Zoltan::Color function requires a non-const SparsityPattern object
+    SparsityPattern graph;
+    graph.copy_from (sparsity_pattern);
+
+    //Set query functions required by coloring function
+    zz->Set_Num_Obj_Fn (get_number_of_objects, &graph);
+    zz->Set_Obj_List_Fn (get_object_list, &graph);
+    zz->Set_Num_Edges_Multi_Fn (get_num_edges_list, &graph);
+    zz->Set_Edge_List_Multi_Fn (get_edge_list, &graph);
+
+    //Variables needed by coloring function
+    int num_gid_entries=1;
+    const int num_objects = graph.n_rows();
+
+    //Preallocate input variables. Element type fixed by ZOLTAN.
+    std::vector<ZOLTAN_ID_TYPE> global_ids (num_objects);
+    std::vector<int> color_exp (num_objects);
+
+    //Set ids for which coloring needs to be done
+    for (int i=0; i < num_objects; i++)
+      global_ids[i] = i;
+
+    //Call ZOLTAN coloring algorithm
+    int rc = zz->Color (num_gid_entries,
+                        num_objects,
+                        global_ids.data(),
+                        color_exp.data());
+
+    //Check for error code
+    Assert (rc == ZOLTAN_OK , ExcInternalError());
+
+    //Allocate and assign color
+    std::copy (color_exp.begin(), color_exp.end(), color_indices.begin());
+
+    unsigned int n_colors = *(std::max_element (color_indices.begin(),
+                                                color_indices.end()));
+    return n_colors;
+#endif
+  }
+
+
   namespace internal
   {
     /**

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.