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