EpetraExt Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
EpetraExt_MapColoring.cpp
Go to the documentation of this file.
1//@HEADER
2// ***********************************************************************
3//
4// EpetraExt: Epetra Extended - Linear Algebra Services Package
5// Copyright (2011) Sandia Corporation
6//
7// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8// the U.S. Government retains certain rights in this software.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38//
39// ***********************************************************************
40//@HEADER
41
42#include <Epetra_ConfigDefs.h>
44
47
48#include <Epetra_CrsGraph.h>
50#include <Epetra_MapColoring.h>
51#include <Epetra_Map.h>
52#include <Epetra_Comm.h>
53#include <Epetra_Util.h>
54#include <Epetra_Import.h>
55#include <Epetra_Export.h>
56
57#include <Epetra_Time.h>
58
59#include <vector>
60#include <set>
61#include <map>
62#include <random>
63
64using std::vector;
65using std::set;
66using std::map;
67
68namespace EpetraExt {
69
70template<typename int_type>
74{
75 Epetra_Time timer( orig.Comm() );
76
77 origObj_ = &orig;
78
79 const Epetra_BlockMap & RowMap = orig.RowMap();
80 int nRows = RowMap.NumMyElements();
81 const Epetra_BlockMap & ColMap = orig.ColMap();
82 int nCols = ColMap.NumMyElements();
83
84 int MyPID = RowMap.Comm().MyPID();
85
86 if( verbosity_ > 1 ) std::cout << "RowMap:\n" << RowMap;
87 if( verbosity_ > 1 ) std::cout << "ColMap:\n" << ColMap;
88
89 Epetra_MapColoring * ColorMap = 0;
90
91 if( algo_ == GREEDY || algo_ == LUBY )
92 {
93
94 Epetra_CrsGraph * base = &( const_cast<Epetra_CrsGraph&>(orig) );
95
96 int NumIndices;
97 int * Indices;
98
99 double wTime1 = timer.WallTime();
100
101 // For parallel applications, add in boundaries to coloring
102 bool distributedGraph = RowMap.DistributedGlobal();
103 if( distributedGraph )
104 {
105 base = new Epetra_CrsGraph( Copy, ColMap, ColMap, 0 );
106
107 for( int i = 0; i < nRows; ++i )
108 {
109 orig.ExtractMyRowView( i, NumIndices, Indices );
110 base->InsertMyIndices( i, NumIndices, Indices );
111
112 //Do this with a single insert
113 //Is this the right thing?
114 for( int j = 0; j < NumIndices; ++j )
115 if( Indices[j] >= nRows )
116 base->InsertMyIndices( Indices[j], 1, &i );
117 }
118 base->FillComplete();
119 }
120
121 if( verbosity_ > 1 ) std::cout << "Base Graph:\n" << *base << std::endl;
122
123 double wTime2 = timer.WallTime();
124 if( verbosity_ > 0 )
125 std::cout << "EpetraExt::MapColoring [INSERT BOUNDARIES] Time: " << wTime2-wTime1 << std::endl;
126
127 //Generate Local Distance-1 Adjacency Graph
128 //(Transpose of orig excluding non-local column indices)
129 EpetraExt::CrsGraph_Transpose transposeTransform( true );
130 Epetra_CrsGraph & Adj1 = transposeTransform( *base );
131 if( verbosity_ > 1 ) std::cout << "Adjacency 1 Graph!\n" << Adj1;
132
133 wTime1 = timer.WallTime();
134 if( verbosity_ > 0 )
135 std::cout << "EpetraExt::MapColoring [TRANSPOSE GRAPH] Time: " << wTime1-wTime2 << std::endl;
136
137 int Delta = Adj1.MaxNumIndices();
138 if( verbosity_ > 0 ) std::cout << std::endl << "Delta: " << Delta << std::endl;
139
140 //Generation of Local Distance-2 Adjacency Graph
141 Epetra_CrsGraph * Adj2 = &Adj1;
142 if( !distance1_ )
143 {
144 Adj2 = new Epetra_CrsGraph( Copy, ColMap, ColMap, 0 );
145 int NumAdj1Indices;
146 int * Adj1Indices;
147 for( int i = 0; i < nCols; ++i )
148 {
149 Adj1.ExtractMyRowView( i, NumAdj1Indices, Adj1Indices );
150
151 set<int> Cols;
152 for( int j = 0; j < NumAdj1Indices; ++j )
153 {
154 base->ExtractMyRowView( Adj1Indices[j], NumIndices, Indices );
155 for( int k = 0; k < NumIndices; ++k )
156 if( Indices[k] < nCols ) Cols.insert( Indices[k] );
157 }
158 int nCols2 = Cols.size();
159 std::vector<int> ColVec( nCols2 );
160 set<int>::iterator iterIS = Cols.begin();
161 set<int>::iterator iendIS = Cols.end();
162 for( int j = 0 ; iterIS != iendIS; ++iterIS, ++j ) ColVec[j] = *iterIS;
163 Adj2->InsertMyIndices( i, nCols2, &ColVec[0] );
164 }
165 Adj2->FillComplete();
166
167 if( verbosity_ > 1 ) std::cout << "Adjacency 2 Graph!\n" << *Adj2;
168
169 wTime2 = timer.WallTime();
170 if( verbosity_ > 0 )
171 std::cout << "EpetraExt::MapColoring [GEN DIST-2 GRAPH] Time: " << wTime2-wTime1 << std::endl;
172 }
173
174 wTime2 = timer.WallTime();
175
176 ColorMap = new Epetra_MapColoring( ColMap );
177
178 std::vector<int> rowOrder( nCols );
179 if( reordering_ == 0 || reordering_ == 1 )
180 {
181 std::multimap<int,int> adjMap;
182 typedef std::multimap<int,int>::value_type adjMapValueType;
183 for( int i = 0; i < nCols; ++i )
184 adjMap.insert( adjMapValueType( Adj2->NumMyIndices(i), i ) );
185 std::multimap<int,int>::iterator iter = adjMap.begin();
186 std::multimap<int,int>::iterator end = adjMap.end();
187 if( reordering_ == 0 ) //largest first (less colors)
188 {
189 for( int i = 1; iter != end; ++iter, ++i )
190 rowOrder[ nCols - i ] = (*iter).second;
191 }
192 else //smallest first (better balance)
193 {
194 for( int i = 0; iter != end; ++iter, ++i )
195 rowOrder[ i ] = (*iter).second;
196 }
197 }
198 else if( reordering_ == 2 ) //random
199 {
200 for( int i = 0; i < nCols; ++i )
201 rowOrder[ i ] = i;
202#ifndef TFLOP
203 std::random_device rd;
204 std::mt19937 g(rd());
205 std::shuffle( rowOrder.begin(), rowOrder.end(), g );
206#endif
207 }
208
209 wTime1 = timer.WallTime();
210 if( verbosity_ > 0 )
211 std::cout << "EpetraExt::MapColoring [REORDERING] Time: " << wTime1-wTime2 << std::endl;
212
213 //Application of Greedy Algorithm to generate Color Map
214 if( algo_ == GREEDY )
215 {
216 for( int col = 0; col < nCols; ++col )
217 {
218 Adj2->ExtractMyRowView( rowOrder[col], NumIndices, Indices );
219
220 set<int> usedColors;
221 int color;
222 for( int i = 0; i < NumIndices; ++i )
223 {
224 color = (*ColorMap)[ Indices[i] ];
225 if( color > 0 ) usedColors.insert( color );
226 color = 0;
227 int testcolor = 1;
228 while( !color )
229 {
230 if( !usedColors.count( testcolor ) ) color = testcolor;
231 ++testcolor;
232 }
233 }
234 (*ColorMap)[ rowOrder[col] ] = color;
235 }
236
237 wTime2 = timer.WallTime();
238 if( verbosity_ > 0 )
239 std::cout << "EpetraExt::MapColoring [GREEDY COLORING] Time: " << wTime2-wTime1 << std::endl;
240 if( verbosity_ > 0 )
241 std::cout << "Num GREEDY Colors: " << ColorMap->NumColors() << std::endl;
242 }
243 else if( algo_ == LUBY )
244 {
245 //Assign Random Keys To Rows
246 Epetra_Util util;
247 std::vector<int> Keys(nCols);
248 std::vector<int> State(nCols,-1);
249
250 for( int col = 0; col < nCols; ++col )
251 Keys[col] = util.RandomInt();
252
253 int NumRemaining = nCols;
254 int CurrentColor = 1;
255
256 while( NumRemaining > 0 )
257 {
258 //maximal independent set
259 while( NumRemaining > 0 )
260 {
261 NumRemaining = 0;
262
263 //zero out everyone less than neighbor
264 for( int i = 0; i < nCols; ++i )
265 {
266 int col = rowOrder[i];
267 if( State[col] < 0 )
268 {
269 Adj2->ExtractMyRowView( col, NumIndices, Indices );
270 int MyKey = Keys[col];
271 for( int j = 0; j < NumIndices; ++j )
272 if( col != Indices[j] && State[ Indices[j] ] < 0 )
273 {
274 if( MyKey > Keys[ Indices[j] ] ) State[ Indices[j] ] = 0;
275 else State[ col ] = 0;
276 }
277 }
278 }
279
280 //assign -1's to current color
281 for( int col = 0; col < nCols; ++col )
282 {
283 if( State[col] < 0 )
284 State[col] = CurrentColor;
285 }
286
287 //reinstate any zero not neighboring current color
288 for( int col = 0; col < nCols; ++col )
289 if( State[col] == 0 )
290 {
291 Adj2->ExtractMyRowView( col, NumIndices, Indices );
292 bool flag = false;
293 for( int i = 0; i < NumIndices; ++i )
294 if( col != Indices[i] && State[ Indices[i] ] == CurrentColor )
295 {
296 flag = true;
297 break;
298 }
299 if( !flag )
300 {
301 State[col] = -1;
302 ++NumRemaining;
303 }
304 }
305 }
306
307 //Reset Status for all non-colored nodes
308 for( int col = 0; col < nCols; ++col )
309 if( State[col] == 0 )
310 {
311 State[col] = -1;
312 ++NumRemaining;
313 }
314
315 if( verbosity_ > 2 )
316 {
317 std::cout << "Finished Color: " << CurrentColor << std::endl;
318 std::cout << "NumRemaining: " << NumRemaining << std::endl;
319 }
320
321 //New color
322 ++CurrentColor;
323 }
324
325 for( int col = 0; col < nCols; ++col )
326 (*ColorMap)[col] = State[col]-1;
327
328 wTime2 = timer.WallTime();
329 if( verbosity_ > 0 )
330 std::cout << "EpetraExt::MapColoring [LUBI COLORING] Time: " << wTime2-wTime1 << std::endl;
331 if( verbosity_ > 0 )
332 std::cout << "Num LUBI Colors: " << ColorMap->NumColors() << std::endl;
333
334 }
335 else
336 abort(); //UNKNOWN ALGORITHM
337
338 if( distributedGraph ) delete base;
339 if( !distance1_ ) delete Adj2;
340 }
341 else
342 {
343 //Generate Parallel Adjacency-2 Graph
344 EpetraExt::CrsGraph_Overlap OverlapTrans(1);
345 Epetra_CrsGraph & OverlapGraph = OverlapTrans( orig );
346
347 if( verbosity_ > 1 ) std::cout << "OverlapGraph:\n" << OverlapGraph;
348
349 Epetra_CrsGraph * Adj2 = &orig;
350
351 int NumIndices;
352 int * Indices;
353 if( !distance1_ )
354 {
355 Adj2 = new Epetra_CrsGraph( Copy, RowMap, OverlapGraph.ColMap(), 0 );
356 int NumAdj1Indices;
357 int * Adj1Indices;
358 for( int i = 0; i < nRows; ++i )
359 {
360 OverlapGraph.ExtractMyRowView( i, NumAdj1Indices, Adj1Indices );
361
362 set<int> Cols;
363 for( int j = 0; j < NumAdj1Indices; ++j )
364 {
365 int GID = OverlapGraph.LRID( OverlapGraph.GCID64( Adj1Indices[j] ) );
366 OverlapGraph.ExtractMyRowView( GID, NumIndices, Indices );
367 for( int k = 0; k < NumIndices; ++k ) Cols.insert( Indices[k] );
368 }
369 int nCols2 = Cols.size();
370 std::vector<int> ColVec( nCols2 );
371 set<int>::iterator iterIS = Cols.begin();
372 set<int>::iterator iendIS = Cols.end();
373 for( int j = 0 ; iterIS != iendIS; ++iterIS, ++j ) ColVec[j] = *iterIS;
374 Adj2->InsertMyIndices( i, nCols2, &ColVec[0] );
375 }
376
377#ifdef NDEBUG
378 (void) Adj2->FillComplete();
379#else
380 // assert() statements go away if NDEBUG is defined. Don't
381 // declare the 'flag' variable if it never gets used.
382 int flag = Adj2->FillComplete();
383 assert( flag == 0 );
384#endif // NDEBUG
385
386 RowMap.Comm().Barrier();
387 if( verbosity_ > 1 ) std::cout << "Adjacency 2 Graph!\n" << *Adj2;
388 }
389
390 //collect GIDs on boundary
391 std::vector<int_type> boundaryGIDs;
392 std::vector<int_type> interiorGIDs;
393 for( int row = 0; row < nRows; ++row )
394 {
395 Adj2->ExtractMyRowView( row, NumIndices, Indices );
396 bool testFlag = false;
397 for( int i = 0; i < NumIndices; ++i )
398 if( Indices[i] >= nRows ) testFlag = true;
399 if( testFlag ) boundaryGIDs.push_back( (int_type) Adj2->GRID64(row) );
400 else interiorGIDs.push_back( (int_type) Adj2->GRID64(row) );
401 }
402
403 int LocalBoundarySize = boundaryGIDs.size();
404
405 Epetra_Map BoundaryMap( -1, boundaryGIDs.size(),
406 LocalBoundarySize ? &boundaryGIDs[0]: 0,
407 0, RowMap.Comm() );
408 if( verbosity_ > 1 ) std::cout << "BoundaryMap:\n" << BoundaryMap;
409
410 int_type BoundarySize = (int_type) BoundaryMap.NumGlobalElements64();
411 Epetra_MapColoring BoundaryColoring( BoundaryMap );
412
413 if( algo_ == PSEUDO_PARALLEL )
414 {
415 Epetra_Map BoundaryIndexMap( BoundarySize, LocalBoundarySize, 0, RowMap.Comm() );
416 if( verbosity_ > 1) std::cout << "BoundaryIndexMap:\n" << BoundaryIndexMap;
417
418 typename Epetra_GIDTypeVector<int_type>::impl bGIDs( View, BoundaryIndexMap, &boundaryGIDs[0] );
419 if( verbosity_ > 1) std::cout << "BoundaryGIDs:\n" << bGIDs;
420
421 int_type NumLocalBs = 0;
422 if( !RowMap.Comm().MyPID() ) NumLocalBs = BoundarySize;
423
424 Epetra_Map LocalBoundaryIndexMap( BoundarySize, NumLocalBs, 0, RowMap.Comm() );
425 if( verbosity_ > 1) std::cout << "LocalBoundaryIndexMap:\n" << LocalBoundaryIndexMap;
426
427 typename Epetra_GIDTypeVector<int_type>::impl lbGIDs( LocalBoundaryIndexMap );
428 Epetra_Import lbImport( LocalBoundaryIndexMap, BoundaryIndexMap );
429 lbGIDs.Import( bGIDs, lbImport, Insert );
430 if( verbosity_ > 1) std::cout << "LocalBoundaryGIDs:\n" << lbGIDs;
431
432 Epetra_Map LocalBoundaryMap( BoundarySize, NumLocalBs, lbGIDs.Values(), 0, RowMap.Comm() );
433 if( verbosity_ > 1) std::cout << "LocalBoundaryMap:\n" << LocalBoundaryMap;
434
435 Epetra_CrsGraph LocalBoundaryGraph( Copy, LocalBoundaryMap, LocalBoundaryMap, 0 );
436 Epetra_Import LocalBoundaryImport( LocalBoundaryMap, Adj2->RowMap() );
437 LocalBoundaryGraph.Import( *Adj2, LocalBoundaryImport, Insert );
438 LocalBoundaryGraph.FillComplete();
439 if( verbosity_ > 1 ) std::cout << "LocalBoundaryGraph:\n " << LocalBoundaryGraph;
440
442 Epetra_MapColoring & LocalBoundaryColoring = BoundaryTrans( LocalBoundaryGraph );
443 if( verbosity_ > 1 ) std::cout << "LocalBoundaryColoring:\n " << LocalBoundaryColoring;
444
445 Epetra_Export BoundaryExport( LocalBoundaryMap, BoundaryMap );
446 BoundaryColoring.Export( LocalBoundaryColoring, BoundaryExport, Insert );
447 }
448 else if( algo_ == JONES_PLASSMAN )
449 {
450 /* Alternative Distrib. Memory Coloring of Boundary based on JonesPlassman(sic) paper
451 * 1.Random number assignment to all boundary nodes using GID as seed to function
452 * (This allows any processor to compute adj. off proc values with a local computation)
453 * 2.Find all nodes greater than any neighbor off processor, color them.
454 * 3.Send colored node info to neighbors
455 * 4.Constrained color all nodes with all off proc neighbors smaller or colored.
456 * 5.Goto 3
457 */
458
459 std::vector<int_type> OverlapBoundaryGIDs( boundaryGIDs );
460 for( int i = nRows; i < Adj2->ColMap().NumMyElements(); ++i )
461 OverlapBoundaryGIDs.push_back( (int_type) Adj2->ColMap().GID64(i) );
462
463 int_type OverlapBoundarySize = OverlapBoundaryGIDs.size();
464 Epetra_Map BoundaryColMap( (int_type) -1, OverlapBoundarySize,
465 OverlapBoundarySize ? &OverlapBoundaryGIDs[0] : 0,
466 0, RowMap.Comm() );
467
468 Epetra_CrsGraph BoundaryGraph( Copy, BoundaryMap, BoundaryColMap, 0 );
469 Epetra_Import BoundaryImport( BoundaryMap, Adj2->RowMap() );
470 BoundaryGraph.Import( *Adj2, BoundaryImport, Insert );
471 BoundaryGraph.FillComplete();
472 if( verbosity_ > 1) std::cout << "BoundaryGraph:\n" << BoundaryGraph;
473
474 Epetra_Import ReverseOverlapBoundaryImport( BoundaryMap, BoundaryColMap );
475 Epetra_Import OverlapBoundaryImport( BoundaryColMap, BoundaryMap );
476
477 int Colored = 0;
478 int GlobalColored = 0;
479 int Level = 0;
480 Epetra_MapColoring OverlapBoundaryColoring( BoundaryColMap );
481
482 //Setup random integers for boundary nodes
483 Epetra_IntVector BoundaryValues( BoundaryMap );
484 Epetra_Util Util;
485 Util.SetSeed( 47954118 * (MyPID+1) );
486 for( int i=0; i < LocalBoundarySize; ++i )
487 {
488 int val = Util.RandomInt();
489 if( val < 0 ) val *= -1;
490 BoundaryValues[i] = val;
491 }
492
493 //Get Random Values for External Boundary
494 Epetra_IntVector OverlapBoundaryValues( BoundaryColMap );
495 OverlapBoundaryValues.Import( BoundaryValues, OverlapBoundaryImport, Insert );
496
497 while( GlobalColored < BoundarySize )
498 {
499 //Find current "Level" of boundary indices to color
500 int theNumIndices;
501 int * theIndices;
502 std::vector<int> LevelIndices;
503 for( int i = 0; i < LocalBoundarySize; ++i )
504 {
505 if( !OverlapBoundaryColoring[i] )
506 {
507 //int MyVal = PRAND(BoundaryColMap.GID(i));
508 int MyVal = OverlapBoundaryValues[i];
509 BoundaryGraph.ExtractMyRowView( i, theNumIndices, theIndices );
510 bool ColorFlag = true;
511 int Loc = 0;
512 while( Loc<theNumIndices && theIndices[Loc]<LocalBoundarySize ) ++Loc;
513 for( int j = Loc; j < theNumIndices; ++j )
514 if( (OverlapBoundaryValues[theIndices[j]]>MyVal)
515 && !OverlapBoundaryColoring[theIndices[j]] )
516 {
517 ColorFlag = false;
518 break;
519 }
520 if( ColorFlag ) LevelIndices.push_back(i);
521 }
522 }
523
524 if( verbosity_ > 1 )
525 {
526 std::cout << MyPID << " Level Indices: ";
527 int Lsize = (int) LevelIndices.size();
528 for( int i = 0; i < Lsize; ++i ) std::cout << LevelIndices[i] << " ";
529 std::cout << std::endl;
530 }
531
532 //Greedy coloring of current level
533 set<int> levelColors;
534 int Lsize = (int) LevelIndices.size();
535 for( int i = 0; i < Lsize; ++i )
536 {
537 BoundaryGraph.ExtractMyRowView( LevelIndices[i], theNumIndices, theIndices );
538
539 set<int> usedColors;
540 int color;
541 for( int j = 0; j < theNumIndices; ++j )
542 {
543 color = OverlapBoundaryColoring[ theIndices[j] ];
544 if( color > 0 ) usedColors.insert( color );
545 color = 0;
546 int testcolor = 1;
547 while( !color )
548 {
549 if( !usedColors.count( testcolor ) ) color = testcolor;
550 ++testcolor;
551 }
552 }
553 OverlapBoundaryColoring[ LevelIndices[i] ] = color;
554 levelColors.insert( color );
555 }
556
557 if( verbosity_ > 2 ) std::cout << MyPID << " Level: " << Level << " Count: " << LevelIndices.size() << " NumColors: " << levelColors.size() << std::endl;
558
559 if( verbosity_ > 2 ) std::cout << "Current Level Boundary Coloring:\n" << OverlapBoundaryColoring;
560
561 //Update off processor coloring info
562 BoundaryColoring.Import( OverlapBoundaryColoring, ReverseOverlapBoundaryImport, Insert );
563 OverlapBoundaryColoring.Import( BoundaryColoring, OverlapBoundaryImport, Insert );
564 Colored += LevelIndices.size();
565 Level++;
566
567 RowMap.Comm().SumAll( &Colored, &GlobalColored, 1 );
568 if( verbosity_ > 2 ) std::cout << "Num Globally Colored: " << GlobalColored << " from Num Global Boundary Nodes: " << BoundarySize << std::endl;
569 }
570 }
571
572 if( verbosity_ > 1 ) std::cout << "BoundaryColoring:\n " << BoundaryColoring;
573
574 Epetra_MapColoring RowColorMap( RowMap );
575
576 //Add Boundary Colors
577 for( int i = 0; i < LocalBoundarySize; ++i )
578 {
579 int_type GID = (int_type) BoundaryMap.GID64(i);
580 RowColorMap(GID) = BoundaryColoring(GID);
581 }
582
583 Epetra_MapColoring Adj2ColColorMap( Adj2->ColMap() );
584 Epetra_Import Adj2Import( Adj2->ColMap(), RowMap );
585 Adj2ColColorMap.Import( RowColorMap, Adj2Import, Insert );
586
587 if( verbosity_ > 1 ) std::cout << "RowColoringMap:\n " << RowColorMap;
588 if( verbosity_ > 1 ) std::cout << "Adj2ColColorMap:\n " << Adj2ColColorMap;
589
590 std::vector<int> rowOrder( nRows );
591 if( reordering_ == 0 || reordering_ == 1 )
592 {
593 std::multimap<int,int> adjMap;
594 typedef std::multimap<int,int>::value_type adjMapValueType;
595 for( int i = 0; i < nRows; ++i )
596 adjMap.insert( adjMapValueType( Adj2->NumMyIndices(i), i ) );
597 std::multimap<int,int>::iterator iter = adjMap.begin();
598 std::multimap<int,int>::iterator end = adjMap.end();
599 if( reordering_ == 0 ) //largest first (less colors)
600 {
601 for( int i = 1; iter != end; ++iter, ++i )
602 rowOrder[nRows-i] = (*iter).second;
603 }
604 else //smallest first (better balance)
605 {
606 for( int i = 0; iter != end; ++iter, ++i )
607 rowOrder[i] = (*iter).second;
608 }
609 }
610 else if( reordering_ == 2 ) //random
611 {
612 for( int i = 0; i < nRows; ++i )
613 rowOrder[i] = i;
614#ifdef TFLOP
615 std::random_device rd;
616 std::mt19937 g(rd());
617 std::shuffle( rowOrder.begin(), rowOrder.end(), g );
618#endif
619 }
620
621 //Constrained greedy coloring of interior
622 set<int> InteriorColors;
623 for( int row = 0; row < nRows; ++row )
624 {
625 if( !RowColorMap[ rowOrder[row] ] )
626 {
627 Adj2->ExtractMyRowView( rowOrder[row], NumIndices, Indices );
628
629 set<int> usedColors;
630 int color;
631 for( int i = 0; i < NumIndices; ++i )
632 {
633 color = Adj2ColColorMap[ Indices[i] ];
634 if( color > 0 ) usedColors.insert( color );
635 color = 0;
636 int testcolor = 1;
637 while( !color )
638 {
639 if( !usedColors.count( testcolor ) ) color = testcolor;
640 ++testcolor;
641 }
642 }
643 Adj2ColColorMap[ rowOrder[row] ] = color;
644 InteriorColors.insert( color );
645 }
646 }
647 if( verbosity_ > 1 ) std::cout << MyPID << " Num Interior Colors: " << InteriorColors.size() << std::endl;
648 if( verbosity_ > 1 ) std::cout << "RowColorMap after Greedy:\n " << RowColorMap;
649
650 ColorMap = new Epetra_MapColoring( ColMap );
651 Epetra_Import ColImport( ColMap, Adj2->ColMap() );
652 ColorMap->Import( Adj2ColColorMap, ColImport, Insert );
653
654 if( !distance1_ ) delete Adj2;
655 }
656
657 if( verbosity_ > 0 ) std::cout << MyPID << " ColorMap Color Count: " << ColorMap->NumColors() << std::endl;
658 if( verbosity_ > 1 ) std::cout << "ColorMap!\n" << *ColorMap;
659
660 newObj_ = ColorMap;
661
662 return *ColorMap;
663}
664
668{
669#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
670 if(orig.RowMap().GlobalIndicesInt()) {
671 return Toperator<int>(orig);
672 }
673 else
674#endif
675#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
676 if(orig.RowMap().GlobalIndicesLongLong()) {
677 return Toperator<long long>(orig);
678 }
679 else
680#endif
681 throw "CrsGraph_MapColoring::operator(): ERROR, GlobalIndices type unknown.";
682}
683
684} // namespace EpetraExt
Insert
View
Copy
Map Coloring of independent columns in a Graph.
CrsGraph_MapColoring::NewTypeRef operator()(CrsGraph_MapColoring::OriginalTypeRef orig)
Generates the Epetra_MapColoring object from an input Epetra_CrsGraph.
CrsGraph_MapColoring::NewTypeRef Toperator(CrsGraph_MapColoring::OriginalTypeRef orig)
Given an input Epetra_CrsGraph, a "overlapped" Epetra_CrsGraph is generated including rows associated...
Transform to generate the explicit transpose of a Epetra_CrsGraph.
bool DistributedGlobal() const
long long GID64(int LID) const
long long NumGlobalElements64() const
const Epetra_Comm & Comm() const
int NumMyElements() const
virtual int MyPID() const=0
virtual void Barrier() const=0
virtual int SumAll(double *PartialSums, double *GlobalSums, int Count) const=0
int ExtractMyRowView(int LocalRow, int &NumIndices, int *&Indices) const
int MaxNumIndices() const
const Epetra_BlockMap & RowMap() const
long long GCID64(int LCID_in) const
int NumMyIndices(int Row) const
int InsertMyIndices(int LocalRow, int NumIndices, int *Indices)
long long GRID64(int LRID_in) const
const Epetra_BlockMap & ColMap() const
int LRID(int GRID_in) const
int Import(const Epetra_SrcDistObject &A, const Epetra_Import &Importer, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor=0)
int Export(const Epetra_SrcDistObject &A, const Epetra_Import &Importer, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor=0)
int NumColors() const
double WallTime(void) const
int SetSeed(unsigned int Seed_in)
unsigned int RandomInt()
EpetraExt::BlockCrsMatrix: A class for constructing a distributed block matrix.