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++) ocolors[i] = maxcolors;
 69:   }

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

166: static PetscErrorCode GreedyColoringLocalDistanceTwo_Private(MatColoring mc, PetscReal *wts, PetscInt *lperm, ISColoringValue *colors)
167: {
168:   MC_Greedy       *gr = (MC_Greedy *)mc->data;
169:   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;
170:   Mat              m   = mc->mat, mt;
171:   Mat_MPIAIJ      *aij = (Mat_MPIAIJ *)m->data;
172:   Mat              md = NULL, mo = NULL;
173:   const PetscInt  *md_i, *mo_i, *md_j, *mo_j;
174:   const PetscInt  *rmd_i, *rmo_i, *rmd_j, *rmo_j;
175:   PetscBool        isMPIAIJ, isSEQAIJ;
176:   PetscInt         pcol, *dcolors, *ocolors;
177:   ISColoringValue *badidx;
178:   const PetscInt  *cidx;
179:   PetscReal       *owts, *colorweights;
180:   PetscInt        *oconf, *conf;
181:   PetscSF          sf;
182:   PetscLayout      layout;

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

479: static PetscErrorCode MatColoringApply_Greedy(MatColoring mc, ISColoring *iscoloring)
480: {
481:   PetscInt         finalcolor, finalcolor_global;
482:   ISColoringValue *colors;
483:   PetscInt         ncolstotal, ncols;
484:   PetscReal       *wts;
485:   PetscInt         i, *lperm;

487:   MatGetSize(mc->mat, NULL, &ncolstotal);
488:   MatGetLocalSize(mc->mat, NULL, &ncols);
489:   if (!mc->user_weights) {
490:     MatColoringCreateWeights(mc, &wts, &lperm);
491:   } else {
492:     wts   = mc->user_weights;
493:     lperm = mc->user_lperm;
494:   }
495:   PetscMalloc1(ncols, &colors);
496:   if (mc->dist == 1) {
497:     GreedyColoringLocalDistanceOne_Private(mc, wts, lperm, colors);
498:   } else if (mc->dist == 2) {
499:     GreedyColoringLocalDistanceTwo_Private(mc, wts, lperm, colors);
500:   } else SETERRQ(PetscObjectComm((PetscObject)mc), PETSC_ERR_ARG_OUTOFRANGE, "Only distance 1 and distance 2 supported by MatColoringGreedy");
501:   finalcolor = 0;
502:   for (i = 0; i < ncols; i++) {
503:     if (colors[i] > finalcolor) finalcolor = colors[i];
504:   }
505:   finalcolor_global = 0;
506:   MPIU_Allreduce(&finalcolor, &finalcolor_global, 1, MPIU_INT, MPI_MAX, PetscObjectComm((PetscObject)mc));
507:   PetscLogEventBegin(MATCOLORING_ISCreate, mc, 0, 0, 0);
508:   ISColoringCreate(PetscObjectComm((PetscObject)mc), finalcolor_global + 1, ncols, colors, PETSC_OWN_POINTER, iscoloring);
509:   PetscLogEventEnd(MATCOLORING_ISCreate, mc, 0, 0, 0);
510:   if (!mc->user_weights) {
511:     PetscFree(wts);
512:     PetscFree(lperm);
513:   }
514:   return 0;
515: }

517: static PetscErrorCode MatColoringSetFromOptions_Greedy(MatColoring mc, PetscOptionItems *PetscOptionsObject)
518: {
519:   MC_Greedy *gr = (MC_Greedy *)mc->data;

521:   PetscOptionsHeadBegin(PetscOptionsObject, "Greedy options");
522:   PetscOptionsBool("-mat_coloring_greedy_symmetric", "Flag for assuming a symmetric matrix", "", gr->symmetric, &gr->symmetric, NULL);
523:   PetscOptionsHeadEnd();
524:   return 0;
525: }

527: /*MC
528:   MATCOLORINGGREEDY - Greedy-with-conflict correction based matrix coloring for distance 1 and 2.

530:    Level: beginner

532:    Notes:

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

543:    Supports both distance one and distance two colorings.

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

549: .seealso: `MatColoringType`, `MatColoringCreate()`, `MatColoring`, `MatColoringSetType()`, `MatColoringType`
550: M*/
551: PETSC_EXTERN PetscErrorCode MatColoringCreate_Greedy(MatColoring mc)
552: {
553:   MC_Greedy *gr;

555:   PetscNew(&gr);
556:   mc->data                = gr;
557:   mc->ops->apply          = MatColoringApply_Greedy;
558:   mc->ops->view           = NULL;
559:   mc->ops->destroy        = MatColoringDestroy_Greedy;
560:   mc->ops->setfromoptions = MatColoringSetFromOptions_Greedy;

562:   gr->symmetric = PETSC_TRUE;
563:   return 0;
564: }