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 CHKERRQ(PetscOptionsBoundedInt("-debug", "The debugging level", "ex4.c", options->debug, &options->debug, NULL,0)); 32 CHKERRQ(PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex4.c", options->dim, &options->dim, NULL,1,3)); 33 CHKERRQ(PetscOptionsBool("-cell_hybrid", "Use a hybrid mesh", "ex4.c", options->cellHybrid, &options->cellHybrid, NULL)); 34 CHKERRQ(PetscOptionsBool("-cell_simplex", "Use simplices if true, otherwise hexes", "ex4.c", options->cellSimplex, &options->cellSimplex, NULL)); 35 CHKERRQ(PetscOptionsBool("-test_partition", "Use a fixed partition for testing", "ex4.c", options->testPartition, &options->testPartition, NULL)); 36 CHKERRQ(PetscOptionsBoundedInt("-test_num", "The particular mesh to test", "ex4.c", options->testNum, &options->testNum, NULL,0)); 37 CHKERRQ(PetscOptionsBool("-uninterpolate", "Uninterpolate the mesh at the end", "ex4.c", options->uninterpolate, &options->uninterpolate, NULL)); 38 CHKERRQ(PetscOptionsBool("-reinterpolate", "Reinterpolate the mesh at the end", "ex4.c", options->reinterpolate, &options->reinterpolate, NULL)); 39 ierr = PetscOptionsEnd();CHKERRQ(ierr); 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 57 PetscFunctionBegin; 58 CHKERRMPI(MPI_Comm_rank(comm, &rank)); 59 if (rank == 0) { 60 PetscInt numPoints[2] = {3, 2}; 61 PetscInt coneSize[5] = {2, 2, 0, 0, 0}; 62 PetscInt cones[4] = {2, 3, 3, 4}; 63 PetscInt coneOrientations[16] = {0, 0, 0, 0}; 64 PetscScalar vertexCoords[3] = {-1.0, 0.0, 1.0}; 65 66 CHKERRQ(DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 67 } else { 68 PetscInt numPoints[2] = {0, 0}; 69 70 CHKERRQ(DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL)); 71 } 72 PetscFunctionReturn(0); 73 } 74 75 /* Two triangles 76 4 77 / | \ 78 8 | 10 79 / | \ 80 2 0 7 1 5 81 \ | / 82 6 | 9 83 \ | / 84 3 85 86 Becomes 87 10 88 / | \ 89 21 | 26 90 / | \ 91 14 2 20 4 16 92 /|\ | /|\ 93 22 | 28 | 32 | 25 94 / | \ | / | 6\ 95 8 29 3 13 7 31 11 96 \0 | / | \ | / 97 17 | 27 | 30 | 24 98 \|/ | \|/ 99 12 1 19 5 15 100 \ | / 101 18 | 23 102 \ | / 103 9 104 */ 105 PetscErrorCode CreateSimplex_2D(MPI_Comm comm, DM *dm) 106 { 107 PetscInt depth = 2; 108 PetscMPIInt rank; 109 110 PetscFunctionBegin; 111 CHKERRMPI(MPI_Comm_rank(comm, &rank)); 112 if (rank == 0) { 113 PetscInt numPoints[3] = {4, 5, 2}; 114 PetscInt coneSize[11] = {3, 3, 0, 0, 0, 0, 2, 2, 2, 2, 2}; 115 PetscInt cones[16] = {6, 7, 8, 7, 9, 10, 2, 3, 3, 4, 4, 2, 3, 5, 5, 4}; 116 PetscInt coneOrientations[16] = {0, 0, 0, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 117 PetscScalar vertexCoords[8] = {-0.5, 0.0, 0.0, -0.5, 0.0, 0.5, 0.5, 0.0}; 118 119 CHKERRQ(DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 120 } else { 121 PetscInt numPoints[3] = {0, 0, 0}; 122 123 CHKERRQ(DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL)); 124 } 125 PetscFunctionReturn(0); 126 } 127 128 /* Two triangles separated by a zero-volume cell with 4 vertices/2 edges 129 5--16--8 130 / | | \ 131 11 | | 12 132 / | | \ 133 3 0 10 2 14 1 6 134 \ | | / 135 9 | | 13 136 \ | | / 137 4--15--7 138 */ 139 PetscErrorCode CreateSimplexHybrid_2D(MPI_Comm comm, PetscInt testNum, DM *dm) 140 { 141 DM idm, hdm = NULL; 142 DMLabel faultLabel, hybridLabel; 143 PetscInt p; 144 PetscMPIInt rank; 145 146 PetscFunctionBegin; 147 CHKERRMPI(MPI_Comm_rank(comm, &rank)); 148 if (rank == 0) { 149 switch (testNum) { 150 case 0: 151 { 152 PetscInt numPoints[2] = {4, 2}; 153 PetscInt coneSize[6] = {3, 3, 0, 0, 0, 0}; 154 PetscInt cones[6] = {2, 3, 4, 5, 4, 3}; 155 PetscInt coneOrientations[6] = {0, 0, 0, 0, 0, 0}; 156 PetscScalar vertexCoords[8] = {-1.0, -0.5, 0.0, -0.5, 0.0, 0.5, 1.0, 0.5}; 157 PetscInt faultPoints[2] = {3, 4}; 158 159 CHKERRQ(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 160 for (p = 0; p < 2; ++p) CHKERRQ(DMSetLabelValue(*dm, "fault", faultPoints[p], 1)); 161 } 162 break; 163 case 1: 164 { 165 PetscInt numPoints[2] = {5, 4}; 166 PetscInt coneSize[9] = {3, 3, 3, 3, 0, 0, 0, 0, 0}; 167 PetscInt cones[12] = {4, 5, 6, 6, 7, 4, 6, 5, 8, 6, 8, 7}; 168 PetscInt coneOrientations[12] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 169 PetscScalar vertexCoords[10] = {-1.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0, 1.0, 1.0, 0.0}; 170 PetscInt faultPoints[3] = {5, 6, 7}; 171 172 CHKERRQ(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 173 for (p = 0; p < 3; ++p) CHKERRQ(DMSetLabelValue(*dm, "fault", faultPoints[p], 1)); 174 } 175 break; 176 default: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum); 177 } 178 CHKERRQ(DMPlexInterpolate(*dm, &idm)); 179 CHKERRQ(PetscObjectSetOptionsPrefix((PetscObject) idm, "in_")); 180 CHKERRQ(DMPlexDistributeSetDefault(idm, PETSC_FALSE)); 181 CHKERRQ(DMSetFromOptions(idm)); 182 CHKERRQ(DMViewFromOptions(idm, NULL, "-dm_view")); 183 CHKERRQ(DMGetLabel(*dm, "fault", &faultLabel)); 184 CHKERRQ(DMPlexCreateHybridMesh(idm, faultLabel, NULL, &hybridLabel, NULL, NULL, &hdm)); 185 CHKERRQ(DMLabelDestroy(&hybridLabel)); 186 CHKERRQ(DMDestroy(&idm)); 187 CHKERRQ(DMDestroy(dm)); 188 *dm = hdm; 189 } else { 190 PetscInt numPoints[2] = {0, 0}; 191 192 CHKERRQ(DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL)); 193 CHKERRQ(DMPlexInterpolate(*dm, &idm)); 194 CHKERRQ(PetscObjectSetOptionsPrefix((PetscObject) idm, "in_")); 195 CHKERRQ(DMPlexDistributeSetDefault(idm, PETSC_FALSE)); 196 CHKERRQ(DMSetFromOptions(idm)); 197 CHKERRQ(DMViewFromOptions(idm, NULL, "-dm_view")); 198 CHKERRQ(DMPlexCreateHybridMesh(idm, NULL, NULL, NULL, NULL, NULL, &hdm)); 199 CHKERRQ(DMDestroy(&idm)); 200 CHKERRQ(DMDestroy(dm)); 201 *dm = hdm; 202 } 203 PetscFunctionReturn(0); 204 } 205 206 /* Two quadrilaterals 207 208 5----10-----4----14-----7 209 | | | 210 | | | 211 | | | 212 11 0 9 1 13 213 | | | 214 | | | 215 | | | 216 2-----8-----3----12-----6 217 */ 218 PetscErrorCode CreateTensorProduct_2D(MPI_Comm comm, PetscInt testNum, DM *dm) 219 { 220 PetscInt depth = 2; 221 PetscMPIInt rank; 222 223 PetscFunctionBegin; 224 CHKERRMPI(MPI_Comm_rank(comm, &rank)); 225 if (rank == 0) { 226 PetscInt numPoints[3] = {6, 7, 2}; 227 PetscInt coneSize[15] = {4, 4, 0, 0, 0, 0, 0, 0, 2, 2, 2, 2, 2, 2, 2}; 228 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}; 229 PetscInt coneOrientations[22] = {0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 230 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}; 231 232 CHKERRQ(DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 233 } else { 234 PetscInt numPoints[3] = {0, 0, 0}; 235 236 CHKERRQ(DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL)); 237 } 238 PetscFunctionReturn(0); 239 } 240 241 PetscErrorCode CreateTensorProductHybrid_2D(MPI_Comm comm, PetscInt testNum, DM *dm) 242 { 243 DM idm, hdm = NULL; 244 DMLabel faultLabel, hybridLabel; 245 PetscInt p; 246 PetscMPIInt rank; 247 248 PetscFunctionBegin; 249 CHKERRMPI(MPI_Comm_rank(comm, &rank)); 250 if (rank == 0) { 251 PetscInt numPoints[2] = {6, 2}; 252 PetscInt coneSize[8] = {4, 4, 0, 0, 0, 0, 0, 0}; 253 PetscInt cones[8] = {2, 3, 4, 5, 3, 6, 7, 4,}; 254 PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0}; 255 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}; 256 PetscInt faultPoints[2] = {3, 4}; 257 258 CHKERRQ(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 259 for (p = 0; p < 2; ++p) CHKERRQ(DMSetLabelValue(*dm, "fault", faultPoints[p], 1)); 260 CHKERRQ(DMPlexInterpolate(*dm, &idm)); 261 CHKERRQ(PetscObjectSetOptionsPrefix((PetscObject) idm, "in_")); 262 CHKERRQ(DMPlexDistributeSetDefault(idm, PETSC_FALSE)); 263 CHKERRQ(DMSetFromOptions(idm)); 264 CHKERRQ(DMViewFromOptions(idm, NULL, "-dm_view")); 265 CHKERRQ(DMGetLabel(*dm, "fault", &faultLabel)); 266 CHKERRQ(DMPlexCreateHybridMesh(idm, faultLabel, NULL, &hybridLabel, NULL, NULL, &hdm)); 267 CHKERRQ(DMLabelDestroy(&hybridLabel)); 268 } else { 269 PetscInt numPoints[3] = {0, 0, 0}; 270 271 CHKERRQ(DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL)); 272 CHKERRQ(DMPlexInterpolate(*dm, &idm)); 273 CHKERRQ(PetscObjectSetOptionsPrefix((PetscObject) idm, "in_")); 274 CHKERRQ(DMPlexDistributeSetDefault(idm, PETSC_FALSE)); 275 CHKERRQ(DMSetFromOptions(idm)); 276 CHKERRQ(DMViewFromOptions(idm, NULL, "-dm_view")); 277 CHKERRQ(DMPlexCreateHybridMesh(idm, NULL, NULL, NULL, NULL, NULL, &hdm)); 278 } 279 CHKERRQ(DMDestroy(&idm)); 280 CHKERRQ(DMDestroy(dm)); 281 *dm = hdm; 282 PetscFunctionReturn(0); 283 } 284 285 /* Two tetrahedrons 286 287 cell 5 5______ cell 288 0 / | \ |\ \ 1 289 17 | 18 | 18 13 21 290 /8 19 10\ 19 \ \ 291 2-14-|----4 | 4--22--6 292 \ 9 | 7 / |10 / / 293 16 | 15 | 15 12 20 294 \ | / |/ / 295 3 3------ 296 */ 297 PetscErrorCode CreateSimplex_3D(MPI_Comm comm, PetscInt testNum, DM *dm) 298 { 299 DM idm; 300 PetscInt depth = 3; 301 PetscMPIInt rank; 302 303 PetscFunctionBegin; 304 CHKERRMPI(MPI_Comm_rank(comm, &rank)); 305 if (rank == 0) { 306 switch (testNum) { 307 case 0: 308 { 309 PetscInt numPoints[4] = {5, 9, 7, 2}; 310 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}; 311 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}; 312 PetscInt coneOrientations[47] = { 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, -1, -1, 0, -1, -1, -1, -1, 0, 0, -1, -1, 0, -1, -1, -1, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 313 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}; 314 315 CHKERRQ(DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 316 } 317 break; 318 case 1: 319 { 320 PetscInt numPoints[2] = {5, 2}; 321 PetscInt coneSize[7] = {4, 4, 0, 0, 0, 0, 0}; 322 PetscInt cones[8] = {4, 3, 5, 2, 5, 3, 4, 6}; 323 PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0}; 324 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}; 325 326 depth = 1; 327 CHKERRQ(DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 328 CHKERRQ(DMPlexInterpolate(*dm, &idm)); 329 CHKERRQ(DMViewFromOptions(idm, NULL, "-in_dm_view")); 330 CHKERRQ(DMDestroy(dm)); 331 *dm = idm; 332 } 333 break; 334 case 2: 335 { 336 PetscInt numPoints[2] = {4, 1}; 337 PetscInt coneSize[5] = {4, 0, 0, 0, 0}; 338 PetscInt cones[4] = {2, 3, 4, 1}; 339 PetscInt coneOrientations[4] = {0, 0, 0, 0}; 340 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}; 341 342 depth = 1; 343 CHKERRQ(DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 344 CHKERRQ(DMPlexInterpolate(*dm, &idm)); 345 CHKERRQ(DMViewFromOptions(idm, NULL, "-in_dm_view")); 346 CHKERRQ(DMDestroy(dm)); 347 *dm = idm; 348 } 349 break; 350 default: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum); 351 } 352 } else { 353 PetscInt numPoints[4] = {0, 0, 0, 0}; 354 355 CHKERRQ(DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL)); 356 switch (testNum) { 357 case 1: 358 CHKERRQ(DMPlexInterpolate(*dm, &idm)); 359 CHKERRQ(DMViewFromOptions(idm, NULL, "-in_dm_view")); 360 CHKERRQ(DMDestroy(dm)); 361 *dm = idm; 362 break; 363 } 364 } 365 PetscFunctionReturn(0); 366 } 367 368 /* Two tetrahedrons separated by a zero-volume cell with 6 vertices 369 370 cell 6 ___33___10______ cell 371 0 / | \ |\ \ 1 372 21 | 23 | 29 27 373 /12 24 14\ 30 \ \ 374 3-20-|----5--32-|---9--26--7 375 \ 13| 11/ |18 / / 376 19 | 22 | 28 25 377 \ | / |/ / 378 4----31----8------ 379 cell 2 380 */ 381 PetscErrorCode CreateSimplexHybrid_3D(MPI_Comm comm, PetscInt testNum, DM *dm) 382 { 383 DM idm, hdm = NULL; 384 DMLabel faultLabel, hybridLabel; 385 PetscInt p; 386 PetscMPIInt rank; 387 388 PetscFunctionBegin; 389 CHKERRMPI(MPI_Comm_rank(comm, &rank)); 390 if (rank == 0) { 391 switch (testNum) { 392 case 0: 393 { 394 PetscInt numPoints[2] = {5, 2}; 395 PetscInt coneSize[7] = {4, 4, 0, 0, 0, 0, 0}; 396 PetscInt cones[8] = {4, 3, 5, 2, 5, 3, 4, 6}; 397 PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0}; 398 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}; 399 PetscInt faultPoints[3] = {3, 4, 5}; 400 401 CHKERRQ(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 402 for (p = 0; p < 3; ++p) CHKERRQ(DMSetLabelValue(*dm, "fault", faultPoints[p], 1)); 403 } 404 break; 405 case 1: 406 { 407 /* Tets 0,3,5 and 1,2,4 */ 408 PetscInt numPoints[2] = {9, 6}; 409 PetscInt coneSize[15] = {4, 4, 4, 4, 4, 4, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 410 PetscInt cones[24] = { 7, 8, 9, 6, 11, 13, 9, 14, 10, 11, 13, 9, 411 9, 10, 11, 7, 9, 14, 13, 12, 7, 8, 11, 9}; 412 PetscInt coneOrientations[24] = { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 413 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 414 PetscScalar vertexCoords[27] = {-2.0, -1.0, 0.0, -2.0, 0.0, 0.0, -2.0, 0.0, 1.0, 415 0.0, -1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 416 2.0, -1.0, 0.0, 2.0, 0.0, 0.0, 2.0, 0.0, 1.0}; 417 PetscInt faultPoints[3] = {9, 10, 11}; 418 419 CHKERRQ(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 420 for (p = 0; p < 3; ++p) CHKERRQ(DMSetLabelValue(*dm, "fault", faultPoints[p], 1)); 421 } 422 break; 423 default: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum); 424 } 425 CHKERRQ(DMPlexInterpolate(*dm, &idm)); 426 CHKERRQ(PetscObjectSetOptionsPrefix((PetscObject) idm, "in_")); 427 CHKERRQ(DMPlexDistributeSetDefault(idm, PETSC_FALSE)); 428 CHKERRQ(DMSetFromOptions(idm)); 429 CHKERRQ(DMViewFromOptions(idm, NULL, "-dm_view")); 430 CHKERRQ(DMGetLabel(*dm, "fault", &faultLabel)); 431 CHKERRQ(DMPlexCreateHybridMesh(idm, faultLabel, NULL, &hybridLabel, NULL, NULL, &hdm)); 432 CHKERRQ(DMLabelDestroy(&hybridLabel)); 433 CHKERRQ(DMDestroy(&idm)); 434 CHKERRQ(DMDestroy(dm)); 435 *dm = hdm; 436 } else { 437 PetscInt numPoints[4] = {0, 0, 0, 0}; 438 439 CHKERRQ(DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL)); 440 CHKERRQ(DMPlexInterpolate(*dm, &idm)); 441 CHKERRQ(PetscObjectSetOptionsPrefix((PetscObject) idm, "in_")); 442 CHKERRQ(DMPlexDistributeSetDefault(idm, PETSC_FALSE)); 443 CHKERRQ(DMSetFromOptions(idm)); 444 CHKERRQ(DMViewFromOptions(idm, NULL, "-dm_view")); 445 CHKERRQ(DMPlexCreateHybridMesh(idm, NULL, NULL, NULL, NULL, NULL, &hdm)); 446 CHKERRQ(DMDestroy(&idm)); 447 CHKERRQ(DMDestroy(dm)); 448 *dm = hdm; 449 } 450 PetscFunctionReturn(0); 451 } 452 453 PetscErrorCode CreateTensorProduct_3D(MPI_Comm comm, PetscInt testNum, DM *dm) 454 { 455 DM idm; 456 PetscMPIInt rank; 457 458 PetscFunctionBegin; 459 CHKERRMPI(MPI_Comm_rank(comm, &rank)); 460 if (rank == 0) { 461 switch (testNum) { 462 case 0: 463 { 464 PetscInt numPoints[2] = {12, 2}; 465 PetscInt coneSize[14] = {8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 466 PetscInt cones[16] = {2, 3, 4, 5, 6, 7, 8, 9, 5, 4, 10, 11, 7, 12, 13, 8}; 467 PetscInt coneOrientations[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 468 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, 469 -1.0, -0.5, 0.5, 0.0, -0.5, 0.5, 0.0, 0.5, 0.5, -1.0, 0.5, 0.5, 470 1.0, 0.5, -0.5, 1.0, -0.5, -0.5, 1.0, -0.5, 0.5, 1.0, 0.5, 0.5}; 471 472 CHKERRQ(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 473 } 474 break; 475 case 1: 476 { 477 PetscInt numPoints[2] = {8, 1}; 478 PetscInt coneSize[9] = {8, 0, 0, 0, 0, 0, 0, 0, 0}; 479 PetscInt cones[8] = {1, 2, 3, 4, 5, 6, 7, 8}; 480 PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0}; 481 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, 482 -1.0, -1.0, 1.0, 1.0, -1.0, 1.0, 1.0, 1.0, 1.0, -1.0, 1.0, 1.0}; 483 484 CHKERRQ(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 485 } 486 break; 487 default: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum); 488 } 489 } else { 490 PetscInt numPoints[4] = {0, 0, 0, 0}; 491 492 CHKERRQ(DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL)); 493 } 494 CHKERRQ(DMPlexInterpolate(*dm, &idm)); 495 CHKERRQ(DMViewFromOptions(idm, NULL, "-in_dm_view")); 496 CHKERRQ(DMDestroy(dm)); 497 *dm = idm; 498 PetscFunctionReturn(0); 499 } 500 501 PetscErrorCode CreateTensorProductHybrid_3D(MPI_Comm comm, PetscInt testNum, DM *dm) 502 { 503 DM idm, hdm = NULL; 504 DMLabel faultLabel; 505 PetscInt p; 506 PetscMPIInt rank; 507 508 PetscFunctionBegin; 509 CHKERRMPI(MPI_Comm_rank(comm, &rank)); 510 if (rank == 0) { 511 switch (testNum) { 512 case 0: 513 { 514 PetscInt numPoints[2] = {12, 2}; 515 PetscInt coneSize[14] = {8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 516 PetscInt cones[16] = {2, 3, 4, 5, 6, 7, 8, 9, 5, 4, 10, 11, 7, 12, 13, 8}; 517 PetscInt coneOrientations[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 518 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, 519 -1.0, -0.5, 0.5, 0.0, -0.5, 0.5, 0.0, 0.5, 0.5, -1.0, 0.5, 0.5, 520 1.0, 0.5, -0.5, 1.0, -0.5, -0.5, 1.0, -0.5, 0.5, 1.0, 0.5, 0.5}; 521 PetscInt faultPoints[4] = {2, 3, 5, 6}; 522 523 CHKERRQ(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 524 for (p = 0; p < 4; ++p) CHKERRQ(DMSetLabelValue(*dm, "fault", faultPoints[p], 1)); 525 } 526 break; 527 case 1: 528 { 529 PetscInt numPoints[2] = {30, 7}; 530 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}; 531 PetscInt cones[56] = { 8, 21, 20, 7, 13, 12, 23, 24, 532 14, 15, 10, 9, 13, 8, 21, 24, 533 15, 16, 11, 10, 24, 21, 22, 25, 534 30, 29, 28, 21, 35, 24, 33, 34, 535 24, 21, 30, 35, 25, 36, 31, 22, 536 27, 20, 21, 28, 32, 33, 24, 23, 537 15, 24, 13, 14, 19, 18, 17, 26}; 538 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}; 539 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, 540 -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, 541 -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, 542 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, 543 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}; 544 PetscInt faultPoints[6] = {20, 21, 22, 23, 24, 25}; 545 546 CHKERRQ(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 547 for (p = 0; p < 6; ++p) CHKERRQ(DMSetLabelValue(*dm, "fault", faultPoints[p], 1)); 548 } 549 break; 550 default: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum); 551 } 552 CHKERRQ(DMPlexInterpolate(*dm, &idm)); 553 CHKERRQ(PetscObjectSetOptionsPrefix((PetscObject) idm, "in_")); 554 CHKERRQ(DMPlexDistributeSetDefault(idm, PETSC_FALSE)); 555 CHKERRQ(DMSetFromOptions(idm)); 556 CHKERRQ(DMViewFromOptions(idm, NULL, "-dm_view")); 557 CHKERRQ(DMGetLabel(*dm, "fault", &faultLabel)); 558 CHKERRQ(DMPlexCreateHybridMesh(idm, faultLabel, NULL, NULL, NULL, NULL, &hdm)); 559 CHKERRQ(DMDestroy(&idm)); 560 CHKERRQ(DMDestroy(dm)); 561 *dm = hdm; 562 } else { 563 PetscInt numPoints[4] = {0, 0, 0, 0}; 564 565 CHKERRQ(DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL)); 566 CHKERRQ(DMPlexInterpolate(*dm, &idm)); 567 CHKERRQ(PetscObjectSetOptionsPrefix((PetscObject) idm, "in_")); 568 CHKERRQ(DMPlexDistributeSetDefault(idm, PETSC_FALSE)); 569 CHKERRQ(DMSetFromOptions(idm)); 570 CHKERRQ(DMViewFromOptions(idm, NULL, "-dm_view")); 571 CHKERRQ(DMPlexCreateHybridMesh(idm, NULL, NULL, NULL, NULL, NULL, &hdm)); 572 CHKERRQ(DMDestroy(&idm)); 573 CHKERRQ(DMDestroy(dm)); 574 *dm = hdm; 575 } 576 PetscFunctionReturn(0); 577 } 578 579 PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm) 580 { 581 PetscInt dim = user->dim; 582 PetscBool cellHybrid = user->cellHybrid; 583 PetscBool cellSimplex = user->cellSimplex; 584 PetscMPIInt rank, size; 585 586 PetscFunctionBegin; 587 CHKERRMPI(MPI_Comm_rank(comm, &rank)); 588 CHKERRMPI(MPI_Comm_size(comm, &size)); 589 CHKERRQ(DMCreate(comm, dm)); 590 CHKERRQ(DMSetType(*dm, DMPLEX)); 591 CHKERRQ(DMSetDimension(*dm, dim)); 592 switch (dim) { 593 case 1: 594 PetscCheck(!cellHybrid,comm, PETSC_ERR_ARG_OUTOFRANGE, "Cannot make hybrid meshes for dimension %d", dim); 595 CHKERRQ(CreateSimplex_1D(comm, dm)); 596 break; 597 case 2: 598 if (cellSimplex) { 599 if (cellHybrid) { 600 CHKERRQ(CreateSimplexHybrid_2D(comm, user->testNum, dm)); 601 } else { 602 CHKERRQ(CreateSimplex_2D(comm, dm)); 603 } 604 } else { 605 if (cellHybrid) { 606 CHKERRQ(CreateTensorProductHybrid_2D(comm, user->testNum, dm)); 607 } else { 608 CHKERRQ(CreateTensorProduct_2D(comm, user->testNum, dm)); 609 } 610 } 611 break; 612 case 3: 613 if (cellSimplex) { 614 if (cellHybrid) { 615 CHKERRQ(CreateSimplexHybrid_3D(comm, user->testNum, dm)); 616 } else { 617 CHKERRQ(CreateSimplex_3D(comm, user->testNum, dm)); 618 } 619 } else { 620 if (cellHybrid) { 621 CHKERRQ(CreateTensorProductHybrid_3D(comm, user->testNum, dm)); 622 } else { 623 CHKERRQ(CreateTensorProduct_3D(comm, user->testNum, dm)); 624 } 625 } 626 break; 627 default: 628 SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "Cannot make meshes for dimension %d", dim); 629 } 630 if (user->testPartition && size > 1) { 631 PetscPartitioner part; 632 PetscInt *sizes = NULL; 633 PetscInt *points = NULL; 634 635 if (rank == 0) { 636 if (dim == 2 && cellSimplex && !cellHybrid && size == 2) { 637 switch (user->testNum) { 638 case 0: { 639 PetscInt triSizes_p2[2] = {1, 1}; 640 PetscInt triPoints_p2[2] = {0, 1}; 641 642 CHKERRQ(PetscMalloc2(2, &sizes, 2, &points)); 643 CHKERRQ(PetscArraycpy(sizes, triSizes_p2, 2)); 644 CHKERRQ(PetscArraycpy(points, triPoints_p2, 2));break;} 645 default: 646 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for triangular mesh on 2 procs", user->testNum); 647 } 648 } else if (dim == 2 && cellSimplex && cellHybrid && size == 2) { 649 switch (user->testNum) { 650 case 0: { 651 PetscInt triSizes_p2[2] = {1, 2}; 652 PetscInt triPoints_p2[3] = {0, 1, 2}; 653 654 CHKERRQ(PetscMalloc2(2, &sizes, 3, &points)); 655 CHKERRQ(PetscArraycpy(sizes, triSizes_p2, 2)); 656 CHKERRQ(PetscArraycpy(points, triPoints_p2, 3));break;} 657 default: 658 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for triangular hybrid mesh on 2 procs", user->testNum); 659 } 660 } else if (dim == 2 && !cellSimplex && !cellHybrid && size == 2) { 661 switch (user->testNum) { 662 case 0: { 663 PetscInt quadSizes_p2[2] = {1, 1}; 664 PetscInt quadPoints_p2[2] = {0, 1}; 665 666 CHKERRQ(PetscMalloc2(2, &sizes, 2, &points)); 667 CHKERRQ(PetscArraycpy(sizes, quadSizes_p2, 2)); 668 CHKERRQ(PetscArraycpy(points, quadPoints_p2, 2));break;} 669 default: 670 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for quadrilateral mesh on 2 procs", user->testNum); 671 } 672 } else if (dim == 2 && !cellSimplex && cellHybrid && size == 2) { 673 switch (user->testNum) { 674 case 0: { 675 PetscInt quadSizes_p2[2] = {1, 2}; 676 PetscInt quadPoints_p2[3] = {0, 1, 2}; 677 678 CHKERRQ(PetscMalloc2(2, &sizes, 3, &points)); 679 CHKERRQ(PetscArraycpy(sizes, quadSizes_p2, 2)); 680 CHKERRQ(PetscArraycpy(points, quadPoints_p2, 3));break;} 681 default: 682 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for quadrilateral hybrid mesh on 2 procs", user->testNum); 683 } 684 } else if (dim == 3 && cellSimplex && !cellHybrid && size == 2) { 685 switch (user->testNum) { 686 case 0: { 687 PetscInt tetSizes_p2[2] = {1, 1}; 688 PetscInt tetPoints_p2[2] = {0, 1}; 689 690 CHKERRQ(PetscMalloc2(2, &sizes, 2, &points)); 691 CHKERRQ(PetscArraycpy(sizes, tetSizes_p2, 2)); 692 CHKERRQ(PetscArraycpy(points, tetPoints_p2, 2));break;} 693 case 1: { 694 PetscInt tetSizes_p2[2] = {1, 1}; 695 PetscInt tetPoints_p2[2] = {0, 1}; 696 697 CHKERRQ(PetscMalloc2(2, &sizes, 2, &points)); 698 CHKERRQ(PetscArraycpy(sizes, tetSizes_p2, 2)); 699 CHKERRQ(PetscArraycpy(points, tetPoints_p2, 2));break;} 700 default: 701 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for tetrahedral mesh on 2 procs", user->testNum); 702 } 703 } else if (dim == 3 && cellSimplex && cellHybrid && size == 2) { 704 switch (user->testNum) { 705 case 0: { 706 PetscInt tetSizes_p2[2] = {1, 2}; 707 PetscInt tetPoints_p2[3] = {0, 1, 2}; 708 709 CHKERRQ(PetscMalloc2(2, &sizes, 3, &points)); 710 CHKERRQ(PetscArraycpy(sizes, tetSizes_p2, 2)); 711 CHKERRQ(PetscArraycpy(points, tetPoints_p2, 3));break;} 712 case 1: { 713 PetscInt tetSizes_p2[2] = {3, 4}; 714 PetscInt tetPoints_p2[7] = {0, 3, 5, 1, 2, 4, 6}; 715 716 CHKERRQ(PetscMalloc2(2, &sizes, 7, &points)); 717 CHKERRQ(PetscArraycpy(sizes, tetSizes_p2, 2)); 718 CHKERRQ(PetscArraycpy(points, tetPoints_p2, 7));break;} 719 default: 720 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for tetrahedral hybrid mesh on 2 procs", user->testNum); 721 } 722 } else if (dim == 3 && !cellSimplex && !cellHybrid && size == 2) { 723 switch (user->testNum) { 724 case 0: { 725 PetscInt hexSizes_p2[2] = {1, 1}; 726 PetscInt hexPoints_p2[2] = {0, 1}; 727 728 CHKERRQ(PetscMalloc2(2, &sizes, 2, &points)); 729 CHKERRQ(PetscArraycpy(sizes, hexSizes_p2, 2)); 730 CHKERRQ(PetscArraycpy(points, hexPoints_p2, 2));break;} 731 default: 732 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for hexahedral mesh on 2 procs", user->testNum); 733 } 734 } else if (dim == 3 && !cellSimplex && cellHybrid && size == 2) { 735 switch (user->testNum) { 736 case 0: { 737 PetscInt hexSizes_p2[2] = {1, 1}; 738 PetscInt hexPoints_p2[2] = {0, 1}; 739 740 CHKERRQ(PetscMalloc2(2, &sizes, 2, &points)); 741 CHKERRQ(PetscArraycpy(sizes, hexSizes_p2, 2)); 742 CHKERRQ(PetscArraycpy(points, hexPoints_p2, 2));break;} 743 case 1: { 744 PetscInt hexSizes_p2[2] = {5, 4}; 745 PetscInt hexPoints_p2[9] = {3, 4, 5, 7, 8, 0, 1, 2, 6}; 746 747 CHKERRQ(PetscMalloc2(2, &sizes, 9, &points)); 748 CHKERRQ(PetscArraycpy(sizes, hexSizes_p2, 2)); 749 CHKERRQ(PetscArraycpy(points, hexPoints_p2, 9));break;} 750 default: 751 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for hexahedral hybrid mesh on 2 procs", user->testNum); 752 } 753 } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test partition"); 754 } 755 CHKERRQ(DMPlexGetPartitioner(*dm, &part)); 756 CHKERRQ(PetscPartitionerSetType(part, PETSCPARTITIONERSHELL)); 757 CHKERRQ(PetscPartitionerShellSetPartition(part, size, sizes, points)); 758 CHKERRQ(PetscFree2(sizes, points)); 759 } else { 760 PetscPartitioner part; 761 762 CHKERRQ(DMPlexGetPartitioner(*dm,&part)); 763 CHKERRQ(PetscPartitionerSetFromOptions(part)); 764 } 765 { 766 DM pdm = NULL; 767 768 CHKERRQ(DMPlexDistribute(*dm, 0, NULL, &pdm)); 769 if (pdm) { 770 CHKERRQ(DMViewFromOptions(pdm, NULL, "-dm_view")); 771 CHKERRQ(DMDestroy(dm)); 772 *dm = pdm; 773 } 774 } 775 CHKERRQ(DMPlexDistributeSetDefault(*dm, PETSC_FALSE)); 776 CHKERRQ(DMViewFromOptions(*dm, NULL, "-dm_view_pre")); 777 CHKERRQ(DMSetFromOptions(*dm)); 778 if (user->uninterpolate || user->reinterpolate) { 779 DM udm = NULL; 780 781 CHKERRQ(DMPlexUninterpolate(*dm, &udm)); 782 CHKERRQ(DMPlexCopyCoordinates(*dm, udm)); 783 CHKERRQ(DMDestroy(dm)); 784 *dm = udm; 785 } 786 if (user->reinterpolate) { 787 DM idm = NULL; 788 789 CHKERRQ(DMPlexInterpolate(*dm, &idm)); 790 CHKERRQ(DMPlexCopyCoordinates(*dm, idm)); 791 CHKERRQ(DMDestroy(dm)); 792 *dm = idm; 793 } 794 CHKERRQ(DMPlexDistributeSetDefault(*dm, PETSC_FALSE)); 795 CHKERRQ(PetscObjectSetName((PetscObject) *dm, "Hybrid Mesh")); 796 CHKERRQ(DMViewFromOptions(*dm, NULL, "-dm_view")); 797 CHKERRQ(PetscObjectSetOptionsPrefix((PetscObject) *dm, "hyb_")); 798 CHKERRQ(DMSetFromOptions(*dm)); 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 CHKERRQ(ProcessOptions(PETSC_COMM_WORLD, &user)); 810 CHKERRQ(CreateMesh(PETSC_COMM_WORLD, &user, &dm)); 811 CHKERRQ(DMDestroy(&dm)); 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