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