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