Actual source code: ex10.c
1: static char help[] = "Test for mesh reordering\n\n";
3: #include <petscdmplex.h>
5: typedef struct {
6: PetscInt dim; /* The topological mesh dimension */
7: PetscReal refinementLimit; /* Maximum volume of a refined cell */
8: PetscInt numFields; /* The number of section fields */
9: PetscInt *numComponents; /* The number of field components */
10: PetscInt *numDof; /* The dof signature for the section */
11: PetscInt numGroups; /* If greater than 1, use grouping in test */
12: } AppCtx;
14: PetscErrorCode ProcessOptions(AppCtx *options)
15: {
16: PetscInt len;
17: PetscBool flg;
19: PetscFunctionBegin;
20: options->numFields = 1;
21: options->numComponents = NULL;
22: options->numDof = NULL;
23: options->numGroups = 0;
25: PetscOptionsBegin(PETSC_COMM_SELF, "", "Meshing Problem Options", "DMPLEX");
26: PetscCall(PetscOptionsBoundedInt("-num_fields", "The number of section fields", "ex10.c", options->numFields, &options->numFields, NULL, 1));
27: if (options->numFields) {
28: len = options->numFields;
29: PetscCall(PetscCalloc1(len, &options->numComponents));
30: PetscCall(PetscOptionsIntArray("-num_components", "The number of components per field", "ex10.c", options->numComponents, &len, &flg));
31: PetscCheck(!flg || !(len != options->numFields), PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Length of components array is %" PetscInt_FMT " should be %" PetscInt_FMT, len, options->numFields);
32: }
33: PetscCall(PetscOptionsBoundedInt("-num_groups", "Group permutation by this many label values", "ex10.c", options->numGroups, &options->numGroups, NULL, 0));
34: PetscOptionsEnd();
35: PetscFunctionReturn(PETSC_SUCCESS);
36: }
38: PetscErrorCode CleanupContext(AppCtx *user)
39: {
40: PetscFunctionBegin;
41: PetscCall(PetscFree(user->numComponents));
42: PetscCall(PetscFree(user->numDof));
43: PetscFunctionReturn(PETSC_SUCCESS);
44: }
46: /* This mesh comes from~\cite{saad2003}, Fig. 2.10, p. 70. */
47: PetscErrorCode CreateTestMesh(MPI_Comm comm, DM *dm, AppCtx *options)
48: {
49: const PetscInt cells[16 * 3] = {6, 7, 8, 7, 9, 10, 10, 11, 12, 11, 13, 14, 0, 6, 8, 6, 2, 7, 1, 8, 7, 1, 7, 10, 2, 9, 7, 10, 9, 4, 1, 10, 12, 10, 4, 11, 12, 11, 3, 3, 11, 14, 11, 4, 13, 14, 13, 5};
50: const PetscReal coords[15 * 2] = {0, -3, 0, -1, 2, -1, 0, 1, 2, 1, 0, 3, 1, -2, 1, -1, 0, -2, 2, 0, 1, 0, 1, 1, 0, 0, 1, 2, 0, 2};
52: PetscFunctionBegin;
53: PetscCall(DMPlexCreateFromCellListPetsc(comm, 2, 16, 15, 3, PETSC_FALSE, cells, 2, coords, dm));
54: PetscFunctionReturn(PETSC_SUCCESS);
55: }
57: PetscErrorCode TestReordering(DM dm, AppCtx *user)
58: {
59: DM pdm;
60: IS perm;
61: Mat A, pA;
62: PetscInt bw, pbw;
63: MatOrderingType order = MATORDERINGRCM;
65: PetscFunctionBegin;
66: PetscCall(DMPlexGetOrdering(dm, order, NULL, &perm));
67: PetscCall(DMPlexPermute(dm, perm, &pdm));
68: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)pdm, "perm_"));
69: PetscCall(DMSetFromOptions(pdm));
70: PetscCall(ISDestroy(&perm));
71: PetscCall(DMViewFromOptions(dm, NULL, "-orig_dm_view"));
72: PetscCall(DMViewFromOptions(pdm, NULL, "-dm_view"));
73: PetscCall(DMCreateMatrix(dm, &A));
74: PetscCall(DMCreateMatrix(pdm, &pA));
75: PetscCall(MatComputeBandwidth(A, 0.0, &bw));
76: PetscCall(MatComputeBandwidth(pA, 0.0, &pbw));
77: PetscCall(MatViewFromOptions(A, NULL, "-orig_mat_view"));
78: PetscCall(MatViewFromOptions(pA, NULL, "-perm_mat_view"));
79: PetscCall(MatDestroy(&A));
80: PetscCall(MatDestroy(&pA));
81: PetscCall(DMDestroy(&pdm));
82: if (pbw > bw) {
83: PetscCall(PetscPrintf(PetscObjectComm((PetscObject)dm), "Ordering method %s increased bandwidth from %" PetscInt_FMT " to %" PetscInt_FMT "\n", order, bw, pbw));
84: } else {
85: PetscCall(PetscPrintf(PetscObjectComm((PetscObject)dm), "Ordering method %s reduced bandwidth from %" PetscInt_FMT " to %" PetscInt_FMT "\n", order, bw, pbw));
86: }
87: PetscFunctionReturn(PETSC_SUCCESS);
88: }
90: PetscErrorCode CreateGroupLabel(DM dm, PetscInt numGroups, DMLabel *label, AppCtx *options)
91: {
92: const PetscInt groupA[10] = {15, 3, 13, 12, 2, 10, 7, 6, 0, 4};
93: const PetscInt groupB[6] = {14, 11, 9, 1, 8, 5};
94: PetscInt c;
96: PetscFunctionBegin;
97: if (numGroups < 2) {
98: *label = NULL;
99: PetscFunctionReturn(PETSC_SUCCESS);
100: }
101: PetscCheck(numGroups == 2, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Test only coded for 2 groups, not %" PetscInt_FMT, numGroups);
102: PetscCall(DMLabelCreate(PETSC_COMM_SELF, "groups", label));
103: for (c = 0; c < 10; ++c) PetscCall(DMLabelSetValue(*label, groupA[c], 101));
104: for (c = 0; c < 6; ++c) PetscCall(DMLabelSetValue(*label, groupB[c], 1001));
105: PetscFunctionReturn(PETSC_SUCCESS);
106: }
108: PetscErrorCode TestReorderingByGroup(DM dm, AppCtx *user)
109: {
110: DM pdm;
111: DMLabel label;
112: Mat A, pA;
113: MatOrderingType order = MATORDERINGRCM;
114: IS perm;
116: PetscFunctionBegin;
117: PetscCall(CreateGroupLabel(dm, user->numGroups, &label, user));
118: PetscCall(DMPlexGetOrdering(dm, order, label, &perm));
119: PetscCall(DMLabelDestroy(&label));
120: PetscCall(DMPlexPermute(dm, perm, &pdm));
121: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)pdm, "perm_"));
122: PetscCall(DMSetFromOptions(pdm));
123: PetscCall(DMViewFromOptions(dm, NULL, "-orig_dm_view"));
124: PetscCall(DMViewFromOptions(pdm, NULL, "-perm_dm_view"));
125: PetscCall(ISDestroy(&perm));
126: PetscCall(DMCreateMatrix(dm, &A));
127: PetscCall(DMCreateMatrix(pdm, &pA));
128: PetscCall(MatViewFromOptions(A, NULL, "-orig_mat_view"));
129: PetscCall(MatViewFromOptions(pA, NULL, "-perm_mat_view"));
130: PetscCall(MatDestroy(&A));
131: PetscCall(MatDestroy(&pA));
132: PetscCall(DMDestroy(&pdm));
133: PetscFunctionReturn(PETSC_SUCCESS);
134: }
136: int main(int argc, char **argv)
137: {
138: DM dm;
139: PetscSection s;
140: AppCtx user;
141: PetscInt dim;
143: PetscFunctionBeginUser;
144: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
145: PetscCall(ProcessOptions(&user));
146: if (user.numGroups < 1) {
147: PetscCall(DMCreate(PETSC_COMM_WORLD, &dm));
148: PetscCall(DMSetType(dm, DMPLEX));
149: } else {
150: PetscCall(CreateTestMesh(PETSC_COMM_WORLD, &dm, &user));
151: }
152: PetscCall(DMSetFromOptions(dm));
153: PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
154: PetscCall(DMGetDimension(dm, &dim));
155: {
156: PetscInt len = (dim + 1) * PetscMax(1, user.numFields);
157: PetscBool flg;
159: PetscCall(PetscCalloc1(len, &user.numDof));
160: PetscOptionsBegin(PETSC_COMM_SELF, "", "Meshing Problem Options", "DMPLEX");
161: PetscCall(PetscOptionsIntArray("-num_dof", "The dof signature for the section", "ex10.c", user.numDof, &len, &flg));
162: if (flg) PetscCheck(len == ((dim + 1) * PetscMax(1, user.numFields)), PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Length of dof array is %" PetscInt_FMT " should be %" PetscInt_FMT, len, (dim + 1) * PetscMax(1, user.numFields));
163: PetscOptionsEnd();
164: }
165: if (user.numGroups < 1) {
166: PetscCall(DMSetNumFields(dm, user.numFields));
167: PetscCall(DMCreateDS(dm));
168: PetscCall(DMPlexCreateSection(dm, NULL, user.numComponents, user.numDof, 0, NULL, NULL, NULL, NULL, &s));
169: PetscCall(DMSetLocalSection(dm, s));
170: PetscCall(PetscSectionDestroy(&s));
171: PetscCall(TestReordering(dm, &user));
172: } else {
173: PetscCall(DMSetNumFields(dm, user.numFields));
174: PetscCall(DMCreateDS(dm));
175: PetscCall(DMPlexCreateSection(dm, NULL, user.numComponents, user.numDof, 0, NULL, NULL, NULL, NULL, &s));
176: PetscCall(DMSetLocalSection(dm, s));
177: PetscCall(PetscSectionDestroy(&s));
178: PetscCall(TestReorderingByGroup(dm, &user));
179: }
180: PetscCall(DMDestroy(&dm));
181: PetscCall(CleanupContext(&user));
182: PetscCall(PetscFinalize());
183: return 0;
184: }
186: /*TEST
188: # Two cell tests 0-3
189: test:
190: suffix: 0
191: requires: triangle
192: args: -dm_plex_simplex 1 -num_dof 1,0,0 -mat_view -dm_coord_space 0
193: test:
194: suffix: 1
195: args: -dm_plex_simplex 0 -num_dof 1,0,0 -mat_view -dm_coord_space 0
196: test:
197: suffix: 2
198: requires: ctetgen
199: args: -dm_plex_dim 3 -dm_plex_simplex 1 -num_dof 1,0,0,0 -mat_view -dm_coord_space 0
200: test:
201: suffix: 3
202: args: -dm_plex_dim 3 -dm_plex_simplex 0 -num_dof 1,0,0,0 -mat_view -dm_coord_space 0
203: # Refined tests 4-7
204: test:
205: suffix: 4
206: requires: triangle
207: args: -dm_plex_simplex 1 -dm_refine_volume_limit_pre 0.00625 -num_dof 1,0,0
208: test:
209: suffix: 5
210: args: -dm_plex_simplex 0 -dm_refine 1 -num_dof 1,0,0
211: test:
212: suffix: 6
213: requires: ctetgen
214: args: -dm_plex_dim 3 -dm_plex_simplex 1 -dm_refine_volume_limit_pre 0.00625 -num_dof 1,0,0,0
215: test:
216: suffix: 7
217: args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_refine 1 -num_dof 1,0,0,0
218: # Parallel tests
219: # Grouping tests
220: test:
221: suffix: group_1
222: args: -num_groups 1 -num_dof 1,0,0 -is_view -orig_mat_view -perm_mat_view
223: test:
224: suffix: group_2
225: args: -num_groups 2 -num_dof 1,0,0 -is_view -perm_mat_view
227: TEST*/