Actual source code: greedy.c

  1: #include <petsc/private/matimpl.h>
  2: #include <../src/mat/impls/aij/seq/aij.h>
  3: #include <../src/mat/impls/aij/mpi/mpiaij.h>
  4: #include <petscsf.h>

  6: typedef struct {
  7:   PetscBool symmetric;
  8: } MC_Greedy;

 10: static PetscErrorCode MatColoringDestroy_Greedy(MatColoring mc)
 11: {
 12:   PetscFree(mc->data);
 13:   return 0;
 14: }

 16: static PetscErrorCode GreedyColoringLocalDistanceOne_Private(MatColoring mc,PetscReal *wts,PetscInt *lperm,ISColoringValue *colors)
 17: {
 18:   PetscInt        i,j,k,s,e,n,no,nd,nd_global,n_global,idx,ncols,maxcolors,masksize,ccol,*mask;
 19:   Mat             m=mc->mat;
 20:   Mat_MPIAIJ      *aij = (Mat_MPIAIJ*)m->data;
 21:   Mat             md=NULL,mo=NULL;
 22:   const PetscInt  *md_i,*mo_i,*md_j,*mo_j;
 23:   PetscBool       isMPIAIJ,isSEQAIJ;
 24:   ISColoringValue pcol;
 25:   const PetscInt  *cidx;
 26:   PetscInt        *lcolors,*ocolors;
 27:   PetscReal       *owts=NULL;
 28:   PetscSF         sf;
 29:   PetscLayout     layout;

 31:   MatGetSize(m,&n_global,NULL);
 32:   MatGetOwnershipRange(m,&s,&e);
 33:   n=e-s;
 34:   masksize=20;
 35:   nd_global = 0;
 36:   /* get the matrix communication structures */
 37:   PetscObjectTypeCompare((PetscObject)m, MATMPIAIJ, &isMPIAIJ);
 38:   PetscObjectBaseTypeCompare((PetscObject)m, MATSEQAIJ, &isSEQAIJ);
 39:   if (isMPIAIJ) {
 40:     /* get the CSR data for on and off diagonal portions of m */
 41:     Mat_SeqAIJ *dseq;
 42:     Mat_SeqAIJ *oseq;
 43:     md=aij->A;
 44:     dseq = (Mat_SeqAIJ*)md->data;
 45:     mo=aij->B;
 46:     oseq = (Mat_SeqAIJ*)mo->data;
 47:     md_i = dseq->i;
 48:     md_j = dseq->j;
 49:     mo_i = oseq->i;
 50:     mo_j = oseq->j;
 51:   } else if (isSEQAIJ) {
 52:     /* get the CSR data for m */
 53:     Mat_SeqAIJ *dseq;
 54:     /* no off-processor nodes */
 55:     md=m;
 56:     dseq = (Mat_SeqAIJ*)md->data;
 57:     mo=NULL;
 58:     no=0;
 59:     md_i = dseq->i;
 60:     md_j = dseq->j;
 61:     mo_i = NULL;
 62:     mo_j = NULL;
 63:   } else SETERRQ(PetscObjectComm((PetscObject)mc),PETSC_ERR_ARG_WRONG,"Matrix must be AIJ for greedy coloring");
 64:   MatColoringGetMaxColors(mc,&maxcolors);
 65:   if (mo) {
 66:     VecGetSize(aij->lvec,&no);
 67:     PetscMalloc2(no,&ocolors,no,&owts);
 68:     for (i=0;i<no;i++) {
 69:       ocolors[i]=maxcolors;
 70:     }
 71:   }

 73:   PetscMalloc1(masksize,&mask);
 74:   PetscMalloc1(n,&lcolors);
 75:   for (i=0;i<n;i++) {
 76:     lcolors[i]=maxcolors;
 77:   }
 78:   for (i=0;i<masksize;i++) {
 79:     mask[i]=-1;
 80:   }
 81:   if (mo) {
 82:     /* transfer neighbor weights */
 83:     PetscSFCreate(PetscObjectComm((PetscObject)m),&sf);
 84:     MatGetLayouts(m,&layout,NULL);
 85:     PetscSFSetGraphLayout(sf,layout,no,NULL,PETSC_COPY_VALUES,aij->garray);
 86:     PetscSFBcastBegin(sf,MPIU_REAL,wts,owts,MPI_REPLACE);
 87:     PetscSFBcastEnd(sf,MPIU_REAL,wts,owts,MPI_REPLACE);
 88:   }
 89:   while (nd_global < n_global) {
 90:     nd=n;
 91:     /* assign lowest possible color to each local vertex */
 92:     PetscLogEventBegin(MATCOLORING_Local,mc,0,0,0);
 93:     for (i=0;i<n;i++) {
 94:       idx=lperm[i];
 95:       if (lcolors[idx] == maxcolors) {
 96:         ncols = md_i[idx+1]-md_i[idx];
 97:         cidx = &(md_j[md_i[idx]]);
 98:         for (j=0;j<ncols;j++) {
 99:           if (lcolors[cidx[j]] != maxcolors) {
100:             ccol=lcolors[cidx[j]];
101:             if (ccol>=masksize) {
102:               PetscInt *newmask;
103:               PetscMalloc1(masksize*2,&newmask);
104:               for (k=0;k<2*masksize;k++) {
105:                 newmask[k]=-1;
106:               }
107:               for (k=0;k<masksize;k++) {
108:                 newmask[k]=mask[k];
109:               }
110:               PetscFree(mask);
111:               mask=newmask;
112:               masksize*=2;
113:             }
114:             mask[ccol]=idx;
115:           }
116:         }
117:         if (mo) {
118:           ncols = mo_i[idx+1]-mo_i[idx];
119:           cidx = &(mo_j[mo_i[idx]]);
120:           for (j=0;j<ncols;j++) {
121:             if (ocolors[cidx[j]] != maxcolors) {
122:               ccol=ocolors[cidx[j]];
123:               if (ccol>=masksize) {
124:                 PetscInt *newmask;
125:                 PetscMalloc1(masksize*2,&newmask);
126:                 for (k=0;k<2*masksize;k++) {
127:                   newmask[k]=-1;
128:                 }
129:                 for (k=0;k<masksize;k++) {
130:                   newmask[k]=mask[k];
131:                 }
132:                 PetscFree(mask);
133:                 mask=newmask;
134:                 masksize*=2;
135:               }
136:               mask[ccol]=idx;
137:             }
138:           }
139:         }
140:         for (j=0;j<masksize;j++) {
141:           if (mask[j]!=idx) {
142:             break;
143:           }
144:         }
145:         pcol=j;
146:         if (pcol>maxcolors)pcol=maxcolors;
147:         lcolors[idx]=pcol;
148:       }
149:     }
150:     PetscLogEventEnd(MATCOLORING_Local,mc,0,0,0);
151:     if (mo) {
152:       /* transfer neighbor colors */
153:       PetscLogEventBegin(MATCOLORING_Comm,mc,0,0,0);
154:       PetscSFBcastBegin(sf,MPIU_INT,lcolors,ocolors,MPI_REPLACE);
155:       PetscSFBcastEnd(sf,MPIU_INT,lcolors,ocolors,MPI_REPLACE);
156:       /* check for conflicts -- this is merely checking if any adjacent off-processor rows have the same color and marking the ones that are lower weight locally for changing */
157:       for (i=0;i<n;i++) {
158:         ncols = mo_i[i+1]-mo_i[i];
159:         cidx = &(mo_j[mo_i[i]]);
160:         for (j=0;j<ncols;j++) {
161:           /* in the case of conflicts, the highest weight one stays and the others go */
162:           if ((ocolors[cidx[j]] == lcolors[i]) && (owts[cidx[j]] > wts[i]) && lcolors[i] < maxcolors) {
163:             lcolors[i]=maxcolors;
164:             nd--;
165:           }
166:         }
167:       }
168:       nd_global=0;
169:     }
170:     MPIU_Allreduce(&nd,&nd_global,1,MPIU_INT,MPI_SUM,PetscObjectComm((PetscObject)mc));
171:   }
172:   for (i=0;i<n;i++) {
173:     colors[i] = (ISColoringValue)lcolors[i];
174:   }
175:   PetscFree(mask);
176:   PetscFree(lcolors);
177:   if (mo) {
178:     PetscFree2(ocolors,owts);
179:     PetscSFDestroy(&sf);
180:   }
181:   return 0;
182: }

184: static PetscErrorCode GreedyColoringLocalDistanceTwo_Private(MatColoring mc,PetscReal *wts,PetscInt *lperm,ISColoringValue *colors)
185: {
186:   MC_Greedy       *gr = (MC_Greedy *) mc->data;
187:   PetscInt        i,j,k,l,s,e,n,nd,nd_global,n_global,idx,ncols,maxcolors,mcol,mcol_global,nd1cols,*mask,masksize,*d1cols,*bad,*badnext,nbad,badsize,ccol,no,cbad;
188:   Mat             m = mc->mat, mt;
189:   Mat_MPIAIJ      *aij = (Mat_MPIAIJ*)m->data;
190:   Mat             md=NULL,mo=NULL;
191:   const PetscInt  *md_i,*mo_i,*md_j,*mo_j;
192:   const PetscInt  *rmd_i,*rmo_i,*rmd_j,*rmo_j;
193:   PetscBool       isMPIAIJ,isSEQAIJ;
194:   PetscInt        pcol,*dcolors,*ocolors;
195:   ISColoringValue *badidx;
196:   const PetscInt  *cidx;
197:   PetscReal       *owts,*colorweights;
198:   PetscInt        *oconf,*conf;
199:   PetscSF         sf;
200:   PetscLayout     layout;

202:   MatGetSize(m,&n_global,NULL);
203:   MatGetOwnershipRange(m,&s,&e);
204:   n=e-s;
205:   nd_global = 0;
206:   /* get the matrix communication structures */
207:   PetscObjectBaseTypeCompare((PetscObject)m, MATMPIAIJ, &isMPIAIJ);
208:   PetscObjectBaseTypeCompare((PetscObject)m, MATSEQAIJ, &isSEQAIJ);
209:   if (isMPIAIJ) {
210:     Mat_SeqAIJ *dseq;
211:     Mat_SeqAIJ *oseq;
212:     md=aij->A;
213:     dseq = (Mat_SeqAIJ*)md->data;
214:     mo=aij->B;
215:     oseq = (Mat_SeqAIJ*)mo->data;
216:     md_i = dseq->i;
217:     md_j = dseq->j;
218:     mo_i = oseq->i;
219:     mo_j = oseq->j;
220:     rmd_i = dseq->i;
221:     rmd_j = dseq->j;
222:     rmo_i = oseq->i;
223:     rmo_j = oseq->j;
224:   } else if (isSEQAIJ) {
225:     Mat_SeqAIJ *dseq;
226:     /* no off-processor nodes */
227:     md=m;
228:     dseq = (Mat_SeqAIJ*)md->data;
229:     md_i = dseq->i;
230:     md_j = dseq->j;
231:     mo_i = NULL;
232:     mo_j = NULL;
233:     rmd_i = dseq->i;
234:     rmd_j = dseq->j;
235:     rmo_i = NULL;
236:     rmo_j = NULL;
237:   } else SETERRQ(PetscObjectComm((PetscObject)mc),PETSC_ERR_ARG_WRONG,"Matrix must be AIJ for greedy coloring");
238:   if (!gr->symmetric) {
239:     MatTranspose(m, MAT_INITIAL_MATRIX, &mt);
240:     if (isSEQAIJ) {
241:       Mat_SeqAIJ *dseq = (Mat_SeqAIJ*) mt->data;
242:       rmd_i = dseq->i;
243:       rmd_j = dseq->j;
244:       rmo_i = NULL;
245:       rmo_j = NULL;
246:     } else SETERRQ(PetscObjectComm((PetscObject) mc), PETSC_ERR_SUP, "Nonsymmetric greedy coloring only works in serial");
247:   }
248:   /* create the vectors and communication structures if necessary */
249:   no=0;
250:   if (mo) {
251:     VecGetLocalSize(aij->lvec,&no);
252:     PetscSFCreate(PetscObjectComm((PetscObject)m),&sf);
253:     MatGetLayouts(m,&layout,NULL);
254:     PetscSFSetGraphLayout(sf,layout,no,NULL,PETSC_COPY_VALUES,aij->garray);
255:   }
256:   MatColoringGetMaxColors(mc,&maxcolors);
257:   masksize=n;
258:   nbad=0;
259:   badsize=n;
260:   PetscMalloc1(masksize,&mask);
261:   PetscMalloc4(n,&d1cols,n,&dcolors,n,&conf,n,&bad);
262:   PetscMalloc2(badsize,&badidx,badsize,&badnext);
263:   for (i=0;i<masksize;i++) {
264:     mask[i]=-1;
265:   }
266:   for (i=0;i<n;i++) {
267:     dcolors[i]=maxcolors;
268:     bad[i]=-1;
269:   }
270:   for (i=0;i<badsize;i++) {
271:     badnext[i]=-1;
272:   }
273:   if (mo) {
274:     PetscMalloc3(no,&owts,no,&oconf,no,&ocolors);
275:     PetscSFBcastBegin(sf,MPIU_REAL,wts,owts,MPI_REPLACE);
276:     PetscSFBcastEnd(sf,MPIU_REAL,wts,owts,MPI_REPLACE);
277:     for (i=0;i<no;i++) {
278:       ocolors[i]=maxcolors;
279:     }
280:   } else {                      /* Appease overzealous -Wmaybe-initialized */
281:     owts = NULL;
282:     oconf = NULL;
283:     ocolors = NULL;
284:   }
285:   mcol=0;
286:   while (nd_global < n_global) {
287:     nd=n;
288:     /* assign lowest possible color to each local vertex */
289:     mcol_global=0;
290:     PetscLogEventBegin(MATCOLORING_Local,mc,0,0,0);
291:     for (i=0;i<n;i++) {
292:       idx=lperm[i];
293:       if (dcolors[idx] == maxcolors) {
294:         /* entries in bad */
295:         cbad=bad[idx];
296:         while (cbad>=0) {
297:           ccol=badidx[cbad];
298:           if (ccol>=masksize) {
299:             PetscInt *newmask;
300:             PetscMalloc1(masksize*2,&newmask);
301:             for (k=0;k<2*masksize;k++) {
302:               newmask[k]=-1;
303:             }
304:             for (k=0;k<masksize;k++) {
305:               newmask[k]=mask[k];
306:             }
307:             PetscFree(mask);
308:             mask=newmask;
309:             masksize*=2;
310:           }
311:           mask[ccol]=idx;
312:           cbad=badnext[cbad];
313:         }
314:         /* diagonal distance-one rows */
315:         nd1cols=0;
316:         ncols = rmd_i[idx+1]-rmd_i[idx];
317:         cidx = &(rmd_j[rmd_i[idx]]);
318:         for (j=0;j<ncols;j++) {
319:           d1cols[nd1cols] = cidx[j];
320:           nd1cols++;
321:           ccol=dcolors[cidx[j]];
322:           if (ccol != maxcolors) {
323:             if (ccol>=masksize) {
324:               PetscInt *newmask;
325:               PetscMalloc1(masksize*2,&newmask);
326:               for (k=0;k<2*masksize;k++) {
327:                 newmask[k]=-1;
328:               }
329:               for (k=0;k<masksize;k++) {
330:                 newmask[k]=mask[k];
331:               }
332:               PetscFree(mask);
333:               mask=newmask;
334:               masksize*=2;
335:             }
336:             mask[ccol]=idx;
337:           }
338:         }
339:         /* off-diagonal distance-one rows */
340:         if (mo) {
341:           ncols = rmo_i[idx+1]-rmo_i[idx];
342:           cidx = &(rmo_j[rmo_i[idx]]);
343:           for (j=0;j<ncols;j++) {
344:             ccol=ocolors[cidx[j]];
345:             if (ccol != maxcolors) {
346:               if (ccol>=masksize) {
347:                 PetscInt *newmask;
348:                 PetscMalloc1(masksize*2,&newmask);
349:                 for (k=0;k<2*masksize;k++) {
350:                   newmask[k]=-1;
351:                 }
352:                 for (k=0;k<masksize;k++) {
353:                   newmask[k]=mask[k];
354:                 }
355:                 PetscFree(mask);
356:                 mask=newmask;
357:                 masksize*=2;
358:               }
359:               mask[ccol]=idx;
360:             }
361:           }
362:         }
363:         /* diagonal distance-two rows */
364:         for (j=0;j<nd1cols;j++) {
365:           ncols = md_i[d1cols[j]+1]-md_i[d1cols[j]];
366:           cidx = &(md_j[md_i[d1cols[j]]]);
367:           for (l=0;l<ncols;l++) {
368:             ccol=dcolors[cidx[l]];
369:             if (ccol != maxcolors) {
370:               if (ccol>=masksize) {
371:                 PetscInt *newmask;
372:                 PetscMalloc1(masksize*2,&newmask);
373:                 for (k=0;k<2*masksize;k++) {
374:                   newmask[k]=-1;
375:                 }
376:                 for (k=0;k<masksize;k++) {
377:                   newmask[k]=mask[k];
378:                 }
379:                 PetscFree(mask);
380:                 mask=newmask;
381:                 masksize*=2;
382:               }
383:               mask[ccol]=idx;
384:             }
385:           }
386:         }
387:         /* off-diagonal distance-two rows */
388:         if (mo) {
389:           for (j=0;j<nd1cols;j++) {
390:             ncols = mo_i[d1cols[j]+1]-mo_i[d1cols[j]];
391:             cidx = &(mo_j[mo_i[d1cols[j]]]);
392:             for (l=0;l<ncols;l++) {
393:               ccol=ocolors[cidx[l]];
394:               if (ccol != maxcolors) {
395:                 if (ccol>=masksize) {
396:                   PetscInt *newmask;
397:                   PetscMalloc1(masksize*2,&newmask);
398:                   for (k=0;k<2*masksize;k++) {
399:                     newmask[k]=-1;
400:                   }
401:                   for (k=0;k<masksize;k++) {
402:                     newmask[k]=mask[k];
403:                   }
404:                   PetscFree(mask);
405:                   mask=newmask;
406:                   masksize*=2;
407:                 }
408:                 mask[ccol]=idx;
409:               }
410:             }
411:           }
412:         }
413:         /* assign this one the lowest color possible by seeing if there's a gap in the sequence of sorted neighbor colors */
414:         for (j=0;j<masksize;j++) {
415:           if (mask[j]!=idx) {
416:             break;
417:           }
418:         }
419:         pcol=j;
420:         if (pcol>maxcolors) pcol=maxcolors;
421:         dcolors[idx]=pcol;
422:         if (pcol>mcol) mcol=pcol;
423:       }
424:     }
425:     PetscLogEventEnd(MATCOLORING_Local,mc,0,0,0);
426:     if (mo) {
427:       /* transfer neighbor colors */
428:       PetscSFBcastBegin(sf,MPIU_INT,dcolors,ocolors,MPI_REPLACE);
429:       PetscSFBcastEnd(sf,MPIU_INT,dcolors,ocolors,MPI_REPLACE);
430:       /* find the maximum color assigned locally and allocate a mask */
431:       MPIU_Allreduce(&mcol,&mcol_global,1,MPIU_INT,MPI_MAX,PetscObjectComm((PetscObject)mc));
432:       PetscMalloc1(mcol_global+1,&colorweights);
433:       /* check for conflicts */
434:       for (i=0;i<n;i++) {
435:         conf[i]=PETSC_FALSE;
436:       }
437:       for (i=0;i<no;i++) {
438:         oconf[i]=PETSC_FALSE;
439:       }
440:       for (i=0;i<n;i++) {
441:         ncols = mo_i[i+1]-mo_i[i];
442:         cidx = &(mo_j[mo_i[i]]);
443:         if (ncols > 0) {
444:           /* fill in the mask */
445:           for (j=0;j<mcol_global+1;j++) {
446:             colorweights[j]=0;
447:           }
448:           colorweights[dcolors[i]]=wts[i];
449:           /* fill in the off-diagonal part of the mask */
450:           for (j=0;j<ncols;j++) {
451:             ccol=ocolors[cidx[j]];
452:             if (ccol < maxcolors) {
453:               if (colorweights[ccol] < owts[cidx[j]]) {
454:                 colorweights[ccol] = owts[cidx[j]];
455:               }
456:             }
457:           }
458:           /* fill in the on-diagonal part of the mask */
459:           ncols = md_i[i+1]-md_i[i];
460:           cidx = &(md_j[md_i[i]]);
461:           for (j=0;j<ncols;j++) {
462:             ccol=dcolors[cidx[j]];
463:             if (ccol < maxcolors) {
464:               if (colorweights[ccol] < wts[cidx[j]]) {
465:                 colorweights[ccol] = wts[cidx[j]];
466:               }
467:             }
468:           }
469:           /* go back through and set up on and off-diagonal conflict vectors */
470:           ncols = md_i[i+1]-md_i[i];
471:           cidx = &(md_j[md_i[i]]);
472:           for (j=0;j<ncols;j++) {
473:             ccol=dcolors[cidx[j]];
474:             if (ccol < maxcolors) {
475:               if (colorweights[ccol] > wts[cidx[j]]) {
476:                 conf[cidx[j]]=PETSC_TRUE;
477:               }
478:             }
479:           }
480:           ncols = mo_i[i+1]-mo_i[i];
481:           cidx = &(mo_j[mo_i[i]]);
482:           for (j=0;j<ncols;j++) {
483:             ccol=ocolors[cidx[j]];
484:             if (ccol < maxcolors) {
485:               if (colorweights[ccol] > owts[cidx[j]]) {
486:                 oconf[cidx[j]]=PETSC_TRUE;
487:               }
488:             }
489:           }
490:         }
491:       }
492:       nd_global=0;
493:       PetscFree(colorweights);
494:       PetscLogEventBegin(MATCOLORING_Comm,mc,0,0,0);
495:       PetscSFReduceBegin(sf,MPIU_INT,oconf,conf,MPIU_SUM);
496:       PetscSFReduceEnd(sf,MPIU_INT,oconf,conf,MPIU_SUM);
497:       PetscLogEventEnd(MATCOLORING_Comm,mc,0,0,0);
498:       /* go through and unset local colors that have conflicts */
499:       for (i=0;i<n;i++) {
500:         if (conf[i]>0) {
501:           /* push this color onto the bad stack */
502:           badidx[nbad]=dcolors[i];
503:           badnext[nbad]=bad[i];
504:           bad[i]=nbad;
505:           nbad++;
506:           if (nbad>=badsize) {
507:             PetscInt *newbadnext;
508:             ISColoringValue *newbadidx;
509:             PetscMalloc2(badsize*2,&newbadidx,badsize*2,&newbadnext);
510:             for (k=0;k<2*badsize;k++) {
511:               newbadnext[k]=-1;
512:             }
513:             for (k=0;k<badsize;k++) {
514:               newbadidx[k]=badidx[k];
515:               newbadnext[k]=badnext[k];
516:             }
517:             PetscFree2(badidx,badnext);
518:             badidx=newbadidx;
519:             badnext=newbadnext;
520:             badsize*=2;
521:           }
522:           dcolors[i] = maxcolors;
523:           nd--;
524:         }
525:       }
526:     }
527:     MPIU_Allreduce(&nd,&nd_global,1,MPIU_INT,MPI_SUM,PetscObjectComm((PetscObject)mc));
528:   }
529:   if (mo) {
530:     PetscSFDestroy(&sf);
531:     PetscFree3(owts,oconf,ocolors);
532:   }
533:   for (i=0;i<n;i++) {
534:     colors[i]=dcolors[i];
535:   }
536:   PetscFree(mask);
537:   PetscFree4(d1cols,dcolors,conf,bad);
538:   PetscFree2(badidx,badnext);
539:   if (!gr->symmetric) MatDestroy(&mt);
540:   return 0;
541: }

543: static PetscErrorCode MatColoringApply_Greedy(MatColoring mc,ISColoring *iscoloring)
544: {
545:   PetscInt        finalcolor,finalcolor_global;
546:   ISColoringValue *colors;
547:   PetscInt        ncolstotal,ncols;
548:   PetscReal       *wts;
549:   PetscInt        i,*lperm;

551:   MatGetSize(mc->mat,NULL,&ncolstotal);
552:   MatGetLocalSize(mc->mat,NULL,&ncols);
553:   if (!mc->user_weights) {
554:     MatColoringCreateWeights(mc,&wts,&lperm);
555:   } else {
556:     wts = mc->user_weights;
557:     lperm = mc->user_lperm;
558:   }
559:   PetscMalloc1(ncols,&colors);
560:   if (mc->dist == 1) {
561:     GreedyColoringLocalDistanceOne_Private(mc,wts,lperm,colors);
562:   } else if (mc->dist == 2) {
563:     GreedyColoringLocalDistanceTwo_Private(mc,wts,lperm,colors);
564:   } else SETERRQ(PetscObjectComm((PetscObject)mc),PETSC_ERR_ARG_OUTOFRANGE,"Only distance 1 and distance 2 supported by MatColoringGreedy");
565:   finalcolor=0;
566:   for (i=0;i<ncols;i++) {
567:     if (colors[i] > finalcolor) finalcolor=colors[i];
568:   }
569:   finalcolor_global=0;
570:   MPIU_Allreduce(&finalcolor,&finalcolor_global,1,MPIU_INT,MPI_MAX,PetscObjectComm((PetscObject)mc));
571:   PetscLogEventBegin(MATCOLORING_ISCreate,mc,0,0,0);
572:   ISColoringCreate(PetscObjectComm((PetscObject)mc),finalcolor_global+1,ncols,colors,PETSC_OWN_POINTER,iscoloring);
573:   PetscLogEventEnd(MATCOLORING_ISCreate,mc,0,0,0);
574:   if (!mc->user_weights) {
575:     PetscFree(wts);
576:     PetscFree(lperm);
577:   }
578:   return 0;
579: }

581: static PetscErrorCode MatColoringSetFromOptions_Greedy(PetscOptionItems *PetscOptionsObject, MatColoring mc)
582: {
583:   MC_Greedy     *gr = (MC_Greedy *) mc->data;

585:   PetscOptionsHead(PetscOptionsObject, "Greedy options");
586:   PetscOptionsBool("-mat_coloring_greedy_symmetric", "Flag for assuming a symmetric matrix", "", gr->symmetric, &gr->symmetric, NULL);
587:   PetscOptionsTail();
588:   return 0;
589: }

591: /*MC
592:   MATCOLORINGGREEDY - Greedy-with-conflict correction based Matrix Coloring for distance 1 and 2.

594:    Level: beginner

596:    Notes:

598:    These algorithms proceed in two phases -- local coloring and conflict resolution.  The local coloring
599:    tentatively colors all vertices at the distance required given what's known of the global coloring.  Then,
600:    the updated colors are transferred to different processors at distance one.  In the distance one case, each
601:    vertex with nonlocal neighbors is then checked to see if it conforms, with the vertex being
602:    marked for recoloring if its lower weight than its same colored neighbor.  In the distance two case,
603:    each boundary vertex's immediate star is checked for validity of the coloring.  Lower-weight conflict
604:    vertices are marked, and then the conflicts are gathered back on owning processors.  In both cases
605:    this is done until each column has received a valid color.

607:    Supports both distance one and distance two colorings.

609:    References:
610: .  * - Bozdag et al. "A Parallel Distance 2 Graph Coloring Algorithm for Distributed Memory Computers"
611:    HPCC'05 Proceedings of the First international conference on High Performance Computing and Communications

613: .seealso: MatColoringCreate(), MatColoring, MatColoringSetType(), MatColoringType
614: M*/
615: PETSC_EXTERN PetscErrorCode MatColoringCreate_Greedy(MatColoring mc)
616: {
617:   MC_Greedy      *gr;

619:   PetscNewLog(mc,&gr);
620:   mc->data                = gr;
621:   mc->ops->apply          = MatColoringApply_Greedy;
622:   mc->ops->view           = NULL;
623:   mc->ops->destroy        = MatColoringDestroy_Greedy;
624:   mc->ops->setfromoptions = MatColoringSetFromOptions_Greedy;

626:   gr->symmetric = PETSC_TRUE;
627:   return 0;
628: }