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