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: }