Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
graph_coloring.h
1 
2 // ---------------------------------------------------------------------
3 // @f$Id: graph_coloring.h 31932 2013-12-08 02:15:54Z heister @f$
4 //
5 // Copyright (C) 2013 by the deal.II authors
6 //
7 // This file is part of the deal.II library.
8 //
9 // The deal.II library is free software; you can use it, redistribute
10 // it, and/or modify it under the terms of the GNU Lesser General
11 // Public License as published by the Free Software Foundation; either
12 // version 2.1 of the License, or (at your option) any later version.
13 // The full text of the license can be found in the file LICENSE at
14 // the top level of the deal.II distribution.
15 //
16 // ---------------------------------------------------------------------
17 
18 #ifndef __deal2__graph_coloring_h
19 #define __deal2__graph_coloring_h
20 
21 
22 #include <deal.II/base/config.h>
23 #include <deal.II/base/thread_management.h>
24 #include <deal.II/base/std_cxx1x/function.h>
25 #include <deal.II/dofs/dof_handler.h>
26 #include <boost/unordered_map.hpp>
27 #include <boost/unordered_set.hpp>
28 
29 #include <set>
30 #include <vector>
31 
32 
34 
38 namespace GraphColoring
39 {
40  namespace internal
41  {
50  inline
51  bool
52  have_nonempty_intersection (const std::vector<types::global_dof_index> &indices1,
53  const std::vector<types::global_dof_index> &indices2)
54  {
55  // we assume that both arrays are sorted, so we can walk
56  // them in lockstep and see if we encounter an index that's
57  // in both arrays. once we reach the end of either array,
58  // we know that there is no intersection
59  std::vector<types::global_dof_index>::const_iterator
60  p = indices1.begin(),
61  q = indices2.begin();
62  while ((p != indices1.end()) && (q != indices2.end()))
63  {
64  if (*p < *q)
65  ++p;
66  else if (*p > *q)
67  ++q;
68  else
69  // conflict found!
70  return true;
71  }
72 
73  // no conflict found!
74  return false;
75  }
76 
77 
108  template <typename Iterator>
109  std::vector<std::vector<Iterator> >
110  create_partitioning(const Iterator &begin,
111  const typename identity<Iterator>::type &end,
112  const std_cxx1x::function<std::vector<types::global_dof_index> (const Iterator &)> &get_conflict_indices)
113  {
114  // Number of iterators.
115  unsigned int n_iterators = 0;
116 
117  // Create a map from conflict indices to iterators
118  boost::unordered_map<types::global_dof_index,std::vector<Iterator> > indices_to_iterators;
119  for (Iterator it=begin; it!=end; ++it)
120  {
121  const std::vector<types::global_dof_index> conflict_indices = get_conflict_indices(it);
122  const unsigned int n_conflict_indices = conflict_indices.size();
123  for (unsigned int i=0; i<n_conflict_indices; ++i)
124  indices_to_iterators[conflict_indices[i]].push_back(it);
125  ++n_iterators;
126  }
127 
128  // create the very first zone which contains only the first
129  // iterator. then create the other zones. keep track of all the
130  // iterators that have already been assigned to a zone
131  std::vector<std::vector<Iterator> > zones(1,std::vector<Iterator> (1,begin));
132  std::set<Iterator> used_it;
133  used_it.insert(begin);
134  while (used_it.size()!=n_iterators)
135  {
136  // loop over the elements of the previous zone. for each element of
137  // the previous zone, get the conflict indices and from there get
138  // those iterators that are conflicting with the current element
139  typename std::vector<Iterator>::iterator previous_zone_it(zones.back().begin());
140  typename std::vector<Iterator>::iterator previous_zone_end(zones.back().end());
141  std::vector<Iterator> new_zone;
142  for (; previous_zone_it!=previous_zone_end; ++previous_zone_it)
143  {
144  const std::vector<types::global_dof_index>
145  conflict_indices = get_conflict_indices(*previous_zone_it);
146 
147  const unsigned int n_conflict_indices(conflict_indices.size());
148  for (unsigned int i=0; i<n_conflict_indices; ++i)
149  {
150  const std::vector<Iterator> &conflicting_elements
151  = indices_to_iterators[conflict_indices[i]];
152  for (unsigned int j=0; j<conflicting_elements.size(); ++j)
153  {
154  // check that the iterator conflicting with the current one is not
155  // associated to a zone yet and if so, assign it to the current
156  // zone. mark it as used
157  //
158  // we can shortcut this test if the conflicting iterator is the
159  // current iterator
160  if ((conflicting_elements[j] != *previous_zone_it)
161  &&
162  (used_it.count(conflicting_elements[j])==0))
163  {
164  new_zone.push_back(conflicting_elements[j]);
165  used_it.insert(conflicting_elements[j]);
166  }
167  }
168  }
169  }
170 
171  // If there are iterators in the new zone, then the zone is added to the
172  // partition. Otherwise, the graph is disconnected and we need to find
173  // an iterator on the other part of the graph. start the whole process again
174  // with the first iterator that hasn't been assigned to a zone yet
175  if (new_zone.size()!=0)
176  zones.push_back(new_zone);
177  else
178  for (Iterator it=begin; it!=end; ++it)
179  if (used_it.count(it)==0)
180  {
181  zones.push_back(std::vector<Iterator> (1,it));
182  used_it.insert(it);
183  break;
184  }
185  }
186 
187  return zones;
188  }
189 
190 
191 
212  template <typename Iterator>
213  void
214  make_dsatur_coloring(std::vector<Iterator> &partition,
215  const std_cxx1x::function<std::vector<types::global_dof_index> (const Iterator &)> &get_conflict_indices,
216  std::vector<std::vector<Iterator> > &partition_coloring)
217  {
218  partition_coloring.clear ();
219 
220  // Number of zones composing the partitioning.
221  const unsigned int partition_size(partition.size());
222  std::vector<unsigned int> sorted_vertices(partition_size);
223  std::vector<int> degrees(partition_size);
224  std::vector<std::vector<types::global_dof_index> > conflict_indices(partition_size);
225  std::vector<std::vector<unsigned int> > graph(partition_size);
226 
227  // Get the conflict indices associated to each iterator. The conflict_indices have to
228  // be sorted so we can more easily find conflicts later on
229  for (unsigned int i=0; i<partition_size; ++i)
230  {
231  conflict_indices[i] = get_conflict_indices(partition[i]);
232  std::sort(conflict_indices[i].begin(), conflict_indices[i].end());
233  }
234 
235  // Compute the degree of each vertex of the graph using the
236  // intersection of the conflict indices.
237  for (unsigned int i=0; i<partition_size; ++i)
238  for (unsigned int j=i+1; j<partition_size; ++j)
239  // If the two iterators share indices then we increase the degree of the
240  // vertices and create an ''edge'' in the graph.
241  if (have_nonempty_intersection (conflict_indices[i], conflict_indices[j]))
242  {
243  ++degrees[i];
244  ++degrees[j];
245  graph[i].push_back(j);
246  graph[j].push_back(i);
247  }
248 
249  // Sort the vertices by decreasing degree.
250  std::vector<int>::iterator degrees_it;
251  for (unsigned int i=0; i<partition_size; ++i)
252  {
253  // Find the largest element.
254  degrees_it = std::max_element(degrees.begin(),degrees.end());
255  sorted_vertices[i] = degrees_it-degrees.begin();
256  // Put the largest element to -1 so it cannot be chosen again.
257  *degrees_it = -1;
258  }
259 
260  // Color the graph.
261  std::vector<boost::unordered_set<unsigned int> > colors_used;
262  for (unsigned int i=0; i<partition_size; ++i)
263  {
264  const unsigned int current_vertex(sorted_vertices[i]);
265  bool new_color(true);
266  // Try to use an existing color, i.e., try to find a color which is not
267  // associated to one of the vertices linked to current_vertex.
268  // Loop over the color.
269  for (unsigned int j=0; j<partition_coloring.size(); ++j)
270  {
271  // Loop on the vertices linked to current_vertex. If one vertex linked
272  // to current_vertex is already using the color j, this color cannot
273  // be used anymore.
274  bool unused_color(true);
275  for (unsigned int k=0; k<graph[current_vertex].size(); ++k)
276  if (colors_used[j].count(graph[current_vertex][k])==1)
277  {
278  unused_color = false;
279  break;
280  }
281  if (unused_color)
282  {
283  partition_coloring[j].push_back(partition[current_vertex]);
284  colors_used[j].insert(current_vertex);
285  new_color = false;
286  break;
287  }
288  }
289  // Add a new color.
290  if (new_color)
291  {
292  partition_coloring.push_back(std::vector<Iterator> (1,
293  partition[current_vertex]));
294  boost::unordered_set<unsigned int> tmp;
295  tmp.insert(current_vertex);
296  colors_used.push_back(tmp);
297  }
298  }
299  }
300 
301 
302 
311  template <typename Iterator>
312  std::vector<std::vector<Iterator> >
313  gather_colors(const std::vector<std::vector<std::vector<Iterator> > > &partition_coloring)
314  {
315  std::vector<std::vector<Iterator> > coloring;
316 
317  // Count the number of iterators in each color.
318  const unsigned int partition_size(partition_coloring.size());
319  std::vector<std::vector<unsigned int> > colors_counter(partition_size);
320  for (unsigned int i=0; i<partition_size; ++i)
321  {
322  const unsigned int n_colors(partition_coloring[i].size());
323  colors_counter[i].resize(n_colors);
324  for (unsigned int j=0; j<n_colors; ++j)
325  colors_counter[i][j] = partition_coloring[i][j].size();
326  }
327 
328  // Find the partition with the largest number of colors for the even partition.
329  unsigned int i_color(0);
330  unsigned int max_even_n_colors(0);
331  const unsigned int colors_size(colors_counter.size());
332  for (unsigned int i=0; i<colors_size; i+=2)
333  {
334  if (max_even_n_colors<colors_counter[i].size())
335  {
336  max_even_n_colors = colors_counter[i].size();
337  i_color = i;
338  }
339  }
340  coloring.resize(max_even_n_colors);
341  for (unsigned int j=0; j<colors_counter[i_color].size(); ++j)
342  coloring[j] = partition_coloring[i_color][j];
343 
344  for (unsigned int i=0; i<partition_size; i+=2)
345  {
346  if (i!=i_color)
347  {
348  boost::unordered_set<unsigned int> used_k;
349  for (unsigned int j=0; j<colors_counter[i].size(); ++j)
350  {
351  // Find the color in the current partition with the largest number of
352  // iterators.
353  std::vector<unsigned int>::iterator it;
354  it = std::max_element(colors_counter[i].begin(),colors_counter[i].end());
355  unsigned int min_iterators(-1);
356  unsigned int pos(0);
357  // Find the color of coloring with the least number of colors among
358  // the colors that have not been used yet.
359  for (unsigned int k=0; k<max_even_n_colors; ++k)
360  if (used_k.count(k)==0)
361  if (colors_counter[i_color][k]<min_iterators)
362  {
363  min_iterators = colors_counter[i_color][k];
364  pos = k;
365  }
366  colors_counter[i_color][pos] += *it;
367  // Concatenate the current color with the existing coloring.
368  coloring[pos].insert(coloring[pos].end(),
369  partition_coloring[i][it-colors_counter[i].begin()].begin(),
370  partition_coloring[i][it-colors_counter[i].begin()].end());
371  used_k.insert(pos);
372  // Put the number of iterators to the current color to zero.
373  *it = 0;
374  }
375  }
376  }
377 
378  // If there is more than one partition, do the same thing that we did for the even partitions
379  // to the odd partitions
380  if (partition_size>1)
381  {
382  unsigned int max_odd_n_colors(0);
383  for (unsigned int i=1; i<partition_size; i+=2)
384  {
385  if (max_odd_n_colors<colors_counter[i].size())
386  {
387  max_odd_n_colors = colors_counter[i].size();
388  i_color = i;
389  }
390  }
391  coloring.resize(max_even_n_colors+max_odd_n_colors);
392  for (unsigned int j=0; j<colors_counter[i_color].size(); ++j)
393  coloring[max_even_n_colors+j] = partition_coloring[i_color][j];
394 
395  for (unsigned int i=1; i<partition_size; i+=2)
396  {
397  if (i!=i_color)
398  {
399  boost::unordered_set<unsigned int> used_k;
400  for (unsigned int j=0; j<colors_counter[i].size(); ++j)
401  {
402  // Find the color in the current partition with the largest number of
403  // iterators.
404  std::vector<unsigned int>::iterator it;
405  it = std::max_element(colors_counter[i].begin(),colors_counter[i].end());
406  unsigned int min_iterators(-1);
407  unsigned int pos(0);
408  // Find the color of coloring with the least number of colors among
409  // the colors that have not been used yet.
410  for (unsigned int k=0; k<max_odd_n_colors; ++k)
411  if (used_k.count(k)==0)
412  if (colors_counter[i_color][k]<min_iterators)
413  {
414  min_iterators = colors_counter[i_color][k];
415  pos = k;
416  }
417  colors_counter[i_color][pos] += *it;
418  // Concatenate the current color with the existing coloring.
419  coloring[max_even_n_colors+pos].insert(coloring[max_even_n_colors+pos].end(),
420  partition_coloring[i][it-colors_counter[i].begin()].begin(),
421  partition_coloring[i][it-colors_counter[i].begin()].end());
422  used_k.insert(pos);
423  // Put the number of iterators to the current color to zero.
424  *it = 0;
425  }
426  }
427  }
428  }
429 
430  return coloring;
431  }
432  }
433 
434 
511  template <typename Iterator>
512  std::vector<std::vector<Iterator> >
513  make_graph_coloring(const Iterator &begin,
514  const typename identity<Iterator>::type &end,
515  const std_cxx1x::function<std::vector<types::global_dof_index> (const typename identity<Iterator>::type &)> &get_conflict_indices)
516  {
517  Assert (begin != end, ExcMessage ("GraphColoring is not prepared to deal with empty ranges!"));
518 
519  // Create the partitioning.
520  std::vector<std::vector<Iterator> >
521  partitioning = internal::create_partitioning (begin,
522  end,
523  get_conflict_indices);
524 
525  // Color the iterators within each partition.
526  // Run the coloring algorithm on each zone in parallel
527  const unsigned int partitioning_size(partitioning.size());
528  std::vector<std::vector<std::vector<Iterator> > >
529  partition_coloring(partitioning_size);
530 
531  Threads::TaskGroup<> tasks;
532  for (unsigned int i=0; i<partitioning_size; ++i)
533  tasks += Threads::new_task (&internal::make_dsatur_coloring<Iterator>,
534  partitioning[i],
535  get_conflict_indices,
536  partition_coloring[i]);
537  tasks.join_all();
538 
539  // Gather the colors together.
540  return internal::gather_colors(partition_coloring);
541  }
542 
543 } // End graph_coloring namespace
544 
545 DEAL_II_NAMESPACE_CLOSE
546 
547 
548 //---------------------------- graph_coloring.h ---------------------------
549 // end of #ifndef __deal2__graph_coloring_h
550 #endif
551 //---------------------------- graph_coloring.h ---------------------------
::ExceptionBase & ExcMessage(std::string arg1)
std::vector< std::vector< Iterator > > make_graph_coloring(const Iterator &begin, const typename identity< Iterator >::type &end, const std_cxx1x::function< std::vector< types::global_dof_index >(const typename identity< Iterator >::type &)> &get_conflict_indices)
#define Assert(cond, exc)
Definition: exceptions.h:299
Task< RT > new_task(const std_cxx1x::function< RT()> &function)