Actual source code: matcoloring.c

  1: #include <petsc/private/matimpl.h>

  3: PetscFunctionList MatColoringList              = NULL;
  4: PetscBool         MatColoringRegisterAllCalled = PETSC_FALSE;
  5: const char *const MatColoringWeightTypes[]     = {"RANDOM", "LEXICAL", "LF", "SL", "MatColoringWeightType", "MAT_COLORING_WEIGHT_", NULL};

  7: /*@C
  8:   MatColoringRegister - Adds a new sparse matrix coloring to the  matrix package.

 10:   Not Collective

 12:   Input Parameters:
 13: + sname    - name of Coloring (for example `MATCOLORINGSL`)
 14: - function - function pointer that creates the coloring

 16:   Level: developer

 18:   Example Usage:
 19: .vb
 20:    MatColoringRegister("my_color", MyColor);
 21: .ve

 23:   Then, your partitioner can be chosen with the procedural interface via
 24: $     MatColoringSetType(part, "my_color")
 25:   or at runtime via the option
 26: $     -mat_coloring_type my_color

 28: .seealso: `MatColoringType`, `MatColoringRegisterDestroy()`, `MatColoringRegisterAll()`
 29: @*/
 30: PetscErrorCode MatColoringRegister(const char sname[], PetscErrorCode (*function)(MatColoring))
 31: {
 32:   PetscFunctionBegin;
 33:   PetscCall(MatInitializePackage());
 34:   PetscCall(PetscFunctionListAdd(&MatColoringList, sname, function));
 35:   PetscFunctionReturn(PETSC_SUCCESS);
 36: }

 38: /*@
 39:   MatColoringCreate - Creates a matrix coloring context.

 41:   Collective

 43:   Input Parameter:
 44: . m - MPI communicator

 46:   Output Parameter:
 47: . mcptr - the new `MatColoring` context

 49:   Options Database Keys:
 50: + -mat_coloring_type      - the type of coloring algorithm used. See `MatColoringType`.
 51: . -mat_coloring_maxcolors - the maximum number of relevant colors, all nodes not in a color are in maxcolors+1
 52: . -mat_coloring_distance  - compute a distance 1,2,... coloring.
 53: . -mat_coloring_view      - print information about the coloring and the produced index sets
 54: . -mat_coloring_test      - debugging option that prints all coloring incompatibilities
 55: - -mat_is_coloring_test   - debugging option that throws an error if MatColoringApply() generates an incorrect iscoloring

 57:   Level: beginner

 59:   Notes:
 60:   A distance one coloring is useful, for example, multi-color SOR.

 62:   A distance two coloring is for the finite difference computation of Jacobians (see `MatFDColoringCreate()`).

 64:   Coloring of matrices can be computed directly from the sparse matrix nonzero structure via the `MatColoring` object or from the mesh from which the
 65:   matrix comes from with `DMCreateColoring()`. In general using the mesh produces a more optimal coloring (fewer colors).

 67:   Some coloring types only support distance two colorings

 69: .seealso: `MatColoringSetFromOptions()`, `MatColoring`, `MatColoringApply()`, `MatFDColoringCreate()`, `DMCreateColoring()`, `MatColoringType`
 70: @*/
 71: PetscErrorCode MatColoringCreate(Mat m, MatColoring *mcptr)
 72: {
 73:   MatColoring mc;

 75:   PetscFunctionBegin;
 77:   PetscAssertPointer(mcptr, 2);
 78:   *mcptr = NULL;

 80:   PetscCall(MatInitializePackage());
 81:   PetscCall(PetscHeaderCreate(mc, MAT_COLORING_CLASSID, "MatColoring", "Matrix coloring", "MatColoring", PetscObjectComm((PetscObject)m), MatColoringDestroy, MatColoringView));
 82:   PetscCall(PetscObjectReference((PetscObject)m));
 83:   mc->mat          = m;
 84:   mc->dist         = 2; /* default to Jacobian computation case */
 85:   mc->maxcolors    = IS_COLORING_MAX;
 86:   *mcptr           = mc;
 87:   mc->valid        = PETSC_FALSE;
 88:   mc->weight_type  = MAT_COLORING_WEIGHT_RANDOM;
 89:   mc->user_weights = NULL;
 90:   mc->user_lperm   = NULL;
 91:   PetscFunctionReturn(PETSC_SUCCESS);
 92: }

 94: /*@
 95:   MatColoringDestroy - Destroys the matrix coloring context

 97:   Collective

 99:   Input Parameter:
100: . mc - the `MatColoring` context

102:   Level: beginner

104: .seealso: `MatColoring`, `MatColoringCreate()`, `MatColoringApply()`
105: @*/
106: PetscErrorCode MatColoringDestroy(MatColoring *mc)
107: {
108:   PetscFunctionBegin;
109:   if (--((PetscObject)*mc)->refct > 0) {
110:     *mc = NULL;
111:     PetscFunctionReturn(PETSC_SUCCESS);
112:   }
113:   PetscCall(MatDestroy(&(*mc)->mat));
114:   PetscTryTypeMethod(*mc, destroy);
115:   PetscCall(PetscFree((*mc)->user_weights));
116:   PetscCall(PetscFree((*mc)->user_lperm));
117:   PetscCall(PetscHeaderDestroy(mc));
118:   PetscFunctionReturn(PETSC_SUCCESS);
119: }

121: /*@C
122:   MatColoringSetType - Sets the type of coloring algorithm used

124:   Collective

126:   Input Parameters:
127: + mc   - the `MatColoring` context
128: - type - the type of coloring

130:   Options Database Key:
131: . -mat_coloring_type type - the name of the type

133:   Level: beginner

135:   Note:
136:   Possible types include the sequential types `MATCOLORINGLF`,
137:   `MATCOLORINGSL`, and `MATCOLORINGID` from the MINPACK package as well
138:   as a parallel `MATCOLORINGGREEDY` algorithm.

140: .seealso: `MatColoring`, `MatColoringSetFromOptions()`, `MatColoringType`, `MatColoringCreate()`, `MatColoringApply()`
141: @*/
142: PetscErrorCode MatColoringSetType(MatColoring mc, MatColoringType type)
143: {
144:   PetscBool match;
145:   PetscErrorCode (*r)(MatColoring);

147:   PetscFunctionBegin;
149:   PetscAssertPointer(type, 2);
150:   PetscCall(PetscObjectTypeCompare((PetscObject)mc, type, &match));
151:   if (match) PetscFunctionReturn(PETSC_SUCCESS);
152:   PetscCall(PetscFunctionListFind(MatColoringList, type, &r));
153:   PetscCheck(r, PetscObjectComm((PetscObject)mc), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unable to find requested MatColoring type %s", type);

155:   PetscTryTypeMethod(mc, destroy);
156:   mc->ops->destroy        = NULL;
157:   mc->ops->apply          = NULL;
158:   mc->ops->view           = NULL;
159:   mc->ops->setfromoptions = NULL;
160:   mc->ops->destroy        = NULL;

162:   PetscCall(PetscObjectChangeTypeName((PetscObject)mc, type));
163:   PetscCall((*r)(mc));
164:   PetscFunctionReturn(PETSC_SUCCESS);
165: }

167: /*@
168:   MatColoringSetFromOptions - Sets `MatColoring` options from options database

170:   Collective

172:   Input Parameter:
173: . mc - `MatColoring` context

175:   Options Database Keys:
176: + -mat_coloring_type      - the type of coloring algorithm used. See `MatColoringType`.
177: . -mat_coloring_maxcolors - the maximum number of relevant colors, all nodes not in a color are in maxcolors+1
178: . -mat_coloring_distance  - compute a distance 1,2,... coloring.
179: . -mat_coloring_view      - print information about the coloring and the produced index sets
180: . -snes_fd_color          - instruct SNES to using coloring and then `MatFDColoring` to compute the Jacobians
181: - -snes_fd_color_use_mat  - instruct `SNES` to color the matrix directly instead of the `DM` from which the matrix comes (the default)

183:   Level: beginner

185: .seealso: `MatColoring`, `MatColoringApply()`, `MatColoringSetDistance()`, `MatColoringSetType()`, `SNESComputeJacobianDefaultColor()`, `MatColoringType`
186: @*/
187: PetscErrorCode MatColoringSetFromOptions(MatColoring mc)
188: {
189:   PetscBool       flg;
190:   MatColoringType deft = MATCOLORINGSL;
191:   char            type[256];
192:   PetscInt        dist, maxcolors;

194:   PetscFunctionBegin;
196:   PetscCall(MatColoringGetDistance(mc, &dist));
197:   if (dist == 2) deft = MATCOLORINGSL;
198:   else deft = MATCOLORINGGREEDY;
199:   PetscCall(MatColoringGetMaxColors(mc, &maxcolors));
200:   PetscCall(MatColoringRegisterAll());
201:   PetscObjectOptionsBegin((PetscObject)mc);
202:   if (((PetscObject)mc)->type_name) deft = ((PetscObject)mc)->type_name;
203:   PetscCall(PetscOptionsFList("-mat_coloring_type", "The coloring method used", "MatColoringSetType", MatColoringList, deft, type, 256, &flg));
204:   if (flg) {
205:     PetscCall(MatColoringSetType(mc, type));
206:   } else if (!((PetscObject)mc)->type_name) {
207:     PetscCall(MatColoringSetType(mc, deft));
208:   }
209:   PetscCall(PetscOptionsInt("-mat_coloring_distance", "Distance of the coloring", "MatColoringSetDistance", dist, &dist, &flg));
210:   if (flg) PetscCall(MatColoringSetDistance(mc, dist));
211:   PetscCall(PetscOptionsInt("-mat_coloring_maxcolors", "Maximum colors returned at the end. 1 returns an independent set", "MatColoringSetMaxColors", maxcolors, &maxcolors, &flg));
212:   if (flg) PetscCall(MatColoringSetMaxColors(mc, maxcolors));
213:   PetscTryTypeMethod(mc, setfromoptions, PetscOptionsObject);
214:   PetscCall(PetscOptionsBool("-mat_coloring_test", "Check that a valid coloring has been produced", "", mc->valid, &mc->valid, NULL));
215:   PetscCall(PetscOptionsBool("-mat_is_coloring_test", "Check that a valid iscoloring has been produced", "", mc->valid_iscoloring, &mc->valid_iscoloring, NULL));
216:   PetscCall(PetscOptionsEnum("-mat_coloring_weight_type", "Sets the type of vertex weighting used", "MatColoringSetWeightType", MatColoringWeightTypes, (PetscEnum)mc->weight_type, (PetscEnum *)&mc->weight_type, NULL));
217:   PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)mc, PetscOptionsObject));
218:   PetscOptionsEnd();
219:   PetscFunctionReturn(PETSC_SUCCESS);
220: }

222: /*@
223:   MatColoringSetDistance - Sets the distance of the coloring

225:   Logically Collective

227:   Input Parameters:
228: + mc   - the `MatColoring` context
229: - dist - the distance the coloring should compute

231:   Options Database Key:
232: . -mat_coloring_type - the type of coloring algorithm used. See `MatColoringType`.

234:   Level: beginner

236:   Note:
237:   The distance of the coloring denotes the minimum number
238:   of edges in the graph induced by the matrix any two vertices
239:   of the same color are.  Distance-1 colorings are the classical
240:   coloring, where no two vertices of the same color are adjacent.
241:   distance-2 colorings are useful for the computation of Jacobians.

243: .seealso: `MatColoring`, `MatColoringSetFromOptions()`, `MatColoringGetDistance()`, `MatColoringApply()`
244: @*/
245: PetscErrorCode MatColoringSetDistance(MatColoring mc, PetscInt dist)
246: {
247:   PetscFunctionBegin;
249:   mc->dist = dist;
250:   PetscFunctionReturn(PETSC_SUCCESS);
251: }

253: /*@
254:   MatColoringGetDistance - Gets the distance of the coloring

256:   Logically Collective

258:   Input Parameter:
259: . mc - the `MatColoring` context

261:   Output Parameter:
262: . dist - the current distance being used for the coloring.

264:   Level: beginner

266:   Note:
267:   The distance of the coloring denotes the minimum number
268:   of edges in the graph induced by the matrix any two vertices
269:   of the same color are.  Distance-1 colorings are the classical
270:   coloring, where no two vertices of the same color are adjacent.
271:   distance-2 colorings are useful for the computation of Jacobians.

273: .seealso: `MatColoring`, `MatColoringSetDistance()`, `MatColoringApply()`
274: @*/
275: PetscErrorCode MatColoringGetDistance(MatColoring mc, PetscInt *dist)
276: {
277:   PetscFunctionBegin;
279:   if (dist) *dist = mc->dist;
280:   PetscFunctionReturn(PETSC_SUCCESS);
281: }

283: /*@
284:   MatColoringSetMaxColors - Sets the maximum number of colors to produce

286:   Logically Collective

288:   Input Parameters:
289: + mc        - the `MatColoring` context
290: - maxcolors - the maximum number of colors to produce

292:   Level: beginner

294:   Notes:
295:   Vertices not in an available color are set to have color maxcolors+1, which is not
296:   a valid color as they may be adjacent.

298:   This works only for  `MATCOLORINGGREEDY` and `MATCOLORINGJP`

300:   This may be used to compute a certain number of
301:   independent sets from the graph.  For instance, while using
302:   `MATCOLORINGGREEDY` and maxcolors = 1, one gets out an MIS.

304: .seealso: `MatColoring`, `MatColoringGetMaxColors()`, `MatColoringApply()`
305: @*/
306: PetscErrorCode MatColoringSetMaxColors(MatColoring mc, PetscInt maxcolors)
307: {
308:   PetscFunctionBegin;
310:   mc->maxcolors = maxcolors;
311:   PetscFunctionReturn(PETSC_SUCCESS);
312: }

314: /*@
315:   MatColoringGetMaxColors - Gets the maximum number of colors

317:   Logically Collective

319:   Input Parameter:
320: . mc - the `MatColoring` context

322:   Output Parameter:
323: . maxcolors - the current maximum number of colors to produce

325:   Level: beginner

327: .seealso: `MatColoring`, `MatColoringSetMaxColors()`, `MatColoringApply()`
328: @*/
329: PetscErrorCode MatColoringGetMaxColors(MatColoring mc, PetscInt *maxcolors)
330: {
331:   PetscFunctionBegin;
333:   if (maxcolors) *maxcolors = mc->maxcolors;
334:   PetscFunctionReturn(PETSC_SUCCESS);
335: }

337: /*@
338:   MatColoringApply - Apply the coloring to the matrix, producing index
339:   sets corresponding to a number of independent sets in the induced
340:   graph.

342:   Collective

344:   Input Parameter:
345: . mc - the `MatColoring` context

347:   Output Parameter:
348: . coloring - the `ISColoring` instance containing the coloring

350:   Level: beginner

352: .seealso: `ISColoring`, `MatColoring`, `MatColoringCreate()`
353: @*/
354: PetscErrorCode MatColoringApply(MatColoring mc, ISColoring *coloring)
355: {
356:   PetscBool         flg;
357:   PetscViewerFormat format;
358:   PetscViewer       viewer;
359:   PetscInt          nc, ncolors;

361:   PetscFunctionBegin;
363:   PetscAssertPointer(coloring, 2);
364:   PetscCall(PetscLogEventBegin(MATCOLORING_Apply, mc, 0, 0, 0));
365:   PetscUseTypeMethod(mc, apply, coloring);
366:   PetscCall(PetscLogEventEnd(MATCOLORING_Apply, mc, 0, 0, 0));

368:   /* valid */
369:   if (mc->valid) PetscCall(MatColoringTest(mc, *coloring));
370:   if (mc->valid_iscoloring) PetscCall(MatISColoringTest(mc->mat, *coloring));

372:   /* view */
373:   PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)mc), ((PetscObject)mc)->options, ((PetscObject)mc)->prefix, "-mat_coloring_view", &viewer, &format, &flg));
374:   if (flg && !PetscPreLoadingOn) {
375:     PetscCall(PetscViewerPushFormat(viewer, format));
376:     PetscCall(MatColoringView(mc, viewer));
377:     PetscCall(MatGetSize(mc->mat, NULL, &nc));
378:     PetscCall(ISColoringGetIS(*coloring, PETSC_USE_POINTER, &ncolors, NULL));
379:     PetscCall(PetscViewerASCIIPrintf(viewer, "  Number of colors %" PetscInt_FMT "\n", ncolors));
380:     PetscCall(PetscViewerASCIIPrintf(viewer, "  Number of total columns %" PetscInt_FMT "\n", nc));
381:     if (nc <= 1000) PetscCall(ISColoringView(*coloring, viewer));
382:     PetscCall(PetscViewerPopFormat(viewer));
383:     PetscCall(PetscOptionsRestoreViewer(&viewer));
384:   }
385:   PetscFunctionReturn(PETSC_SUCCESS);
386: }

388: /*@
389:   MatColoringView - Output details about the `MatColoring`.

391:   Collective

393:   Input Parameters:
394: + mc     - the `MatColoring` context
395: - viewer - the Viewer context

397:   Level: beginner

399: .seealso: `PetscViewer`, `MatColoring`, `MatColoringApply()`
400: @*/
401: PetscErrorCode MatColoringView(MatColoring mc, PetscViewer viewer)
402: {
403:   PetscBool iascii;

405:   PetscFunctionBegin;
407:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)mc), &viewer));
409:   PetscCheckSameComm(mc, 1, viewer, 2);

411:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
412:   if (iascii) {
413:     PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)mc, viewer));
414:     PetscCall(PetscViewerASCIIPrintf(viewer, "  Weight type: %s\n", MatColoringWeightTypes[mc->weight_type]));
415:     if (mc->maxcolors > 0) {
416:       PetscCall(PetscViewerASCIIPrintf(viewer, "  Distance %" PetscInt_FMT ", Max. Colors %" PetscInt_FMT "\n", mc->dist, mc->maxcolors));
417:     } else {
418:       PetscCall(PetscViewerASCIIPrintf(viewer, "  Distance %" PetscInt_FMT "\n", mc->dist));
419:     }
420:   }
421:   PetscFunctionReturn(PETSC_SUCCESS);
422: }

424: /*@
425:   MatColoringSetWeightType - Set the type of weight computation used while computing the coloring

427:   Logically Collective

429:   Input Parameters:
430: + mc - the `MatColoring` context
431: - wt - the weight type

433:   Level: beginner

435: .seealso: `MatColoring`, `MatColoringWeightType`, `MatColoringApply()`
436: @*/
437: PetscErrorCode MatColoringSetWeightType(MatColoring mc, MatColoringWeightType wt)
438: {
439:   PetscFunctionBegin;
440:   mc->weight_type = wt;
441:   PetscFunctionReturn(PETSC_SUCCESS);
442: }