Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
graph_coloring.h
1 
2 // ---------------------------------------------------------------------
3 //
4 // Copyright (C) 2013 - 2019 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE.md at
13 // the top level directory of deal.II.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef dealii_graph_coloring_h
18 # define dealii_graph_coloring_h
19 
20 
21 # include <deal.II/base/config.h>
22 
23 # include <deal.II/base/thread_management.h>
24 
25 # include <deal.II/lac/sparsity_tools.h>
26 
27 DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
28 # include <boost/unordered_map.hpp>
29 # include <boost/unordered_set.hpp>
30 DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
31 
32 # include <functional>
33 # include <set>
34 # include <vector>
35 
36 
37 DEAL_II_NAMESPACE_OPEN
38 
42 namespace GraphColoring
43 {
44  namespace internal
45  {
55  inline bool
56  have_nonempty_intersection(
57  const std::vector<types::global_dof_index> &indices1,
58  const std::vector<types::global_dof_index> &indices2)
59  {
60  // we assume that both arrays are sorted, so we can walk
61  // them in lockstep and see if we encounter an index that's
62  // in both arrays. once we reach the end of either array,
63  // we know that there is no intersection
64  std::vector<types::global_dof_index>::const_iterator p = indices1.begin(),
65  q = indices2.begin();
66  while ((p != indices1.end()) && (q != indices2.end()))
67  {
68  if (*p < *q)
69  ++p;
70  else if (*p > *q)
71  ++q;
72  else
73  // conflict found!
74  return true;
75  }
76 
77  // no conflict found!
78  return false;
79  }
80 
81 
112  template <typename Iterator>
113  std::vector<std::vector<Iterator>>
114  create_partitioning(
115  const Iterator & begin,
116  const typename identity<Iterator>::type &end,
117  const std::function<std::vector<types::global_dof_index>(
118  const Iterator &)> & get_conflict_indices)
119  {
120  // Number of iterators.
121  unsigned int n_iterators = 0;
122 
123  // Create a map from conflict indices to iterators
124  boost::unordered_map<types::global_dof_index, std::vector<Iterator>>
125  indices_to_iterators;
126  for (Iterator it = begin; it != end; ++it)
127  {
128  const std::vector<types::global_dof_index> conflict_indices =
129  get_conflict_indices(it);
130  const unsigned int n_conflict_indices = conflict_indices.size();
131  for (unsigned int i = 0; i < n_conflict_indices; ++i)
132  indices_to_iterators[conflict_indices[i]].push_back(it);
133  ++n_iterators;
134  }
135 
136  // create the very first zone which contains only the first
137  // iterator. then create the other zones. keep track of all the
138  // iterators that have already been assigned to a zone
139  std::vector<std::vector<Iterator>> zones(1,
140  std::vector<Iterator>(1, begin));
141  std::set<Iterator> used_it;
142  used_it.insert(begin);
143  while (used_it.size() != n_iterators)
144  {
145  // loop over the elements of the previous zone. for each element of
146  // the previous zone, get the conflict indices and from there get
147  // those iterators that are conflicting with the current element
148  typename std::vector<Iterator>::iterator previous_zone_it(
149  zones.back().begin());
150  typename std::vector<Iterator>::iterator previous_zone_end(
151  zones.back().end());
152  std::vector<Iterator> new_zone;
153  for (; previous_zone_it != previous_zone_end; ++previous_zone_it)
154  {
155  const std::vector<types::global_dof_index> conflict_indices =
156  get_conflict_indices(*previous_zone_it);
157 
158  const unsigned int n_conflict_indices(conflict_indices.size());
159  for (unsigned int i = 0; i < n_conflict_indices; ++i)
160  {
161  const std::vector<Iterator> &conflicting_elements =
162  indices_to_iterators[conflict_indices[i]];
163  for (unsigned int j = 0; j < conflicting_elements.size(); ++j)
164  {
165  // check that the iterator conflicting with the current
166  // one is not associated to a zone yet and if so, assign
167  // it to the current zone. mark it as used
168  //
169  // we can shortcut this test if the conflicting iterator
170  // is the current iterator
171  if ((conflicting_elements[j] != *previous_zone_it) &&
172  (used_it.count(conflicting_elements[j]) == 0))
173  {
174  new_zone.push_back(conflicting_elements[j]);
175  used_it.insert(conflicting_elements[j]);
176  }
177  }
178  }
179  }
180 
181  // If there are iterators in the new zone, then the zone is added to
182  // the partition. Otherwise, the graph is disconnected and we need to
183  // find an iterator on the other part of the graph. start the whole
184  // process again with the first iterator that hasn't been assigned to
185  // a zone yet
186  if (new_zone.size() != 0)
187  zones.push_back(new_zone);
188  else
189  for (Iterator it = begin; it != end; ++it)
190  if (used_it.count(it) == 0)
191  {
192  zones.push_back(std::vector<Iterator>(1, it));
193  used_it.insert(it);
194  break;
195  }
196  }
197 
198  return zones;
199  }
200 
201 
202 
222  template <typename Iterator>
223  void
224  make_dsatur_coloring(
225  std::vector<Iterator> & partition,
226  const std::function<std::vector<types::global_dof_index>(
227  const Iterator &)> & get_conflict_indices,
228  std::vector<std::vector<Iterator>> &partition_coloring)
229  {
230  partition_coloring.clear();
231 
232  // Number of zones composing the partitioning.
233  const unsigned int partition_size(partition.size());
234  std::vector<unsigned int> sorted_vertices(partition_size);
235  std::vector<int> degrees(partition_size);
236  std::vector<std::vector<types::global_dof_index>> conflict_indices(
237  partition_size);
238  std::vector<std::vector<unsigned int>> graph(partition_size);
239 
240  // Get the conflict indices associated to each iterator. The
241  // conflict_indices have to be sorted so we can more easily find conflicts
242  // later on
243  for (unsigned int i = 0; i < partition_size; ++i)
244  {
245  conflict_indices[i] = get_conflict_indices(partition[i]);
246  std::sort(conflict_indices[i].begin(), conflict_indices[i].end());
247  }
248 
249  // Compute the degree of each vertex of the graph using the
250  // intersection of the conflict indices.
251  for (unsigned int i = 0; i < partition_size; ++i)
252  for (unsigned int j = i + 1; j < partition_size; ++j)
253  // If the two iterators share indices then we increase the degree of
254  // the vertices and create an ''edge'' in the graph.
255  if (have_nonempty_intersection(conflict_indices[i],
256  conflict_indices[j]))
257  {
258  ++degrees[i];
259  ++degrees[j];
260  graph[i].push_back(j);
261  graph[j].push_back(i);
262  }
263 
264  // Sort the vertices by decreasing degree.
265  std::vector<int>::iterator degrees_it;
266  for (unsigned int i = 0; i < partition_size; ++i)
267  {
268  // Find the largest element.
269  degrees_it = std::max_element(degrees.begin(), degrees.end());
270  sorted_vertices[i] = degrees_it - degrees.begin();
271  // Put the largest element to -1 so it cannot be chosen again.
272  *degrees_it = -1;
273  }
274 
275  // Color the graph.
276  std::vector<boost::unordered_set<unsigned int>> colors_used;
277  for (unsigned int i = 0; i < partition_size; ++i)
278  {
279  const unsigned int current_vertex(sorted_vertices[i]);
280  bool new_color(true);
281  // Try to use an existing color, i.e., try to find a color which is
282  // not associated to one of the vertices linked to current_vertex.
283  // Loop over the color.
284  for (unsigned int j = 0; j < partition_coloring.size(); ++j)
285  {
286  // Loop on the vertices linked to current_vertex. If one vertex
287  // linked to current_vertex is already using the color j, this
288  // color cannot be used anymore.
289  bool unused_color(true);
290  for (const auto adjacent_vertex : graph[current_vertex])
291  if (colors_used[j].count(adjacent_vertex) == 1)
292  {
293  unused_color = false;
294  break;
295  }
296  if (unused_color)
297  {
298  partition_coloring[j].push_back(partition[current_vertex]);
299  colors_used[j].insert(current_vertex);
300  new_color = false;
301  break;
302  }
303  }
304  // Add a new color.
305  if (new_color)
306  {
307  partition_coloring.push_back(
308  std::vector<Iterator>(1, partition[current_vertex]));
309  boost::unordered_set<unsigned int> tmp;
310  tmp.insert(current_vertex);
311  colors_used.push_back(tmp);
312  }
313  }
314  }
315 
316 
317 
327  template <typename Iterator>
328  std::vector<std::vector<Iterator>>
329  gather_colors(
330  const std::vector<std::vector<std::vector<Iterator>>> &partition_coloring)
331  {
332  std::vector<std::vector<Iterator>> coloring;
333 
334  // Count the number of iterators in each color.
335  const unsigned int partition_size(partition_coloring.size());
336  std::vector<std::vector<unsigned int>> colors_counter(partition_size);
337  for (unsigned int i = 0; i < partition_size; ++i)
338  {
339  const unsigned int n_colors(partition_coloring[i].size());
340  colors_counter[i].resize(n_colors);
341  for (unsigned int j = 0; j < n_colors; ++j)
342  colors_counter[i][j] = partition_coloring[i][j].size();
343  }
344 
345  // Find the partition with the largest number of colors for the even
346  // partition.
347  unsigned int i_color(0);
348  unsigned int max_even_n_colors(0);
349  const unsigned int colors_size(colors_counter.size());
350  for (unsigned int i = 0; i < colors_size; i += 2)
351  {
352  if (max_even_n_colors < colors_counter[i].size())
353  {
354  max_even_n_colors = colors_counter[i].size();
355  i_color = i;
356  }
357  }
358  coloring.resize(max_even_n_colors);
359  for (unsigned int j = 0; j < colors_counter[i_color].size(); ++j)
360  coloring[j] = partition_coloring[i_color][j];
361 
362  for (unsigned int i = 0; i < partition_size; i += 2)
363  {
364  if (i != i_color)
365  {
366  boost::unordered_set<unsigned int> used_k;
367  for (unsigned int j = 0; j < colors_counter[i].size(); ++j)
368  {
369  // Find the color in the current partition with the largest
370  // number of iterators.
371  std::vector<unsigned int>::iterator it;
372  it = std::max_element(colors_counter[i].begin(),
373  colors_counter[i].end());
374  unsigned int min_iterators(static_cast<unsigned int>(-1));
375  unsigned int pos(0);
376  // Find the color of coloring with the least number of colors
377  // among the colors that have not been used yet.
378  for (unsigned int k = 0; k < max_even_n_colors; ++k)
379  if (used_k.count(k) == 0)
380  if (colors_counter[i_color][k] < min_iterators)
381  {
382  min_iterators = colors_counter[i_color][k];
383  pos = k;
384  }
385  colors_counter[i_color][pos] += *it;
386  // Concatenate the current color with the existing coloring.
387  coloring[pos].insert(
388  coloring[pos].end(),
389  partition_coloring[i][it - colors_counter[i].begin()]
390  .begin(),
391  partition_coloring[i][it - colors_counter[i].begin()]
392  .end());
393  used_k.insert(pos);
394  // Put the number of iterators to the current color to zero.
395  *it = 0;
396  }
397  }
398  }
399 
400  // If there is more than one partition, do the same thing that we did for
401  // the even partitions to the odd partitions
402  if (partition_size > 1)
403  {
404  unsigned int max_odd_n_colors(0);
405  for (unsigned int i = 1; i < partition_size; i += 2)
406  {
407  if (max_odd_n_colors < colors_counter[i].size())
408  {
409  max_odd_n_colors = colors_counter[i].size();
410  i_color = i;
411  }
412  }
413  coloring.resize(max_even_n_colors + max_odd_n_colors);
414  for (unsigned int j = 0; j < colors_counter[i_color].size(); ++j)
415  coloring[max_even_n_colors + j] = partition_coloring[i_color][j];
416 
417  for (unsigned int i = 1; i < partition_size; i += 2)
418  {
419  if (i != i_color)
420  {
421  boost::unordered_set<unsigned int> used_k;
422  for (unsigned int j = 0; j < colors_counter[i].size(); ++j)
423  {
424  // Find the color in the current partition with the
425  // largest number of iterators.
426  std::vector<unsigned int>::iterator it;
427  it = std::max_element(colors_counter[i].begin(),
428  colors_counter[i].end());
429  unsigned int min_iterators(static_cast<unsigned int>(-1));
430  unsigned int pos(0);
431  // Find the color of coloring with the least number of
432  // colors among the colors that have not been used yet.
433  for (unsigned int k = 0; k < max_odd_n_colors; ++k)
434  if (used_k.count(k) == 0)
435  if (colors_counter[i_color][k] < min_iterators)
436  {
437  min_iterators = colors_counter[i_color][k];
438  pos = k;
439  }
440  colors_counter[i_color][pos] += *it;
441  // Concatenate the current color with the existing
442  // coloring.
443  coloring[max_even_n_colors + pos].insert(
444  coloring[max_even_n_colors + pos].end(),
445  partition_coloring[i][it - colors_counter[i].begin()]
446  .begin(),
447  partition_coloring[i][it - colors_counter[i].begin()]
448  .end());
449  used_k.insert(pos);
450  // Put the number of iterators to the current color to
451  // zero.
452  *it = 0;
453  }
454  }
455  }
456  }
457 
458  return coloring;
459  }
460  } // namespace internal
461 
462 
542  template <typename Iterator>
543  std::vector<std::vector<Iterator>>
545  const Iterator & begin,
546  const typename identity<Iterator>::type & end,
547  const std::function<std::vector<types::global_dof_index>(
548  const typename identity<Iterator>::type &)> &get_conflict_indices)
549  {
550  Assert(begin != end,
551  ExcMessage(
552  "GraphColoring is not prepared to deal with empty ranges!"));
553 
554  // Create the partitioning.
555  std::vector<std::vector<Iterator>> partitioning =
556  internal::create_partitioning(begin, end, get_conflict_indices);
557 
558  // Color the iterators within each partition.
559  // Run the coloring algorithm on each zone in parallel
560  const unsigned int partitioning_size(partitioning.size());
561  std::vector<std::vector<std::vector<Iterator>>> partition_coloring(
562  partitioning_size);
563 
564  Threads::TaskGroup<> tasks;
565  for (unsigned int i = 0; i < partitioning_size; ++i)
566  tasks += Threads::new_task(&internal::make_dsatur_coloring<Iterator>,
567  partitioning[i],
568  get_conflict_indices,
569  partition_coloring[i]);
570  tasks.join_all();
571 
572  // Gather the colors together.
573  return internal::gather_colors(partition_coloring);
574  }
575 
582  unsigned int
583  color_sparsity_pattern(const SparsityPattern & sparsity_pattern,
584  std::vector<unsigned int> &color_indices);
585 
586 } // namespace GraphColoring
587 
588 DEAL_II_NAMESPACE_CLOSE
589 
590 
591 //---------------------------- graph_coloring.h ---------------------------
592 // end of #ifndef dealii_graph_coloring_h
593 #endif
594 //---------------------------- graph_coloring.h ---------------------------
Task< RT > new_task(const std::function< RT()> &function)
std::vector< std::vector< Iterator > > make_graph_coloring(const Iterator &begin, const typename identity< Iterator >::type &end, const std::function< std::vector< types::global_dof_index >(const typename identity< Iterator >::type &)> &get_conflict_indices)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1407
unsigned int color_sparsity_pattern(const SparsityPattern &sparsity_pattern, std::vector< unsigned int > &color_indices)