Actual source code: ex104.c
1: static char help[] = "Tests DMPlexCreateColoring() and DMPlexCreateColoringLabel().\n\n";
3: #include <petscdmplex.h>
4: #include <petsc/private/hashseti.h>
6: typedef struct {
7: PetscInt depth;
8: PetscInt distance;
9: PetscInt markCells;
10: PetscBool femAdjacency;
11: } AppCtx;
13: PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
14: {
15: PetscFunctionBegin;
16: options->depth = 0;
17: options->distance = 1;
18: options->markCells = 0;
19: options->femAdjacency = PETSC_FALSE;
20: PetscOptionsBegin(comm, "", "DMPlexCreateColoring() Test Options", "DMPLEX");
21: PetscCall(PetscOptionsInt("-depth", "Stratum depth defining the nodes in the connectivity graph", "ex104.c", options->depth, &options->depth, NULL));
22: PetscCall(PetscOptionsInt("-distance", "How far through the mesh a point reaches", "ex104.c", options->distance, &options->distance, NULL));
23: PetscCall(PetscOptionsInt("-mark_cells", "Color only the points in the closure of this many cells, instead of the whole stratum", "ex104.c", options->markCells, &options->markCells, NULL));
24: PetscCall(PetscOptionsBool("-fem_adjacency", "Use the finite-element adjacency at the cell stratum too, as a patch coloring needs", "ex104.c", options->femAdjacency, &options->femAdjacency, NULL));
25: PetscOptionsEnd();
26: PetscFunctionReturn(PETSC_SUCCESS);
27: }
29: /* Select the closure of the first markCells owned cells. */
30: static PetscErrorCode CreateActiveLabel(DM dm, AppCtx *user, DMLabel *label)
31: {
32: PetscInt cStart, cEnd;
34: PetscFunctionBeginUser;
35: PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd));
36: PetscCall(DMLabelCreate(PetscObjectComm((PetscObject)dm), "active", label));
37: for (PetscInt c = cStart; c < PetscMin(cStart + user->markCells, cEnd); ++c) PetscCall(DMLabelSetValue(*label, c, 1));
38: PetscCall(DMPlexLabelComplete(dm, *label));
39: PetscFunctionReturn(PETSC_SUCCESS);
40: }
42: /*
43: Check that no two points of a color are within distance hops in the mesh. Cross-rank conflicts are not checked.
44: */
45: static PetscErrorCode CheckColoring(DM dm, PetscInt distance, PetscInt ncolors, IS iscolors[])
46: {
47: PetscHSetI ht, nbr;
48: PetscInt *pts = NULL;
49: PetscInt npts = 0, index = 0;
51: PetscFunctionBeginUser;
52: PetscCall(PetscHSetICreate(&ht));
53: PetscCall(PetscHSetICreate(&nbr));
54: for (PetscInt c = 0; c < ncolors; ++c) {
55: const PetscInt *color;
56: PetscInt n;
58: PetscCall(PetscHSetIClear(ht));
59: PetscCall(ISGetLocalSize(iscolors[c], &n));
60: PetscCall(ISGetIndices(iscolors[c], &color));
61: for (PetscInt k = 0; k < n; ++k) PetscCall(PetscHSetIAdd(ht, color[k]));
62: for (PetscInt k = 0; k < n; ++k) {
63: /* Grow the neighborhood of this point one hop at a time */
64: PetscCall(PetscHSetIClear(nbr));
65: PetscCall(PetscHSetIAdd(nbr, color[k]));
66: for (PetscInt r = 0; r < distance; ++r) {
67: PetscCall(PetscHSetIGetSize(nbr, &npts));
68: PetscCall(PetscMalloc1(npts, &pts));
69: index = 0;
70: PetscCall(PetscHSetIGetElems(nbr, &index, pts));
71: for (PetscInt m = 0; m < npts; ++m) {
72: PetscInt nadj = PETSC_DETERMINE;
73: PetscInt *adj = NULL;
75: PetscCall(DMPlexGetAdjacency(dm, pts[m], &nadj, &adj));
76: for (PetscInt a = 0; a < nadj; ++a) PetscCall(PetscHSetIAdd(nbr, adj[a]));
77: PetscCall(PetscFree(adj));
78: }
79: PetscCall(PetscFree(pts));
80: }
81: PetscCall(PetscHSetIGetSize(nbr, &npts));
82: PetscCall(PetscMalloc1(npts, &pts));
83: index = 0;
84: PetscCall(PetscHSetIGetElems(nbr, &index, pts));
85: for (PetscInt m = 0; m < npts; ++m) {
86: PetscBool has;
88: if (pts[m] == color[k]) continue;
89: PetscCall(PetscHSetIHas(ht, pts[m], &has));
90: PetscCheck(has == PETSC_FALSE, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Points %" PetscInt_FMT " and %" PetscInt_FMT " are within %" PetscInt_FMT " of each other but share color %" PetscInt_FMT, color[k], pts[m], distance, c);
91: }
92: PetscCall(PetscFree(pts));
93: }
94: PetscCall(ISRestoreIndices(iscolors[c], &color));
95: }
96: PetscCall(PetscHSetIDestroy(&nbr));
97: PetscCall(PetscHSetIDestroy(&ht));
98: PetscFunctionReturn(PETSC_SUCCESS);
99: }
101: PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
102: {
103: DM pdm = NULL;
104: PetscInt overlap = user->distance;
105: PetscInt dim;
107: PetscFunctionBegin;
108: PetscCall(DMCreate(comm, dm));
109: PetscCall(DMSetType(*dm, DMPLEX));
110: PetscCall(DMPlexDistributeSetDefault(*dm, PETSC_TRUE));
111: PetscCall(DMSetFromOptions(*dm));
112: PetscCall(DMGetDimension(*dm, &dim));
113: /* A cell's finite-element adjacency does not reach another cell; use cone adjacency unless requested. */
114: if (user->depth == dim && user->femAdjacency == PETSC_FALSE) PetscCall(DMSetBasicAdjacency(*dm, PETSC_TRUE, PETSC_FALSE));
115: else PetscCall(DMSetBasicAdjacency(*dm, PETSC_FALSE, PETSC_TRUE));
116: {
117: PetscPartitioner part;
118: PetscCall(DMPlexSetOptionsPrefix(*dm, "lb_"));
119: PetscCall(DMPlexGetPartitioner(*dm, &part));
120: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)part, "lb_"));
121: PetscCall(PetscPartitionerSetFromOptions(part));
122: }
123: PetscCall(DMPlexDistribute(*dm, overlap, NULL, &pdm));
124: if (pdm) {
125: PetscCall(DMDestroy(dm));
126: *dm = pdm;
127: }
128: PetscFunctionReturn(PETSC_SUCCESS);
129: }
131: int main(int argc, char **argv)
132: {
133: DM dm;
134: DMLabel active = NULL;
135: AppCtx user;
136: PetscInt ncolors = 0, maxcolors = 0;
137: IS *iscolors = NULL;
138: ISColoring coloring = NULL;
140: PetscFunctionBeginUser;
141: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
142: PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user));
143: PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm));
144: if (user.markCells > 0) PetscCall(CreateActiveLabel(dm, &user, &active));
145: if (active == NULL) PetscCall(DMPlexCreateColoring(dm, user.depth, user.distance, &coloring));
146: else PetscCall(DMPlexCreateColoringLabel(dm, user.depth, user.distance, active, 1, &coloring));
147: PetscCall(ISColoringGetIS(coloring, PETSC_USE_POINTER, &ncolors, &iscolors));
148: /* Report the largest color count across processes. */
149: PetscCallMPI(MPIU_Allreduce(&ncolors, &maxcolors, 1, MPIU_INT, MPI_MAX, PETSC_COMM_WORLD));
150: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Number of colors: %" PetscInt_FMT "\n", maxcolors));
151: PetscCall(CheckColoring(dm, user.distance, ncolors, iscolors));
152: for (PetscInt c = 0; c < ncolors; c++) {
153: PetscCall(ISViewFromOptions(iscolors[c], NULL, "-iscoloring_view"));
154: }
155: PetscCall(ISColoringRestoreIS(coloring, PETSC_USE_POINTER, &iscolors));
156: PetscCall(ISColoringDestroy(&coloring));
157: PetscCall(DMLabelDestroy(&active));
158: PetscCall(DMDestroy(&dm));
159: PetscCall(PetscFinalize());
160: return 0;
161: }
163: /*TEST
165: test:
166: nsize: {{1 2}separate output}
167: args: -depth {{0 1 2}separate output} -distance {{1 2}separate output} -iscoloring_view -dm_coord_space 0 -dm_plex_simplex 0 -dm_plex_box_faces 4,4 -petscpartitioner_type simple
169: # Lexical weighting uses the four-color vertex pattern; the other cases exercise fixed and reordered point order.
170: testset:
171: nsize: 1
172: args: -depth 0 -distance 1 -dm_coord_space 0 -dm_plex_simplex 0 -dm_plex_box_faces 16,16
173: test:
174: suffix: grid_lexical
175: output_file: output/ex104_grid.out
176: test:
177: suffix: grid_natural
178: args: -dm_plex_coloring_ordering_type natural
179: output_file: output/ex104_grid.out
180: test:
181: suffix: grid_ordering
182: args: -dm_plex_coloring_ordering_type {{rcm nd}separate output}
184: # Local coloring tests the induced graph on each rank; the reported count is the maximum across ranks.
185: test:
186: suffix: local
187: nsize: {{1 2}separate output}
188: args: -depth 0 -distance 1 -dm_plex_coloring_local -dm_coord_space 0 -dm_plex_simplex 0 -dm_plex_box_faces 4,4 -petscpartitioner_type simple
190: # Color only the closure of a few cells.
191: test:
192: suffix: label
193: nsize: {{1 2}separate output}
194: args: -depth 0 -distance 1 -mark_cells 3 -iscoloring_view -dm_coord_space 0 -dm_plex_simplex 0 -dm_plex_box_faces 8,8 -petscpartitioner_type simple
196: # Finite-element cell adjacency gives an edgeless distance-one graph; distance two tests Vanka patches.
197: test:
198: suffix: cell_fem
199: nsize: {{1 2}separate output}
200: args: -depth 2 -fem_adjacency -distance {{1 2}separate output} -dm_coord_space 0 -dm_plex_simplex 0 -dm_plex_box_faces 4,4 -petscpartitioner_type simple
202: TEST*/