1 static char help[] = "Test for mesh reordering\n\n"; 2 3 #include <petscdmplex.h> 4 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; 13 14 PetscErrorCode ProcessOptions(AppCtx *options) 15 { 16 PetscInt len; 17 PetscBool flg; 18 19 PetscFunctionBegin; 20 options->numFields = 1; 21 options->numComponents = NULL; 22 options->numDof = NULL; 23 options->numGroups = 0; 24 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(0); 36 } 37 38 PetscErrorCode CleanupContext(AppCtx *user) 39 { 40 PetscFunctionBegin; 41 PetscCall(PetscFree(user->numComponents)); 42 PetscCall(PetscFree(user->numDof)); 43 PetscFunctionReturn(0); 44 } 45 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, 50 2, 9, 7, 10, 9, 4, 1, 10, 12, 10, 4, 11, 12, 11, 3, 3, 11, 14, 11, 4, 13, 14, 13, 5}; 51 const PetscReal coords[15*2] = {0, -3, 0, -1, 2, -1, 0, 1, 2, 1, 52 0, 3, 1, -2, 1, -1, 0, -2, 2, 0, 53 1, 0, 1, 1, 0, 0, 1, 2, 0, 2}; 54 55 PetscFunctionBegin; 56 PetscCall(DMPlexCreateFromCellListPetsc(comm, 2, 16, 15, 3, PETSC_FALSE, cells, 2, coords, dm)); 57 PetscFunctionReturn(0); 58 } 59 60 PetscErrorCode TestReordering(DM dm, AppCtx *user) 61 { 62 DM pdm; 63 IS perm; 64 Mat A, pA; 65 PetscInt bw, pbw; 66 MatOrderingType order = MATORDERINGRCM; 67 68 PetscFunctionBegin; 69 PetscCall(DMPlexGetOrdering(dm, order, NULL, &perm)); 70 PetscCall(DMPlexPermute(dm, perm, &pdm)); 71 PetscCall(PetscObjectSetOptionsPrefix((PetscObject) pdm, "perm_")); 72 PetscCall(DMSetFromOptions(pdm)); 73 PetscCall(ISDestroy(&perm)); 74 PetscCall(DMViewFromOptions(dm, NULL, "-orig_dm_view")); 75 PetscCall(DMViewFromOptions(pdm, NULL, "-dm_view")); 76 PetscCall(DMCreateMatrix(dm, &A)); 77 PetscCall(DMCreateMatrix(pdm, &pA)); 78 PetscCall(MatComputeBandwidth(A, 0.0, &bw)); 79 PetscCall(MatComputeBandwidth(pA, 0.0, &pbw)); 80 PetscCall(MatViewFromOptions(A, NULL, "-orig_mat_view")); 81 PetscCall(MatViewFromOptions(pA, NULL, "-perm_mat_view")); 82 PetscCall(MatDestroy(&A)); 83 PetscCall(MatDestroy(&pA)); 84 PetscCall(DMDestroy(&pdm)); 85 if (pbw > bw) { 86 PetscCall(PetscPrintf(PetscObjectComm((PetscObject) dm), "Ordering method %s increased bandwidth from %" PetscInt_FMT " to %" PetscInt_FMT "\n", order, bw, pbw)); 87 } else { 88 PetscCall(PetscPrintf(PetscObjectComm((PetscObject) dm), "Ordering method %s reduced bandwidth from %" PetscInt_FMT " to %" PetscInt_FMT "\n", order, bw, pbw)); 89 } 90 PetscFunctionReturn(0); 91 } 92 93 PetscErrorCode CreateGroupLabel(DM dm, PetscInt numGroups, DMLabel *label, AppCtx *options) 94 { 95 const PetscInt groupA[10] = {15, 3, 13, 12, 2, 10, 7, 6, 0, 4}; 96 const PetscInt groupB[6] = {14, 11, 9, 1, 8, 5}; 97 PetscInt c; 98 99 PetscFunctionBegin; 100 if (numGroups < 2) {*label = NULL; PetscFunctionReturn(0);} 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(0); 106 } 107 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; 115 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(0); 134 } 135 136 int main(int argc, char **argv) 137 { 138 DM dm; 139 PetscSection s; 140 AppCtx user; 141 PetscInt dim; 142 143 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 144 PetscCall(ProcessOptions(&user)); 145 if (user.numGroups < 1) { 146 PetscCall(DMCreate(PETSC_COMM_WORLD, &dm)); 147 PetscCall(DMSetType(dm, DMPLEX)); 148 } else { 149 PetscCall(CreateTestMesh(PETSC_COMM_WORLD, &dm, &user)); 150 } 151 PetscCall(DMSetFromOptions(dm)); 152 PetscCall(DMViewFromOptions(dm, NULL, "-dm_view")); 153 PetscCall(DMGetDimension(dm, &dim)); 154 { 155 PetscInt len = (dim+1) * PetscMax(1, user.numFields); 156 PetscBool flg; 157 158 PetscCall(PetscCalloc1(len, &user.numDof)); 159 PetscOptionsBegin(PETSC_COMM_SELF, "", "Meshing Problem Options", "DMPLEX"); 160 PetscCall(PetscOptionsIntArray("-num_dof", "The dof signature for the section", "ex10.c", user.numDof, &len, &flg)); 161 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)); 162 PetscOptionsEnd(); 163 } 164 if (user.numGroups < 1) { 165 PetscCall(DMSetNumFields(dm, user.numFields)); 166 PetscCall(DMCreateDS(dm)); 167 PetscCall(DMPlexCreateSection(dm, NULL, user.numComponents, user.numDof, 0, NULL, NULL, NULL, NULL, &s)); 168 PetscCall(DMSetLocalSection(dm, s)); 169 PetscCall(PetscSectionDestroy(&s)); 170 PetscCall(TestReordering(dm, &user)); 171 } else { 172 PetscCall(DMSetNumFields(dm, user.numFields)); 173 PetscCall(DMCreateDS(dm)); 174 PetscCall(DMPlexCreateSection(dm, NULL, user.numComponents, user.numDof, 0, NULL, NULL, NULL, NULL, &s)); 175 PetscCall(DMSetLocalSection(dm, s)); 176 PetscCall(PetscSectionDestroy(&s)); 177 PetscCall(TestReorderingByGroup(dm, &user)); 178 } 179 PetscCall(DMDestroy(&dm)); 180 PetscCall(CleanupContext(&user)); 181 PetscCall(PetscFinalize()); 182 return 0; 183 } 184 185 /*TEST 186 187 # Two cell tests 0-3 188 test: 189 suffix: 0 190 requires: triangle 191 args: -dm_plex_simplex 1 -num_dof 1,0,0 -mat_view -dm_coord_space 0 192 test: 193 suffix: 1 194 args: -dm_plex_simplex 0 -num_dof 1,0,0 -mat_view -dm_coord_space 0 195 test: 196 suffix: 2 197 requires: ctetgen 198 args: -dm_plex_dim 3 -dm_plex_simplex 1 -num_dof 1,0,0,0 -mat_view -dm_coord_space 0 199 test: 200 suffix: 3 201 args: -dm_plex_dim 3 -dm_plex_simplex 0 -num_dof 1,0,0,0 -mat_view -dm_coord_space 0 202 # Refined tests 4-7 203 test: 204 suffix: 4 205 requires: triangle 206 args: -dm_plex_simplex 1 -dm_refine_volume_limit_pre 0.00625 -num_dof 1,0,0 207 test: 208 suffix: 5 209 args: -dm_plex_simplex 0 -dm_refine 1 -num_dof 1,0,0 210 test: 211 suffix: 6 212 requires: ctetgen 213 args: -dm_plex_dim 3 -dm_plex_simplex 1 -dm_refine_volume_limit_pre 0.00625 -num_dof 1,0,0,0 214 test: 215 suffix: 7 216 args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_refine 1 -num_dof 1,0,0,0 217 # Parallel tests 218 # Grouping tests 219 test: 220 suffix: group_1 221 args: -num_groups 1 -num_dof 1,0,0 -is_view -orig_mat_view -perm_mat_view 222 test: 223 suffix: group_2 224 args: -num_groups 2 -num_dof 1,0,0 -is_view -perm_mat_view 225 226 TEST*/ 227