Epetra Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
Epetra_MapColoring.cpp
Go to the documentation of this file.
1
2//@HEADER
3// ************************************************************************
4//
5// Epetra: Linear Algebra Services Package
6// Copyright 2011 Sandia Corporation
7//
8// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9// the U.S. Government retains certain rights in this software.
10//
11// Redistribution and use in source and binary forms, with or without
12// modification, are permitted provided that the following conditions are
13// met:
14//
15// 1. Redistributions of source code must retain the above copyright
16// notice, this list of conditions and the following disclaimer.
17//
18// 2. Redistributions in binary form must reproduce the above copyright
19// notice, this list of conditions and the following disclaimer in the
20// documentation and/or other materials provided with the distribution.
21//
22// 3. Neither the name of the Corporation nor the names of the
23// contributors may be used to endorse or promote products derived from
24// this software without specific prior written permission.
25//
26// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37//
38// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
39//
40// ************************************************************************
41//@HEADER
42
43#include "Epetra_ConfigDefs.h"
44#include "Epetra_HashTable.h"
45#include "Epetra_Comm.h"
46#include "Epetra_Map.h"
47#include "Epetra_MapColoring.h"
48#include "Epetra_Util.h"
49
50//=============================================================================
52 const int defaultColor)
53 : Epetra_DistObject(map, "Epetra::MapColoring"),
54 DefaultColor_(defaultColor),
55 ColorIDs_(0),
56 FirstColor_(0),
57 NumColors_(0),
58 ListOfColors_(0),
59 ColorCount_(0),
60 ElementColors_(0),
61 ColorLists_(0),
62 Allocated_(false),
63 ListsAreGenerated_(false),
64 ListsAreValid_(false)
65{
66 Allocate(elementColors, 1);
67}
68//=============================================================================
70 const int defaultColor)
71 : Epetra_DistObject(map, "Epetra::MapColoring"),
72 DefaultColor_(defaultColor),
73 ColorIDs_(0),
74 FirstColor_(0),
75 NumColors_(0),
76 ListOfColors_(0),
77 ColorCount_(0),
78 ElementColors_(0),
79 ColorLists_(0),
80 Allocated_(false),
81 ListsAreGenerated_(false),
82 ListsAreValid_(false)
83{
85}
86//=============================================================================
88 : Epetra_DistObject(Source),
89 DefaultColor_(Source.DefaultColor_),
90 ColorIDs_(0),
91 FirstColor_(0),
92 NumColors_(0),
93 ListOfColors_(0),
94 ColorCount_(0),
95 ElementColors_(0),
96 ColorLists_(0),
97 Allocated_(false),
98 ListsAreGenerated_(false),
99 ListsAreValid_(false)
100{
101 Allocate(Source.ElementColors_, 1);
102}
103//=========================================================================
105
106
107 if (Allocated_ && Map().NumMyElements()>0) delete [] ElementColors_;
109}
110
111//=========================================================================
113
114
115 if (ListsAreGenerated_) {
116 for (int i=0; i<NumColors_; i++) if (ColorLists_[i]!=0) delete [] ColorLists_[i];
117 delete [] ColorLists_;
118 delete [] ColorCount_;
119 delete [] ListOfColors_;
120 delete ColorIDs_;
121 ListItem * CurItem = FirstColor_;
122 while (CurItem!=0) {
123 ListItem * NextItem = CurItem->NextItem;
124 delete CurItem;
125 CurItem = NextItem;
126 }
127 }
128 ListsAreValid_ = false;
129 return(0);
130}
131
132//=========================================================================
133int Epetra_MapColoring::Allocate(int * elementColors, int Increment)
134{
135
136 if (Allocated_) return(0);
137
138 int NumMyElements = Map().NumMyElements();
139 if (NumMyElements>0) ElementColors_ = new int[NumMyElements];
140 for (int i=0; i< NumMyElements; i++) ElementColors_[i] = elementColors[i*Increment];
141 Allocated_ = true;
142 return(0);
143}
144
145//=========================================================================
147 int NumMyElements = Map().NumMyElements();
148 if (NumMyElements==0) return(0); // Nothing to do
149
150 if (ListsAreValid_) return(0); // Already been here
151
152 if (ListsAreGenerated_) DeleteLists(); // Delete any existing lists
153
154
155 // Scan the ElementColors to determine how many colors we have
156 NumColors_ = 1;
157 FirstColor_ = new ListItem(ElementColors_[0]); // Initialize First color in list
158 for (int i=1; i<NumMyElements; i++) if (!InItemList(ElementColors_[i])) NumColors_++;
159
160 // Create hash table that maps color IDs to the integers 0,...NumColors_
162 ListOfColors_ = new int[NumColors_];
163 ListItem * CurItem = FirstColor_;
164 {for (int i=0; i<NumColors_; i++) {
165 ColorIDs_->Add(CurItem->ItemValue, i); // Create hash table entry
166 ListOfColors_[i] = CurItem->ItemValue; // Put color value in a list of colors
167 CurItem = CurItem->NextItem;
168 }}
169 Epetra_Util util;
170 util.Sort(true, NumColors_, ListOfColors_, 0, 0, 0, 0, 0, 0); // Sort List of colors in ascending order
171 // Count the number of IDs of each color
172 ColorCount_ = new int[NumColors_];
173 {for (int i=0; i<NumColors_; i++) ColorCount_[i] = 0;}
174 {for (int i=0; i<NumMyElements; i++) ColorCount_[ColorIDs_->Get(ElementColors_[i])]++;}
175
176 // Finally build list of IDs grouped by color
177 ColorLists_ = new int *[NumColors_];
178 {for (int i=0; i<NumColors_; i++) ColorLists_[i] = new int[ColorCount_[i]];}
179 {for (int i=0; i<NumColors_; i++) ColorCount_[i] = 0;} // Reset so we can use for counting
180 {for (int i=0; i<NumMyElements; i++) {
181 int j = ColorIDs_->Get(ElementColors_[i]);
182 ColorLists_[j][ColorCount_[j]++] = i;
183 }}
184 ListsAreValid_ = true;
185 ListsAreGenerated_ = true;
186
187 return(0);
188}
189//=========================================================================
190bool Epetra_MapColoring::InItemList(int ColorValue) const {
191 bool ColorFound = false;
192 ListItem * CurColor = 0;
193 ListItem * NextColor = FirstColor_;
194 while (!ColorFound && NextColor!=0) {
195 CurColor = NextColor;
196 NextColor = CurColor->NextItem;
197 if (ColorValue==CurColor->ItemValue) ColorFound = true;
198 }
199
200 if (!ColorFound) CurColor->NextItem = new ListItem(ColorValue);
201
202 return(ColorFound);
203}
204//=========================================================================
207 int arrayIndex = -1;
208 if( ColorIDs_ )
209 arrayIndex = ColorIDs_->Get(Color);
210 if (arrayIndex>-1) return(ColorCount_[arrayIndex]);
211 else return(0);
212}
213//=========================================================================
214int * Epetra_MapColoring::ColorLIDList(int Color) const {
216 int arrayIndex = -1;
217 if( ColorIDs_ )
218 arrayIndex = ColorIDs_->Get(Color);
219 if (arrayIndex>-1) return(ColorLists_[arrayIndex]);
220 else return(0);
221}
222//=========================================================================
223template<typename int_type>
225
227 int arrayIndex = -1;
228 if( ColorIDs_ )
229 arrayIndex = ColorIDs_->Get(Color);
230 int NumElements = 0;
231 int * ColorElementLIDs = 0;
232 int_type * ColorElementGIDs =0;
233 if (arrayIndex>-1) NumElements = ColorCount_[arrayIndex];
234 if (NumElements>0) {
235 ColorElementLIDs = ColorLIDList(Color);
236 ColorElementGIDs = new int_type[NumElements];
237 for (int i=0; i<NumElements; i++) ColorElementGIDs[i] = (int_type) Map().GID64(ColorElementLIDs[i]);
238 }
239 Epetra_Map * map = new Epetra_Map((int_type) -1, NumElements, ColorElementGIDs,
240 (int_type) Map().IndexBase64(), Map().Comm());
241 if (ColorElementGIDs!=0) delete [] ColorElementGIDs;
242 return(map);
243}
244
246#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
247 if(Map().GlobalIndicesInt()) {
248 return TGenerateMap<int>(Color);
249 }
250 else
251#endif
252#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
253 if(Map().GlobalIndicesLongLong()) {
254 return TGenerateMap<long long>(Color);
255 }
256 else
257#endif
258 throw "Epetra_MapColoring::GenerateMap: GlobalIndices type unknown";
259}
260
261//=========================================================================
262template<typename int_type>
264
266 int arrayIndex = -1;
267 if( ColorIDs_ )
268 arrayIndex = ColorIDs_->Get(Color);
269 int NumElements = 0;
270 int * ColorElementLIDs = 0;
271 int * ColorElementSizes = 0;
272 int_type * ColorElementGIDs = 0;
273 if (arrayIndex>-1) NumElements = ColorCount_[arrayIndex];
274 if (NumElements>0) {
275 ColorElementLIDs = ColorLIDList(Color);
276 ColorElementSizes = new int[NumElements];
277 ColorElementGIDs = new int_type[NumElements];
278 for (int i=0; i<NumElements; i++) ColorElementGIDs[i] = (int_type) Map().GID64(ColorElementLIDs[i]);
279 }
280 int * MapElementSizes = Map().ElementSizeList();
281
282 {for (int i=0; i<NumElements; i++)
283 ColorElementSizes[i] = MapElementSizes[ColorElementLIDs[i]];}
284
285 Epetra_BlockMap * map = new Epetra_BlockMap((int_type) -1, NumElements, ColorElementGIDs,
286 ColorElementSizes,
287 (int_type) Map().IndexBase64(), Map().Comm());
288
289 if (ColorElementGIDs!=0) delete [] ColorElementGIDs;
290 if (ColorElementSizes!=0) delete [] ColorElementSizes;
291
292 return(map);
293}
294
296#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
297 if(Map().GlobalIndicesInt()) {
298 return TGenerateBlockMap<int>(Color);
299 }
300 else
301#endif
302#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
303 if(Map().GlobalIndicesLongLong()) {
304 return TGenerateBlockMap<long long>(Color);
305 }
306 else
307#endif
308 throw "Epetra_MapColoring::GenerateBlockMap: GlobalIndices type unknown";
309}
310
311//=========================================================================
312void Epetra_MapColoring::Print(std::ostream& os) const {
313 int MyPID = Map().Comm().MyPID();
314 int NumProc = Map().Comm().NumProc();
315
316 if (MyPID==0) os
317 << std::endl
318 << " *****************************************" << std::endl
319 << " Coloring information arranged map element" << std::endl
320 << " *****************************************" << std::endl
321 << std::endl;
322 for (int iproc=0; iproc < NumProc; iproc++) {
323 if (MyPID==iproc) {
324 int NumMyElements1 =Map(). NumMyElements();
325
326 if (MyPID==0) {
327 os.width(8);
328 os << " MyPID"; os << " ";
329 os.width(12);
330 os << "GID ";
331 os.width(20);
332 os << "Color ";
333 os << std::endl;
334 }
335 for (int i=0; i < NumMyElements1; i++) {
336 os.width(10);
337 os << MyPID; os << " ";
338 os.width(10);
339
340 if(Map().GlobalIndicesInt()) {
341#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
342 int * MyGlobalElements1 = Map().MyGlobalElements();
343 os << MyGlobalElements1[i] << " ";
344#else
345 throw ReportError("Epetra_MapColoring::Print: ERROR, GlobalIndicesInt but no API for it.",-1);
346#endif
347 }
348 else if(Map().GlobalIndicesLongLong())
349 {
350#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
351 long long * MyGlobalElements1 = Map().MyGlobalElements64();
352 os << MyGlobalElements1[i] << " ";
353#else
354 throw ReportError("Epetra_MapColoring::Print: ERROR, GlobalIndicesLongLong but no API for it.",-1);
355#endif
356 }
357 else
358 throw ReportError("Epetra_MapColoring::Print: ERROR, Don't know map global index type.",-1);
359
360 os.width(20);
361 os << ElementColors_[i];
362 os << std::endl;
363 }
364 os << std::flush;
365 }
366
367 // Do a few global ops to give I/O a chance to complete
368 Map().Comm().Barrier();
369 Map().Comm().Barrier();
370 Map().Comm().Barrier();
371 }
372
373 if (MyPID==0) os
374 << std::endl
375 << " **************************************" << std::endl
376 << " Coloring information arranged by color" << std::endl
377 << " **************************************" << std::endl
378 << std::endl;
379 {for (int iproc=0; iproc < NumProc; iproc++) {
380 if (MyPID==iproc) {
381 if (NumColors()==0) os << " No colored elements on processor " << MyPID << std::endl;
382 else {
383 os << "Number of colors in map = " << NumColors() << std::endl
384 << "Default color = " << DefaultColor() << std::endl << std::endl;
385 if (MyPID==0) {
386 os.width(8);
387 os << " MyPID"; os << " ";
388 os.width(12);
389 os << "LID ";
390 os.width(20);
391 os << "Color ";
392 os << std::endl;
393 }
394 int * ColorValues = ListOfColors();
395 for (int ii=0; ii<NumColors(); ii++) {
396 int CV = ColorValues[ii];
397 int ColorCount = NumElementsWithColor(CV);
398 int * LIDList = ColorLIDList(CV);
399
400
401 for (int i=0; i < ColorCount; i++) {
402 os.width(10);
403 os << MyPID; os << " ";
404 os.width(10);
405 os << LIDList[i] << " ";
406 os.width(20);
407 os << CV;
408 os << std::endl;
409 }
410 os << std::flush;
411 }
412 }
413 }
414 // Do a few global ops to give I/O a chance to complete
415 Map().Comm().Barrier();
416 Map().Comm().Barrier();
417 Map().Comm().Barrier();
418 }}
419 return;
420}
421//=========================================================================
423
425 int tmp1 = NumColors_, tmp2;
426 Map().Comm().MaxAll(&tmp1, &tmp2, 1);
427 return(tmp2);
428}
429//=========================================================================
431 (void)Source;//prevents unused variable compiler warning
432 return(0);
433}
434
435//=========================================================================
437 int NumSameIDs,
438 int NumPermuteIDs,
439 int * PermuteToLIDs,
440 int *PermuteFromLIDs,
441 const Epetra_OffsetIndex * Indexor,
442 Epetra_CombineMode CombineMode)
443{
444 (void)Indexor;
445 const Epetra_MapColoring & A = dynamic_cast<const Epetra_MapColoring &>(Source);
446
447 int * From = A.ElementColors();
448 int *To = ElementColors_;
449
450 // Do copy first
451 if (NumSameIDs>0)
452 if (To!=From) {
453 if (CombineMode==Epetra_AddLocalAlso)
454 for (int j=0; j<NumSameIDs; j++) To[j] += From[j]; // Add to existing value
455 else
456 for (int j=0; j<NumSameIDs; j++) To[j] = From[j];
457 }
458 // Do local permutation next
459 if (NumPermuteIDs>0) {
460 if (CombineMode==Epetra_AddLocalAlso)
461 for (int j=0; j<NumPermuteIDs; j++) To[PermuteToLIDs[j]] += From[PermuteFromLIDs[j]]; // Add to existing value
462 else
463 for (int j=0; j<NumPermuteIDs; j++) To[PermuteToLIDs[j]] = From[PermuteFromLIDs[j]];
464 }
465
466 return(0);
467}
468
469//=========================================================================
471 int NumExportIDs,
472 int * ExportLIDs,
473 int & LenExports,
474 char * & Exports,
475 int & SizeOfPacket,
476 int * Sizes,
477 bool & VarSizes,
478 Epetra_Distributor & Distor)
479{
480 (void)Sizes;
481 (void)VarSizes;
482 (void)Distor;
483 const Epetra_MapColoring & A = dynamic_cast<const Epetra_MapColoring &>(Source);
484
485 int * From = A.ElementColors();
486 int * IntExports = 0;
487
488 SizeOfPacket = (int)sizeof(int);
489
490 if (NumExportIDs*SizeOfPacket>LenExports) {
491 if (LenExports>0) delete [] Exports;
492 LenExports = NumExportIDs*SizeOfPacket;
493 IntExports = new int[LenExports];
494 Exports = (char *) IntExports;
495 }
496
497 int * ptr;
498
499 if (NumExportIDs>0) {
500 ptr = (int *) Exports;
501 for (int j=0; j<NumExportIDs; j++) ptr[j] = From[ExportLIDs[j]];
502 }
503
504 return(0);
505}
506
507//=========================================================================
509 int NumImportIDs,
510 int * ImportLIDs,
511 int LenImports,
512 char * Imports,
513 int & SizeOfPacket,
514 Epetra_Distributor & Distor,
515 Epetra_CombineMode CombineMode,
516 const Epetra_OffsetIndex * Indexor )
517{
518 (void)Source;
519 (void)LenImports;
520 (void)Imports;
521 (void)SizeOfPacket;
522 (void)Distor;
523 (void)Indexor;
524 int j;
525
526 if( CombineMode != Add
527 && CombineMode != Zero
528 && CombineMode != Insert
529 && CombineMode != AbsMax )
530 EPETRA_CHK_ERR(-1); //Unsupported CombinedMode, will default to Zero
531
532 if (NumImportIDs<=0) return(0);
533
534 int * To = ElementColors_;
535
536 int * ptr;
537 // Unpack it...
538
539 ptr = (int *) Imports;
540
541 if (CombineMode==Add)
542 for (j=0; j<NumImportIDs; j++) To[ImportLIDs[j]] += ptr[j]; // Add to existing value
543 else if(CombineMode==Insert)
544 for (j=0; j<NumImportIDs; j++) To[ImportLIDs[j]] = ptr[j];
545 else if(CombineMode==AbsMax) {
546 for (j=0; j<NumImportIDs; j++) To[ImportLIDs[j]] = 0;
547 for (j=0; j<NumImportIDs; j++) To[ImportLIDs[j]] = EPETRA_MAX( To[ImportLIDs[j]],std::abs(ptr[j]));
548 }
549
550 return(0);
551}
Epetra_CombineMode
@ Epetra_AddLocalAlso
#define EPETRA_MAX(x, y)
#define EPETRA_CHK_ERR(a)
Epetra_BlockMap: A class for partitioning block element vectors and matrices.
int MyGlobalElements(int *MyGlobalElementList) const
Puts list of global elements on this processor into the user-provided array.
long long * MyGlobalElements64() const
int * ElementSizeList() const
List of the element sizes corresponding to the array MyGlobalElements().
long long GID64(int LID) const
const Epetra_Comm & Comm() const
Access function for Epetra_Comm communicator.
int NumMyElements() const
Number of elements on the calling processor.
virtual int MaxAll(double *PartialMaxs, double *GlobalMaxs, int Count) const =0
Epetra_Comm Global Max function.
virtual int NumProc() const =0
Returns total number of processes.
virtual int MyPID() const =0
Return my process ID.
virtual void Barrier() const =0
Epetra_Comm Barrier function.
Epetra_DistObject: A class for constructing and using dense multi-vectors, vectors and matrices in pa...
const Epetra_BlockMap & Map() const
Returns the address of the Epetra_BlockMap for this multi-vector.
const Epetra_Comm & Comm() const
Returns the address of the Epetra_Comm for this multi-vector.
Epetra_Distributor: The Epetra Gather/Scatter Setup Base Class.
void Add(const long long key, const value_type value)
value_type Get(const long long key)
Epetra_MapColoring: A class for coloring Epetra_Map and Epetra_BlockMap objects.
int NumColors() const
Returns number of colors on the calling processor.
virtual void Print(std::ostream &os) const
Print method.
int CopyAndPermute(const Epetra_SrcDistObject &Source, int NumSameIDs, int NumPermuteIDs, int *PermuteToLIDs, int *PermuteFromLIDs, const Epetra_OffsetIndex *Indexor, Epetra_CombineMode CombineMode=Zero)
Perform ID copies and permutations that are on processor.
int PackAndPrepare(const Epetra_SrcDistObject &Source, int NumExportIDs, int *ExportLIDs, int &LenExports, char *&Exports, int &SizeOfPacket, int *Sizes, bool &VarSizes, Epetra_Distributor &Distor)
Perform any packing or preparation required for call to DoTransfer().
Epetra_Map * TGenerateMap(int Color) const
int Allocate(int *ElementColors, int Increment)
int NumElementsWithColor(int Color) const
Returns number of map elements on calling processor having specified Color.
int * ListOfColors() const
Array of length NumColors() containing List of color values used in this coloring.
int DefaultColor() const
Returns default color.
Epetra_BlockMap * GenerateBlockMap(int Color) const
Generates an Epetra_BlockMap of the GIDs associated with the specified color.
int UnpackAndCombine(const Epetra_SrcDistObject &Source, int NumImportIDs, int *ImportLIDs, int LenImports, char *Imports, int &SizeOfPacket, Epetra_Distributor &Distor, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor)
Perform any unpacking and combining after call to DoTransfer().
Epetra_BlockMap * TGenerateBlockMap(int Color) const
int MaxNumColors() const
Returns maximum over all processors of the number of colors.
virtual ~Epetra_MapColoring()
Epetra_MapColoring destructor.
Epetra_HashTable< int > * ColorIDs_
Epetra_Map * GenerateMap(int Color) const
Generates an Epetra_Map of the GIDs associated with the specified color.
bool InItemList(int ColorValue) const
Epetra_MapColoring(const Epetra_BlockMap &Map, const int DefaultColor=0)
Epetra_MapColoring basic constructor.
int * ElementColors() const
Returns pointer to array of the colors associated with the LIDs on the calling processor.
int CheckSizes(const Epetra_SrcDistObject &A)
Allows the source and target (this) objects to be compared for compatibility, return nonzero if not.
int * ColorLIDList(int Color) const
Returns pointer to array of Map LIDs associated with the specified color.
Epetra_Map: A class for partitioning vectors and matrices.
Definition Epetra_Map.h:119
virtual int ReportError(const std::string Message, int ErrorCode) const
Error reporting method.
Epetra_OffsetIndex: This class builds index for efficient mapping of data from one Epetra_CrsGraph ba...
Epetra_SrcDistObject: A class for supporting flexible source distributed objects for import/export op...
Epetra_Util: The Epetra Util Wrapper Class.
Definition Epetra_Util.h:79
static void EPETRA_LIB_DLL_EXPORT Sort(bool SortAscending, int NumKeys, T *Keys, int NumDoubleCompanions, double **DoubleCompanions, int NumIntCompanions, int **IntCompanions, int NumLongLongCompanions, long long **LongLongCompanions)
Epetra_Util Sort Routine (Shell sort)