1 static char help[] = "Tests for refinement of meshes created by hand\n\n"; 2 3 #include <petscdmplex.h> 4 5 typedef struct { 6 PetscInt debug; /* The debugging level */ 7 PetscInt dim; /* The topological mesh dimension */ 8 PetscBool cellHybrid; /* Use a hybrid mesh */ 9 PetscBool cellSimplex; /* Use simplices or hexes */ 10 PetscBool testPartition; /* Use a fixed partitioning for testing */ 11 PetscInt testNum; /* The particular mesh to test */ 12 PetscBool uninterpolate; /* Uninterpolate the mesh at the end */ 13 PetscBool reinterpolate; /* Reinterpolate the mesh at the end */ 14 } AppCtx; 15 16 PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 17 { 18 PetscErrorCode ierr; 19 20 PetscFunctionBegin; 21 options->debug = 0; 22 options->dim = 2; 23 options->cellHybrid = PETSC_TRUE; 24 options->cellSimplex = PETSC_TRUE; 25 options->testPartition = PETSC_TRUE; 26 options->testNum = 0; 27 options->uninterpolate = PETSC_FALSE; 28 options->reinterpolate = PETSC_FALSE; 29 30 ierr = PetscOptionsBegin(comm, "", "Meshing Problem Options", "DMPLEX");CHKERRQ(ierr); 31 ierr = PetscOptionsBoundedInt("-debug", "The debugging level", "ex4.c", options->debug, &options->debug, NULL,0);CHKERRQ(ierr); 32 ierr = PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex4.c", options->dim, &options->dim, NULL,1,3);CHKERRQ(ierr); 33 ierr = PetscOptionsBool("-cell_hybrid", "Use a hybrid mesh", "ex4.c", options->cellHybrid, &options->cellHybrid, NULL);CHKERRQ(ierr); 34 ierr = PetscOptionsBool("-cell_simplex", "Use simplices if true, otherwise hexes", "ex4.c", options->cellSimplex, &options->cellSimplex, NULL);CHKERRQ(ierr); 35 ierr = PetscOptionsBool("-test_partition", "Use a fixed partition for testing", "ex4.c", options->testPartition, &options->testPartition, NULL);CHKERRQ(ierr); 36 ierr = PetscOptionsBoundedInt("-test_num", "The particular mesh to test", "ex4.c", options->testNum, &options->testNum, NULL,0);CHKERRQ(ierr); 37 ierr = PetscOptionsBool("-uninterpolate", "Uninterpolate the mesh at the end", "ex4.c", options->uninterpolate, &options->uninterpolate, NULL);CHKERRQ(ierr); 38 ierr = PetscOptionsBool("-reinterpolate", "Reinterpolate the mesh at the end", "ex4.c", options->reinterpolate, &options->reinterpolate, NULL);CHKERRQ(ierr); 39 ierr = PetscOptionsEnd(); 40 PetscFunctionReturn(0); 41 } 42 43 /* Two segments 44 45 2-------0-------3-------1-------4 46 47 become 48 49 4---0---7---1---5---2---8---3---6 50 51 */ 52 PetscErrorCode CreateSimplex_1D(MPI_Comm comm, DM *dm) 53 { 54 PetscInt depth = 1; 55 PetscMPIInt rank; 56 PetscErrorCode ierr; 57 58 PetscFunctionBegin; 59 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 60 if (!rank) { 61 PetscInt numPoints[2] = {3, 2}; 62 PetscInt coneSize[5] = {2, 2, 0, 0, 0}; 63 PetscInt cones[4] = {2, 3, 3, 4}; 64 PetscInt coneOrientations[16] = {0, 0, 0, 0}; 65 PetscScalar vertexCoords[3] = {-1.0, 0.0, 1.0}; 66 67 ierr = DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 68 } else { 69 PetscInt numPoints[2] = {0, 0}; 70 71 ierr = DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr); 72 } 73 PetscFunctionReturn(0); 74 } 75 76 77 /* Two triangles 78 4 79 / | \ 80 8 | 10 81 / | \ 82 2 0 7 1 5 83 \ | / 84 6 | 9 85 \ | / 86 3 87 88 Becomes 89 10 90 / | \ 91 21 | 26 92 / | \ 93 14 2 20 4 16 94 /|\ | /|\ 95 22 | 28 | 32 | 25 96 / | \ | / | 6\ 97 8 29 3 13 7 31 11 98 \0 | / | \ | / 99 17 | 27 | 30 | 24 100 \|/ | \|/ 101 12 1 19 5 15 102 \ | / 103 18 | 23 104 \ | / 105 9 106 */ 107 PetscErrorCode CreateSimplex_2D(MPI_Comm comm, DM *dm) 108 { 109 PetscInt depth = 2; 110 PetscMPIInt rank; 111 PetscErrorCode ierr; 112 113 PetscFunctionBegin; 114 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 115 if (!rank) { 116 PetscInt numPoints[3] = {4, 5, 2}; 117 PetscInt coneSize[11] = {3, 3, 0, 0, 0, 0, 2, 2, 2, 2, 2}; 118 PetscInt cones[16] = {6, 7, 8, 7, 9, 10, 2, 3, 3, 4, 4, 2, 3, 5, 5, 4}; 119 PetscInt coneOrientations[16] = {0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 120 PetscScalar vertexCoords[8] = {-0.5, 0.0, 0.0, -0.5, 0.0, 0.5, 0.5, 0.0}; 121 122 ierr = DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 123 } else { 124 PetscInt numPoints[3] = {0, 0, 0}; 125 126 ierr = DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr); 127 } 128 PetscFunctionReturn(0); 129 } 130 131 /* Two triangles separated by a zero-volume cell with 4 vertices/2 edges 132 5--16--8 133 / | | \ 134 11 | | 12 135 / | | \ 136 3 0 10 2 14 1 6 137 \ | | / 138 9 | | 13 139 \ | | / 140 4--15--7 141 */ 142 PetscErrorCode CreateSimplexHybrid_2D(MPI_Comm comm, PetscInt testNum, DM *dm) 143 { 144 DM idm, hdm = NULL; 145 DMLabel faultLabel, hybridLabel; 146 PetscInt p; 147 PetscMPIInt rank; 148 PetscErrorCode ierr; 149 150 PetscFunctionBegin; 151 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 152 if (!rank) { 153 switch (testNum) { 154 case 0: 155 { 156 PetscInt numPoints[2] = {4, 2}; 157 PetscInt coneSize[6] = {3, 3, 0, 0, 0, 0}; 158 PetscInt cones[6] = {2, 3, 4, 5, 4, 3}; 159 PetscInt coneOrientations[6] = {0, 0, 0, 0, 0, 0}; 160 PetscScalar vertexCoords[8] = {-1.0, -0.5, 0.0, -0.5, 0.0, 0.5, 1.0, 0.5}; 161 PetscInt faultPoints[2] = {3, 4}; 162 163 ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 164 for (p = 0; p < 2; ++p) {ierr = DMSetLabelValue(*dm, "fault", faultPoints[p], 1);CHKERRQ(ierr);} 165 } 166 break; 167 case 1: 168 { 169 PetscInt numPoints[2] = {5, 4}; 170 PetscInt coneSize[9] = {3, 3, 3, 3, 0, 0, 0, 0, 0}; 171 PetscInt cones[12] = {4, 5, 6, 6, 7, 4, 6, 5, 8, 6, 8, 7}; 172 PetscInt coneOrientations[12] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 173 PetscScalar vertexCoords[10] = {-1.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0, 1.0, 1.0, 0.0}; 174 PetscInt faultPoints[3] = {5, 6, 7}; 175 176 ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 177 for (p = 0; p < 3; ++p) {ierr = DMSetLabelValue(*dm, "fault", faultPoints[p], 1);CHKERRQ(ierr);} 178 } 179 break; 180 default: SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum); 181 } 182 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 183 ierr = PetscObjectSetOptionsPrefix((PetscObject) idm, "in_");CHKERRQ(ierr); 184 ierr = DMSetFromOptions(idm);CHKERRQ(ierr); 185 ierr = DMViewFromOptions(idm, NULL, "-dm_view");CHKERRQ(ierr); 186 ierr = DMGetLabel(*dm, "fault", &faultLabel);CHKERRQ(ierr); 187 ierr = DMPlexCreateHybridMesh(idm, faultLabel, NULL, &hybridLabel, NULL, NULL, &hdm);CHKERRQ(ierr); 188 ierr = DMLabelDestroy(&hybridLabel);CHKERRQ(ierr); 189 ierr = DMDestroy(&idm);CHKERRQ(ierr); 190 ierr = DMDestroy(dm);CHKERRQ(ierr); 191 *dm = hdm; 192 } else { 193 PetscInt numPoints[2] = {0, 0}; 194 195 ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr); 196 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 197 ierr = PetscObjectSetOptionsPrefix((PetscObject) idm, "in_");CHKERRQ(ierr); 198 ierr = DMSetFromOptions(idm);CHKERRQ(ierr); 199 ierr = DMViewFromOptions(idm, NULL, "-dm_view");CHKERRQ(ierr); 200 ierr = DMPlexCreateHybridMesh(idm, NULL, NULL, NULL, NULL, NULL, &hdm);CHKERRQ(ierr); 201 ierr = DMDestroy(&idm);CHKERRQ(ierr); 202 ierr = DMDestroy(dm);CHKERRQ(ierr); 203 *dm = hdm; 204 } 205 PetscFunctionReturn(0); 206 } 207 208 /* Two quadrilaterals 209 210 5----10-----4----14-----7 211 | | | 212 | | | 213 | | | 214 11 0 9 1 13 215 | | | 216 | | | 217 | | | 218 2-----8-----3----12-----6 219 */ 220 PetscErrorCode CreateTensorProduct_2D(MPI_Comm comm, PetscInt testNum, DM *dm) 221 { 222 PetscInt depth = 2; 223 PetscMPIInt rank; 224 PetscErrorCode ierr; 225 226 PetscFunctionBegin; 227 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 228 if (!rank) { 229 PetscInt numPoints[3] = {6, 7, 2}; 230 PetscInt coneSize[15] = {4, 4, 0, 0, 0, 0, 0, 0, 2, 2, 2, 2, 2, 2, 2}; 231 PetscInt cones[22] = {8, 9, 10, 11, 12, 13, 14, 9, 2, 3, 3, 4, 4, 5, 5, 2, 3, 6, 6, 7, 7, 4}; 232 PetscInt coneOrientations[22] = {0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 233 PetscScalar vertexCoords[12] = {-1.0, -0.5, 0.0, -0.5, 0.0, 0.5, -1.0, 0.5, 1.0, -0.5, 1.0, 0.5}; 234 235 ierr = DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 236 } else { 237 PetscInt numPoints[3] = {0, 0, 0}; 238 239 ierr = DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr); 240 } 241 PetscFunctionReturn(0); 242 } 243 244 PetscErrorCode CreateTensorProductHybrid_2D(MPI_Comm comm, PetscInt testNum, DM *dm) 245 { 246 DM idm, hdm = NULL; 247 DMLabel faultLabel, hybridLabel; 248 PetscInt p; 249 PetscMPIInt rank; 250 PetscErrorCode ierr; 251 252 PetscFunctionBegin; 253 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 254 if (!rank) { 255 PetscInt numPoints[2] = {6, 2}; 256 PetscInt coneSize[8] = {4, 4, 0, 0, 0, 0, 0, 0}; 257 PetscInt cones[8] = {2, 3, 4, 5, 3, 6, 7, 4,}; 258 PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0}; 259 PetscScalar vertexCoords[12] = {-1.0, -0.5, 0.0, -0.5, 0.0, 0.5, -1.0, 0.5, 1.0, -0.5, 1.0, 0.5}; 260 PetscInt faultPoints[2] = {3, 4}; 261 262 ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 263 for (p = 0; p < 2; ++p) {ierr = DMSetLabelValue(*dm, "fault", faultPoints[p], 1);CHKERRQ(ierr);} 264 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 265 ierr = PetscObjectSetOptionsPrefix((PetscObject) idm, "in_");CHKERRQ(ierr); 266 ierr = DMSetFromOptions(idm);CHKERRQ(ierr); 267 ierr = DMViewFromOptions(idm, NULL, "-dm_view");CHKERRQ(ierr); 268 ierr = DMGetLabel(*dm, "fault", &faultLabel);CHKERRQ(ierr); 269 ierr = DMPlexCreateHybridMesh(idm, faultLabel, NULL, &hybridLabel, NULL, NULL, &hdm);CHKERRQ(ierr); 270 ierr = DMLabelDestroy(&hybridLabel);CHKERRQ(ierr); 271 } else { 272 PetscInt numPoints[3] = {0, 0, 0}; 273 274 ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr); 275 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 276 ierr = PetscObjectSetOptionsPrefix((PetscObject) idm, "in_");CHKERRQ(ierr); 277 ierr = DMSetFromOptions(idm);CHKERRQ(ierr); 278 ierr = DMViewFromOptions(idm, NULL, "-dm_view");CHKERRQ(ierr); 279 ierr = DMPlexCreateHybridMesh(idm, NULL, NULL, NULL, NULL, NULL, &hdm);CHKERRQ(ierr); 280 } 281 ierr = DMDestroy(&idm);CHKERRQ(ierr); 282 ierr = DMDestroy(dm);CHKERRQ(ierr); 283 *dm = hdm; 284 PetscFunctionReturn(0); 285 } 286 287 /* Two tetrahedrons 288 289 cell 5 5______ cell 290 0 / | \ |\ \ 1 291 17 | 18 | 18 13 21 292 /8 19 10\ 19 \ \ 293 2-14-|----4 | 4--22--6 294 \ 9 | 7 / |10 / / 295 16 | 15 | 15 12 20 296 \ | / |/ / 297 3 3------ 298 */ 299 PetscErrorCode CreateSimplex_3D(MPI_Comm comm, PetscInt testNum, DM *dm) 300 { 301 DM idm; 302 PetscInt depth = 3; 303 PetscMPIInt rank; 304 PetscErrorCode ierr; 305 306 PetscFunctionBegin; 307 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 308 if (!rank) { 309 switch (testNum) { 310 case 0: 311 { 312 PetscInt numPoints[4] = {5, 9, 7, 2}; 313 PetscInt coneSize[23] = {4, 4, 0, 0, 0, 0, 0, 3, 3, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2}; 314 PetscInt cones[47] = { 7, 8, 9, 10, 10, 11, 12, 13, 14, 15, 16, 17, 18, 14, 16, 19, 17, 15, 18, 19, 20, 21, 19, 15, 22, 20, 18, 21, 22, 2, 4, 4, 3, 3, 2, 2, 5, 5, 4, 3, 5, 3, 6, 6, 5, 4, 6}; 315 PetscInt coneOrientations[47] = { 0, 0, 0, 0, -3, 0, 0, 0, 0, 0, 0, 0, 0, -2, -2, 0, -2, -2, -2, -2, 0, 0, -2, -2, 0, -2, -2, -2, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 316 PetscScalar vertexCoords[15] = {0.0, 0.0, -0.5, 0.0, -0.5, 0.0, 1.0, 0.0, 0.0, 0.0, 0.5, 0.0, 0.0, 0.0, 0.5}; 317 318 ierr = DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 319 } 320 break; 321 case 1: 322 { 323 PetscInt numPoints[2] = {5, 2}; 324 PetscInt coneSize[7] = {4, 4, 0, 0, 0, 0, 0}; 325 PetscInt cones[8] = {4, 3, 5, 2, 5, 3, 4, 6}; 326 PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0}; 327 PetscScalar vertexCoords[15] = {-1.0, 0.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0, 0.0, 0.0}; 328 329 depth = 1; 330 ierr = DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 331 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 332 ierr = DMViewFromOptions(idm, NULL, "-in_dm_view");CHKERRQ(ierr); 333 ierr = DMDestroy(dm);CHKERRQ(ierr); 334 *dm = idm; 335 } 336 break; 337 case 2: 338 { 339 PetscInt numPoints[2] = {4, 1}; 340 PetscInt coneSize[5] = {4, 0, 0, 0, 0}; 341 PetscInt cones[4] = {2, 3, 4, 1}; 342 PetscInt coneOrientations[4] = {0, 0, 0, 0}; 343 PetscScalar vertexCoords[12] = {0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0}; 344 345 depth = 1; 346 ierr = DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 347 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 348 ierr = DMViewFromOptions(idm, NULL, "-in_dm_view");CHKERRQ(ierr); 349 ierr = DMDestroy(dm);CHKERRQ(ierr); 350 *dm = idm; 351 } 352 break; 353 default: SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum); 354 } 355 } else { 356 PetscInt numPoints[4] = {0, 0, 0, 0}; 357 358 ierr = DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr); 359 switch (testNum) { 360 case 1: 361 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 362 ierr = DMViewFromOptions(idm, NULL, "-in_dm_view");CHKERRQ(ierr); 363 ierr = DMDestroy(dm);CHKERRQ(ierr); 364 *dm = idm; 365 break; 366 } 367 } 368 PetscFunctionReturn(0); 369 } 370 371 /* Two tetrahedrons separated by a zero-volume cell with 6 vertices 372 373 cell 6 ___33___10______ cell 374 0 / | \ |\ \ 1 375 21 | 23 | 29 27 376 /12 24 14\ 30 \ \ 377 3-20-|----5--32-|---9--26--7 378 \ 13| 11/ |18 / / 379 19 | 22 | 28 25 380 \ | / |/ / 381 4----31----8------ 382 cell 2 383 */ 384 PetscErrorCode CreateSimplexHybrid_3D(MPI_Comm comm, PetscInt testNum, DM *dm) 385 { 386 DM idm, hdm = NULL; 387 DMLabel faultLabel, hybridLabel; 388 PetscInt p; 389 PetscMPIInt rank; 390 PetscErrorCode ierr; 391 392 PetscFunctionBegin; 393 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 394 if (!rank) { 395 switch (testNum) { 396 case 0: 397 { 398 PetscInt numPoints[2] = {5, 2}; 399 PetscInt coneSize[7] = {4, 4, 0, 0, 0, 0, 0}; 400 PetscInt cones[8] = {4, 3, 5, 2, 5, 3, 4, 6}; 401 PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0}; 402 PetscScalar vertexCoords[15] = {-1.0, 0.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0, 0.0, 0.0}; 403 PetscInt faultPoints[3] = {3, 4, 5}; 404 405 ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 406 for (p = 0; p < 3; ++p) {ierr = DMSetLabelValue(*dm, "fault", faultPoints[p], 1);CHKERRQ(ierr);} 407 } 408 break; 409 case 1: 410 { 411 /* Tets 0,3,5 and 1,2,4 */ 412 PetscInt numPoints[2] = {9, 6}; 413 PetscInt coneSize[15] = {4, 4, 4, 4, 4, 4, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 414 PetscInt cones[24] = { 7, 8, 9, 6, 11, 13, 9, 14, 10, 11, 13, 9, 415 9, 10, 11, 7, 9, 14, 13, 12, 7, 8, 11, 9}; 416 PetscInt coneOrientations[24] = { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 417 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 418 PetscScalar vertexCoords[27] = {-2.0, -1.0, 0.0, -2.0, 0.0, 0.0, -2.0, 0.0, 1.0, 419 0.0, -1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 420 2.0, -1.0, 0.0, 2.0, 0.0, 0.0, 2.0, 0.0, 1.0}; 421 PetscInt faultPoints[3] = {9, 10, 11}; 422 423 ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 424 for (p = 0; p < 3; ++p) {ierr = DMSetLabelValue(*dm, "fault", faultPoints[p], 1);CHKERRQ(ierr);} 425 } 426 break; 427 default: SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum); 428 } 429 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 430 ierr = PetscObjectSetOptionsPrefix((PetscObject) idm, "in_");CHKERRQ(ierr); 431 ierr = DMSetFromOptions(idm);CHKERRQ(ierr); 432 ierr = DMViewFromOptions(idm, NULL, "-dm_view");CHKERRQ(ierr); 433 ierr = DMGetLabel(*dm, "fault", &faultLabel);CHKERRQ(ierr); 434 ierr = DMPlexCreateHybridMesh(idm, faultLabel, NULL, &hybridLabel, NULL, NULL, &hdm);CHKERRQ(ierr); 435 ierr = DMLabelDestroy(&hybridLabel);CHKERRQ(ierr); 436 ierr = DMDestroy(&idm);CHKERRQ(ierr); 437 ierr = DMDestroy(dm);CHKERRQ(ierr); 438 *dm = hdm; 439 } else { 440 PetscInt numPoints[4] = {0, 0, 0, 0}; 441 442 ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr); 443 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 444 ierr = PetscObjectSetOptionsPrefix((PetscObject) idm, "in_");CHKERRQ(ierr); 445 ierr = DMSetFromOptions(idm);CHKERRQ(ierr); 446 ierr = DMViewFromOptions(idm, NULL, "-dm_view");CHKERRQ(ierr); 447 ierr = DMPlexCreateHybridMesh(idm, NULL, NULL, NULL, NULL, NULL, &hdm);CHKERRQ(ierr); 448 ierr = DMDestroy(&idm);CHKERRQ(ierr); 449 ierr = DMDestroy(dm);CHKERRQ(ierr); 450 *dm = hdm; 451 } 452 PetscFunctionReturn(0); 453 } 454 455 PetscErrorCode CreateTensorProduct_3D(MPI_Comm comm, PetscInt testNum, DM *dm) 456 { 457 DM idm; 458 PetscMPIInt rank; 459 PetscErrorCode ierr; 460 461 PetscFunctionBegin; 462 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 463 if (!rank) { 464 switch (testNum) { 465 case 0: 466 { 467 PetscInt numPoints[2] = {12, 2}; 468 PetscInt coneSize[14] = {8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 469 PetscInt cones[16] = {2, 3, 4, 5, 6, 7, 8, 9, 5, 4, 10, 11, 7, 12, 13, 8}; 470 PetscInt coneOrientations[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 471 PetscScalar vertexCoords[36] = {-1.0, -0.5, -0.5, -1.0, 0.5, -0.5, 0.0, 0.5, -0.5, 0.0, -0.5, -0.5, 472 -1.0, -0.5, 0.5, 0.0, -0.5, 0.5, 0.0, 0.5, 0.5, -1.0, 0.5, 0.5, 473 1.0, 0.5, -0.5, 1.0, -0.5, -0.5, 1.0, -0.5, 0.5, 1.0, 0.5, 0.5}; 474 475 ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 476 } 477 break; 478 case 1: 479 { 480 PetscInt numPoints[2] = {8, 1}; 481 PetscInt coneSize[9] = {8, 0, 0, 0, 0, 0, 0, 0, 0}; 482 PetscInt cones[8] = {1, 2, 3, 4, 5, 6, 7, 8}; 483 PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0}; 484 PetscScalar vertexCoords[24] = {-1.0, -1.0, -1.0, -1.0, 1.0, -1.0, 1.0, 1.0, -1.0, 1.0, -1.0, -1.0, 485 -1.0, -1.0, 1.0, 1.0, -1.0, 1.0, 1.0, 1.0, 1.0, -1.0, 1.0, 1.0}; 486 487 ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 488 } 489 break; 490 default: SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum); 491 } 492 } else { 493 PetscInt numPoints[4] = {0, 0, 0, 0}; 494 495 ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr); 496 } 497 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 498 ierr = DMViewFromOptions(idm, NULL, "-in_dm_view");CHKERRQ(ierr); 499 ierr = DMDestroy(dm);CHKERRQ(ierr); 500 *dm = idm; 501 PetscFunctionReturn(0); 502 } 503 504 PetscErrorCode CreateTensorProductHybrid_3D(MPI_Comm comm, PetscInt testNum, DM *dm) 505 { 506 DM idm, hdm = NULL; 507 DMLabel faultLabel; 508 PetscInt p; 509 PetscMPIInt rank; 510 PetscErrorCode ierr; 511 512 PetscFunctionBegin; 513 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 514 if (!rank) { 515 switch (testNum) { 516 case 0: 517 { 518 PetscInt numPoints[2] = {12, 2}; 519 PetscInt coneSize[14] = {8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 520 PetscInt cones[16] = {2, 3, 4, 5, 6, 7, 8, 9, 5, 4, 10, 11, 7, 12, 13, 8}; 521 PetscInt coneOrientations[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 522 PetscScalar vertexCoords[36] = {-1.0, -0.5, -0.5, -1.0, 0.5, -0.5, 0.0, 0.5, -0.5, 0.0, -0.5, -0.5, 523 -1.0, -0.5, 0.5, 0.0, -0.5, 0.5, 0.0, 0.5, 0.5, -1.0, 0.5, 0.5, 524 1.0, 0.5, -0.5, 1.0, -0.5, -0.5, 1.0, -0.5, 0.5, 1.0, 0.5, 0.5}; 525 PetscInt faultPoints[4] = {2, 3, 5, 6}; 526 527 ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 528 for (p = 0; p < 4; ++p) {ierr = DMSetLabelValue(*dm, "fault", faultPoints[p], 1);CHKERRQ(ierr);} 529 } 530 break; 531 case 1: 532 { 533 PetscInt numPoints[2] = {30, 7}; 534 PetscInt coneSize[37] = {8,8,8,8,8,8,8, 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; 535 PetscInt cones[56] = { 8, 21, 20, 7, 13, 12, 23, 24, 536 14, 15, 10, 9, 13, 8, 21, 24, 537 15, 16, 11, 10, 24, 21, 22, 25, 538 30, 29, 28, 21, 35, 24, 33, 34, 539 24, 21, 30, 35, 25, 36, 31, 22, 540 27, 20, 21, 28, 32, 33, 24, 23, 541 15, 24, 13, 14, 19, 18, 17, 26}; 542 PetscInt coneOrientations[56] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; 543 PetscScalar vertexCoords[90] = {-2.0, -2.0, -2.0, -2.0, -1.0, -2.0, -3.0, 0.0, -2.0, -2.0, 1.0, -2.0, -2.0, 2.0, -2.0, -2.0, -2.0, 0.0, 544 -2.0, -1.0, 0.0, -3.0, 0.0, 0.0, -2.0, 1.0, 0.0, -2.0, 2.0, 0.0, -2.0, -1.0, 2.0, -3.0, 0.0, 2.0, 545 -2.0, 1.0, 2.0, 0.0, -2.0, -2.0, 0.0, 0.0, -2.0, 0.0, 2.0, -2.0, 0.0, -2.0, 0.0, 0.0, 0.0, 0.0, 546 0.0, 2.0, 0.0, 0.0, 0.0, 2.0, 2.0, -2.0, -2.0, 2.0, -1.0, -2.0, 3.0, 0.0, -2.0, 2.0, 1.0, -2.0, 547 2.0, 2.0, -2.0, 2.0, -2.0, 0.0, 2.0, -1.0, 0.0, 3.0, 0.0, 0.0, 2.0, 1.0, 0.0, 2.0, 2.0, 0.0}; 548 PetscInt faultPoints[6] = {20, 21, 22, 23, 24, 25}; 549 550 ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 551 for (p = 0; p < 6; ++p) {ierr = DMSetLabelValue(*dm, "fault", faultPoints[p], 1);CHKERRQ(ierr);} 552 } 553 break; 554 default: SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum); 555 } 556 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 557 ierr = PetscObjectSetOptionsPrefix((PetscObject) idm, "in_");CHKERRQ(ierr); 558 ierr = DMSetFromOptions(idm);CHKERRQ(ierr); 559 ierr = DMViewFromOptions(idm, NULL, "-dm_view");CHKERRQ(ierr); 560 ierr = DMGetLabel(*dm, "fault", &faultLabel);CHKERRQ(ierr); 561 ierr = DMPlexCreateHybridMesh(idm, faultLabel, NULL, NULL, NULL, NULL, &hdm);CHKERRQ(ierr); 562 ierr = DMDestroy(&idm);CHKERRQ(ierr); 563 ierr = DMDestroy(dm);CHKERRQ(ierr); 564 *dm = hdm; 565 } else { 566 PetscInt numPoints[4] = {0, 0, 0, 0}; 567 568 ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr); 569 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 570 ierr = PetscObjectSetOptionsPrefix((PetscObject) idm, "in_");CHKERRQ(ierr); 571 ierr = DMSetFromOptions(idm);CHKERRQ(ierr); 572 ierr = DMViewFromOptions(idm, NULL, "-dm_view");CHKERRQ(ierr); 573 ierr = DMPlexCreateHybridMesh(idm, NULL, NULL, NULL, NULL, NULL, &hdm);CHKERRQ(ierr); 574 ierr = DMDestroy(&idm);CHKERRQ(ierr); 575 ierr = DMDestroy(dm);CHKERRQ(ierr); 576 *dm = hdm; 577 } 578 PetscFunctionReturn(0); 579 } 580 581 PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm) 582 { 583 PetscInt dim = user->dim; 584 PetscBool cellHybrid = user->cellHybrid; 585 PetscBool cellSimplex = user->cellSimplex; 586 PetscMPIInt rank, size; 587 PetscErrorCode ierr; 588 589 PetscFunctionBegin; 590 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 591 ierr = MPI_Comm_size(comm, &size);CHKERRQ(ierr); 592 ierr = DMCreate(comm, dm);CHKERRQ(ierr); 593 ierr = DMSetType(*dm, DMPLEX);CHKERRQ(ierr); 594 ierr = DMSetDimension(*dm, dim);CHKERRQ(ierr); 595 switch (dim) { 596 case 1: 597 if (cellHybrid) SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "Cannot make hybrid meshes for dimension %d", dim); 598 ierr = CreateSimplex_1D(comm, dm);CHKERRQ(ierr); 599 break; 600 case 2: 601 if (cellSimplex) { 602 if (cellHybrid) { 603 ierr = CreateSimplexHybrid_2D(comm, user->testNum, dm);CHKERRQ(ierr); 604 } else { 605 ierr = CreateSimplex_2D(comm, dm);CHKERRQ(ierr); 606 } 607 } else { 608 if (cellHybrid) { 609 ierr = CreateTensorProductHybrid_2D(comm, user->testNum, dm);CHKERRQ(ierr); 610 } else { 611 ierr = CreateTensorProduct_2D(comm, user->testNum, dm);CHKERRQ(ierr); 612 } 613 } 614 break; 615 case 3: 616 if (cellSimplex) { 617 if (cellHybrid) { 618 ierr = CreateSimplexHybrid_3D(comm, user->testNum, dm);CHKERRQ(ierr); 619 } else { 620 ierr = CreateSimplex_3D(comm, user->testNum, dm);CHKERRQ(ierr); 621 } 622 } else { 623 if (cellHybrid) { 624 ierr = CreateTensorProductHybrid_3D(comm, user->testNum, dm);CHKERRQ(ierr); 625 } else { 626 ierr = CreateTensorProduct_3D(comm, user->testNum, dm);CHKERRQ(ierr); 627 } 628 } 629 break; 630 default: 631 SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "Cannot make meshes for dimension %d", dim); 632 } 633 if (user->testPartition && size > 1) { 634 PetscPartitioner part; 635 PetscInt *sizes = NULL; 636 PetscInt *points = NULL; 637 638 if (!rank) { 639 if (dim == 2 && cellSimplex && !cellHybrid && size == 2) { 640 switch (user->testNum) { 641 case 0: { 642 PetscInt triSizes_p2[2] = {1, 1}; 643 PetscInt triPoints_p2[2] = {0, 1}; 644 645 ierr = PetscMalloc2(2, &sizes, 2, &points);CHKERRQ(ierr); 646 ierr = PetscArraycpy(sizes, triSizes_p2, 2);CHKERRQ(ierr); 647 ierr = PetscArraycpy(points, triPoints_p2, 2);CHKERRQ(ierr);break;} 648 default: 649 SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for triangular mesh on 2 procs", user->testNum); 650 } 651 } else if (dim == 2 && cellSimplex && cellHybrid && size == 2) { 652 switch (user->testNum) { 653 case 0: { 654 PetscInt triSizes_p2[2] = {1, 2}; 655 PetscInt triPoints_p2[3] = {0, 1, 2}; 656 657 ierr = PetscMalloc2(2, &sizes, 3, &points);CHKERRQ(ierr); 658 ierr = PetscArraycpy(sizes, triSizes_p2, 2);CHKERRQ(ierr); 659 ierr = PetscArraycpy(points, triPoints_p2, 3);CHKERRQ(ierr);break;} 660 default: 661 SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for triangular hybrid mesh on 2 procs", user->testNum); 662 } 663 } else if (dim == 2 && !cellSimplex && !cellHybrid && size == 2) { 664 switch (user->testNum) { 665 case 0: { 666 PetscInt quadSizes_p2[2] = {1, 1}; 667 PetscInt quadPoints_p2[2] = {0, 1}; 668 669 ierr = PetscMalloc2(2, &sizes, 2, &points);CHKERRQ(ierr); 670 ierr = PetscArraycpy(sizes, quadSizes_p2, 2);CHKERRQ(ierr); 671 ierr = PetscArraycpy(points, quadPoints_p2, 2);CHKERRQ(ierr);break;} 672 default: 673 SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for quadrilateral mesh on 2 procs", user->testNum); 674 } 675 } else if (dim == 2 && !cellSimplex && cellHybrid && size == 2) { 676 switch (user->testNum) { 677 case 0: { 678 PetscInt quadSizes_p2[2] = {1, 2}; 679 PetscInt quadPoints_p2[3] = {0, 1, 2}; 680 681 ierr = PetscMalloc2(2, &sizes, 3, &points);CHKERRQ(ierr); 682 ierr = PetscArraycpy(sizes, quadSizes_p2, 2);CHKERRQ(ierr); 683 ierr = PetscArraycpy(points, quadPoints_p2, 3);CHKERRQ(ierr);break;} 684 default: 685 SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for quadrilateral hybrid mesh on 2 procs", user->testNum); 686 } 687 } else if (dim == 3 && cellSimplex && !cellHybrid && size == 2) { 688 switch (user->testNum) { 689 case 0: { 690 PetscInt tetSizes_p2[2] = {1, 1}; 691 PetscInt tetPoints_p2[2] = {0, 1}; 692 693 ierr = PetscMalloc2(2, &sizes, 2, &points);CHKERRQ(ierr); 694 ierr = PetscArraycpy(sizes, tetSizes_p2, 2);CHKERRQ(ierr); 695 ierr = PetscArraycpy(points, tetPoints_p2, 2);CHKERRQ(ierr);break;} 696 case 1: { 697 PetscInt tetSizes_p2[2] = {1, 1}; 698 PetscInt tetPoints_p2[2] = {0, 1}; 699 700 ierr = PetscMalloc2(2, &sizes, 2, &points);CHKERRQ(ierr); 701 ierr = PetscArraycpy(sizes, tetSizes_p2, 2);CHKERRQ(ierr); 702 ierr = PetscArraycpy(points, tetPoints_p2, 2);CHKERRQ(ierr);break;} 703 default: 704 SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for tetrahedral mesh on 2 procs", user->testNum); 705 } 706 } else if (dim == 3 && cellSimplex && cellHybrid && size == 2) { 707 switch (user->testNum) { 708 case 0: { 709 PetscInt tetSizes_p2[2] = {1, 2}; 710 PetscInt tetPoints_p2[3] = {0, 1, 2}; 711 712 ierr = PetscMalloc2(2, &sizes, 3, &points);CHKERRQ(ierr); 713 ierr = PetscArraycpy(sizes, tetSizes_p2, 2);CHKERRQ(ierr); 714 ierr = PetscArraycpy(points, tetPoints_p2, 3);CHKERRQ(ierr);break;} 715 case 1: { 716 PetscInt tetSizes_p2[2] = {3, 4}; 717 PetscInt tetPoints_p2[7] = {0, 3, 5, 1, 2, 4, 6}; 718 719 ierr = PetscMalloc2(2, &sizes, 7, &points);CHKERRQ(ierr); 720 ierr = PetscArraycpy(sizes, tetSizes_p2, 2);CHKERRQ(ierr); 721 ierr = PetscArraycpy(points, tetPoints_p2, 7);CHKERRQ(ierr);break;} 722 default: 723 SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for tetrahedral hybrid mesh on 2 procs", user->testNum); 724 } 725 } else if (dim == 3 && !cellSimplex && !cellHybrid && size == 2) { 726 switch (user->testNum) { 727 case 0: { 728 PetscInt hexSizes_p2[2] = {1, 1}; 729 PetscInt hexPoints_p2[2] = {0, 1}; 730 731 ierr = PetscMalloc2(2, &sizes, 2, &points);CHKERRQ(ierr); 732 ierr = PetscArraycpy(sizes, hexSizes_p2, 2);CHKERRQ(ierr); 733 ierr = PetscArraycpy(points, hexPoints_p2, 2);CHKERRQ(ierr);break;} 734 default: 735 SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for hexahedral mesh on 2 procs", user->testNum); 736 } 737 } else if (dim == 3 && !cellSimplex && cellHybrid && size == 2) { 738 switch (user->testNum) { 739 case 0: { 740 PetscInt hexSizes_p2[2] = {1, 1}; 741 PetscInt hexPoints_p2[2] = {0, 1}; 742 743 ierr = PetscMalloc2(2, &sizes, 2, &points);CHKERRQ(ierr); 744 ierr = PetscArraycpy(sizes, hexSizes_p2, 2);CHKERRQ(ierr); 745 ierr = PetscArraycpy(points, hexPoints_p2, 2);CHKERRQ(ierr);break;} 746 case 1: { 747 PetscInt hexSizes_p2[2] = {5, 4}; 748 PetscInt hexPoints_p2[9] = {3, 4, 5, 7, 8, 0, 1, 2, 6}; 749 750 ierr = PetscMalloc2(2, &sizes, 9, &points);CHKERRQ(ierr); 751 ierr = PetscArraycpy(sizes, hexSizes_p2, 2);CHKERRQ(ierr); 752 ierr = PetscArraycpy(points, hexPoints_p2, 9);CHKERRQ(ierr);break;} 753 default: 754 SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for hexahedral hybrid mesh on 2 procs", user->testNum); 755 } 756 } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test partition"); 757 } 758 ierr = DMPlexGetPartitioner(*dm, &part);CHKERRQ(ierr); 759 ierr = PetscPartitionerSetType(part, PETSCPARTITIONERSHELL);CHKERRQ(ierr); 760 ierr = PetscPartitionerShellSetPartition(part, size, sizes, points);CHKERRQ(ierr); 761 ierr = PetscFree2(sizes, points);CHKERRQ(ierr); 762 } else { 763 PetscPartitioner part; 764 765 ierr = DMPlexGetPartitioner(*dm,&part);CHKERRQ(ierr); 766 ierr = PetscPartitionerSetFromOptions(part);CHKERRQ(ierr); 767 } 768 { 769 DM pdm = NULL; 770 771 ierr = DMPlexDistribute(*dm, 0, NULL, &pdm);CHKERRQ(ierr); 772 if (pdm) { 773 ierr = DMViewFromOptions(pdm, NULL, "-dm_view");CHKERRQ(ierr); 774 ierr = DMDestroy(dm);CHKERRQ(ierr); 775 *dm = pdm; 776 } 777 } 778 ierr = DMSetFromOptions(*dm);CHKERRQ(ierr); 779 if (user->uninterpolate || user->reinterpolate) { 780 DM udm = NULL; 781 782 ierr = DMPlexUninterpolate(*dm, &udm);CHKERRQ(ierr); 783 ierr = DMPlexCopyCoordinates(*dm, udm);CHKERRQ(ierr); 784 ierr = DMDestroy(dm);CHKERRQ(ierr); 785 *dm = udm; 786 } 787 if (user->reinterpolate) { 788 DM idm = NULL; 789 790 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 791 ierr = DMPlexCopyCoordinates(*dm, idm);CHKERRQ(ierr); 792 ierr = DMDestroy(dm);CHKERRQ(ierr); 793 *dm = idm; 794 } 795 ierr = PetscObjectSetName((PetscObject) *dm, "Hybrid Mesh");CHKERRQ(ierr); 796 ierr = DMViewFromOptions(*dm, NULL, "-dm_view");CHKERRQ(ierr); 797 ierr = PetscObjectSetOptionsPrefix((PetscObject) *dm, "hyb_");CHKERRQ(ierr); 798 ierr = DMSetFromOptions(*dm);CHKERRQ(ierr); 799 PetscFunctionReturn(0); 800 } 801 802 int main(int argc, char **argv) 803 { 804 DM dm; 805 AppCtx user; /* user-defined work context */ 806 PetscErrorCode ierr; 807 808 ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr; 809 ierr = ProcessOptions(PETSC_COMM_WORLD, &user);CHKERRQ(ierr); 810 ierr = CreateMesh(PETSC_COMM_WORLD, &user, &dm);CHKERRQ(ierr); 811 ierr = DMDestroy(&dm);CHKERRQ(ierr); 812 ierr = PetscFinalize(); 813 return ierr; 814 } 815 816 /*TEST 817 818 # 1D Simplex 29-31 819 testset: 820 args: -dim 1 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all 821 test: 822 suffix: 29 823 test: 824 suffix: 30 825 args: -dm_refine 1 826 test: 827 suffix: 31 828 args: -dm_refine 5 829 830 # 2D Simplex 0-3 831 testset: 832 args: -dim 2 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all 833 test: 834 suffix: 0 835 args: -hyb_dm_plex_check_faces 836 test: 837 suffix: 1 838 args: -dm_refine 1 -hyb_dm_plex_check_faces 839 test: 840 suffix: 2 841 nsize: 2 842 args: -hyb_dm_plex_check_faces 843 test: 844 suffix: 3 845 nsize: 2 846 args: -dm_refine 1 -hyb_dm_plex_check_faces 847 test: 848 suffix: 32 849 args: -dm_refine 1 -uninterpolate 850 test: 851 suffix: 33 852 nsize: 2 853 args: -dm_refine 1 -uninterpolate 854 test: 855 suffix: 34 856 nsize: 2 857 args: -dm_refine 3 -uninterpolate 858 859 # 2D Hybrid Simplex 4-7 860 testset: 861 args: -dim 2 -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all 862 test: 863 suffix: 4 864 test: 865 suffix: 5 866 args: -dm_refine 1 867 test: 868 suffix: 6 869 nsize: 2 870 test: 871 suffix: 7 872 nsize: 2 873 args: -dm_refine 1 874 test: 875 suffix: 24 876 args: -test_num 1 -dm_refine 1 877 878 # 2D Quad 12-13 879 testset: 880 args: -dim 2 -cell_simplex 0 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all 881 test: 882 suffix: 12 883 args: -dm_refine 1 884 test: 885 suffix: 13 886 nsize: 2 887 args: -dm_refine 1 888 889 # 2D Hybrid Quad 27-28 890 testset: 891 args: -dim 2 -cell_simplex 0 -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all 892 test: 893 suffix: 27 894 args: -dm_refine 1 895 test: 896 suffix: 28 897 nsize: 2 898 args: -dm_refine 1 899 900 # 3D Simplex 8-11 901 testset: 902 args: -dim 3 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all 903 test: 904 suffix: 8 905 args: -dm_refine 1 906 test: 907 suffix: 9 908 nsize: 2 909 args: -dm_refine 1 910 test: 911 suffix: 10 912 args: -test_num 1 -dm_refine 1 913 test: 914 suffix: 11 915 nsize: 2 916 args: -test_num 1 -dm_refine 1 917 test: 918 suffix: 25 919 args: -test_num 2 -dm_refine 2 920 921 # 3D Hybrid Simplex 16-19 922 testset: 923 args: -dim 3 -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all 924 test: 925 suffix: 16 926 args: -dm_refine 1 927 test: 928 suffix: 17 929 nsize: 2 930 args: -dm_refine 1 931 test: 932 suffix: 18 933 args: -test_num 1 -dm_refine 1 934 test: 935 suffix: 19 936 nsize: 2 937 args: -test_num 1 -dm_refine 1 938 939 # 3D Hex 14-15 940 testset: 941 args: -dim 3 -cell_simplex 0 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all 942 test: 943 suffix: 14 944 args: -dm_refine 1 945 test: 946 suffix: 15 947 nsize: 2 948 args: -dm_refine 1 949 test: 950 suffix: 26 951 args: -test_num 1 -dm_refine 2 952 953 # 3D Hybrid Hex 20-23 954 testset: 955 args: -dim 3 -cell_simplex 0 -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all 956 test: 957 suffix: 20 958 args: -dm_refine 1 959 test: 960 suffix: 21 961 nsize: 2 962 args: -dm_refine 1 963 test: 964 suffix: 22 965 args: -test_num 1 -dm_refine 1 966 test: 967 suffix: 23 968 nsize: 2 969 args: -test_num 1 -dm_refine 1 970 971 # Hybrid interpolation 972 # TODO Setup new tests (like -reinterpolate) that interpolate hybrid cells 973 testset: 974 nsize: 2 975 args: -test_partition 0 -petscpartitioner_type simple -dm_view -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all 976 test: 977 suffix: hybint_2d_0 978 args: -dim 2 -dm_refine 2 979 test: 980 suffix: hybint_2d_1 981 args: -dim 2 -dm_refine 2 -test_num 1 982 test: 983 suffix: hybint_3d_0 984 args: -dim 3 -dm_refine 1 985 test: 986 suffix: hybint_3d_1 987 args: -dim 3 -dm_refine 1 -test_num 1 988 989 TEST*/ 990