1 static char help[] = "Tests for creation of hybrid meshes\n\n"; 2 3 /* TODO 4 - Propagate hybridSize with distribution 5 - Test with multiple fault segments 6 - Test with embedded fault 7 - Test with multiple faults 8 - Move over all PyLith tests? 9 */ 10 11 #include <petscdmplex.h> 12 #include <petscds.h> 13 #include <petsc/private/dmpleximpl.h> 14 /* List of test meshes 15 16 Triangle 17 -------- 18 Test 0: 19 Two triangles sharing a face 20 21 4 22 / | \ 23 8 | 9 24 / | \ 25 2 0 7 1 5 26 \ | / 27 6 | 10 28 \ | / 29 3 30 31 should become two triangles separated by a zero-volume cell with 4 vertices 32 33 5--16--8 4--12--6 3 34 / | | \ / | | | \ 35 11 | | 12 9 | | | 4 36 / | | \ / | | | \ 37 3 0 10 2 14 1 6 2 0 8 1 10 6 0 1 38 \ | | / \ | | | / 39 9 | | 13 7 | | | 5 40 \ | | / \ | | | / 41 4--15--7 3--11--5 2 42 43 Test 1: 44 Four triangles sharing two faces which are oriented against each other 45 46 9 47 / \ 48 / \ 49 17 2 16 50 / \ 51 / \ 52 8-----15----5 53 \ /|\ 54 \ / | \ 55 18 3 12 | 14 56 \ / | \ 57 \ / | \ 58 4 0 11 1 7 59 \ | / 60 \ | / 61 10 | 13 62 \ | / 63 \|/ 64 6 65 66 Fault mesh 67 68 0 --> 0 69 1 --> 1 70 2 --> 2 71 3 --> 3 72 4 --> 5 73 5 --> 6 74 6 --> 8 75 7 --> 11 76 8 --> 15 77 78 2 79 | 80 6----8----4 81 | | 82 3 | 83 0-7-1 84 | 85 | 86 5 87 88 should become four triangles separated by two zero-volume cells with 4 vertices 89 90 11 91 / \ 92 / \ 93 / \ 94 22 2 21 95 / \ 96 / \ 97 10-----20------7 98 28 | 5 26/ \ 99 14----25----12 \ 100 \ /| |\ 101 \ / | | \ 102 23 3 17 | | 19 103 \ / | | \ 104 \ / | | \ 105 6 0 24 4 16 1 9 106 \ | | / 107 \ | | / 108 15 | | 18 109 \ | | / 110 \| |/ 111 13---8 112 27 113 114 Test 2: 115 Six triangles sharing one face 116 117 11-----12------13 118 | /|\ | 119 | 1 / | \ 4 | 120 | / | \ | 121 | / | \ | 122 | / | \ | 123 |/ | \| 124 9 2 | 5 10 125 |\ | /| 126 | \ | / | 127 | \ | / | 128 | \ | / | 129 | 0 \ | / 3 | 130 | \|/ | 131 6------7------8 132 133 Test 3: 134 This is Test 2 on two processes. After the fault, we have 135 136 6--12--7 7--20-10--16-8 137 | /| | |\ | 138 | 1 / | | | \ 1 | 139 13 11 | | | 17 15 140 | / | | | \ | 141 | / | | | \ | 142 |/ | | | \| 143 5 2 14 11 3 18 2 6 144 |\ | | | /| 145 | \ | | | / | 146 | \ | | | / | 147 10 9 | | | 14 13 148 | 0 \ | | | / 0 | 149 | \| | |/ | 150 3---8--4 4--19-9--12--5 151 152 Test 4: 153 This is Test 2 on six processes. After the fault, we have 154 155 Test 5: 156 157 Fault only on points 2 and 5: 158 159 6 160 / | \ 161 13 | 17 162 / 15 \ 163 7 0 | 1 9 164 |\ | /| 165 | 14 | 16 | 166 | \ | / | 167 18| 2 8 3 |21 168 | / | \ | 169 | 19 | 20 | 170 |/ | \| 171 10 4 | 5 12 172 \ 23 / 173 22 | 24 174 \ | / 175 11 176 177 Tetrahedron 178 ----------- 179 Test 0: 180 Two tets sharing a face 181 182 cell 5 _______ cell 183 0 / | \ \ 1 184 16 | 18 22 185 /8 19 10\ \ 186 2-15-|----4--21--6 187 \ 9| 7 / / 188 14 | 17 20 189 \ | / / 190 3------- 191 192 should become two tetrahedrons separated by a zero-volume cell with 3 faces/3 edges/6 vertices 193 194 cell 6 ___36___10______ cell 195 0 / | \ |\ \ 1 196 24 | 26 | 32 30 197 /12 27 14\ 33 \ \ 198 3-23-|----5--35-|---9--29--7 199 \ 13| 11/ |18 / / 200 22 | 25 | 31 28 201 \ | / |/ / 202 4----34----8------ 203 cell 2 204 205 In parallel, 206 207 cell 5 ___28____8 4______ cell 208 0 / | \ |\ |\ \ 0 209 19 | 21 | 24 | 13 6 11 210 /10 22 12\ 25 \ |8 \ \ 211 2-18-|----4--27-|---7 14 3--10--1 212 \ 11| 9 / |13 / | / / 213 17 | 20 | 23 | 12 5 9 214 \ | / |/ |/ / 215 3----26----6 2------ 216 cell 1 217 218 Test 1: 219 Four tets sharing two faces 220 221 Cells: 0-3,4-5 222 Vertices: 6-15 223 Faces: 16-29,30-34 224 Edges: 35-52,53-56 225 226 Quadrilateral 227 ------------- 228 Test 0: 229 Two quads sharing a face 230 231 5--10---4--14---7 232 | | | 233 11 0 9 1 13 234 | | | 235 2---8---3--12---6 236 237 should become two quads separated by a zero-volume cell with 4 vertices 238 239 6--13---5-20-10--17---8 5--10---4-14--7 4---7---2 240 | | | | | | | | | 241 14 0 12 2 18 1 16 11 0 9 1 12 8 0 6 242 | | | | | | | | | 243 3--11---4-19--9--15---7 2---8---3-13--6 3---5---1 244 245 Test 1: 246 247 Original mesh with 9 cells, 248 249 9-----10-----11-----12 250 | | || | 251 | | || | 252 | 0 | 1 || 2 | 253 | | || | 254 13-----14-----15-----16 255 | | || | 256 | | || | 257 | 3 | 4 || 5 | 258 | | || | 259 17-----18-----19=====20 260 | | | | 261 | | | | 262 | 6 | 7 | 8 | 263 | | | | 264 21-----22-----23-----24 265 266 After first fault, 267 268 12 ----13 ----14-28 ----15 269 | | | | | 270 | 0 | 1 | 9| 2 | 271 | | | | | 272 | | | | | 273 16 ----17 ----18-29 ----19 274 | | | | | 275 | 3 | 4 |10| 5 | 276 | | | | | 277 | | | | | 278 20 ----21-----22-30 ----23 279 | | | \--11- | 280 | 6 | 7 | 8 | 281 | | | | 282 | | | | 283 24 ----25 ----26--------27 284 285 After second fault, 286 287 14 ----15 ----16-30 ----17 288 | | | | | 289 | 0 | 1 | 9| 2 | 290 | | | | | 291 | | | | | 292 18 ----19 ----20-31 ----21 293 | | | | | 294 | 3 | 4 |10| 5 | 295 | | | | | 296 | | | | | 297 33 ----34-----24-32 ----25 298 | 12 | 13 / | \-11-- | 299 22 ----23---/ | | 300 | | | | 301 | 6 | 7 | 8 | 302 | | | | 303 | | | | 304 26 ----27 ----28--------29 305 306 Test 2: 307 Two quads sharing a face in parallel 308 309 4---7---3 2---8---4 310 | | | | 311 8 0 6 5 0 7 312 | | | | 313 1---5---2 1---6---3 314 315 should become two quads separated by a zero-volume cell with 4 vertices 316 317 4---7---3 3-14--7--11---5 318 | | | | | 319 8 0 6 8 1 12 0 10 320 | | | | | 321 1---5---2 2-13--6---9---4 322 323 Test 3: 324 Like Test 2, but with different partition 325 326 5--10---4-14--7 2---8---4 327 | | | | | 328 11 0 9 1 12 5 0 7 329 | | | | | 330 2---8---3-13--6 1---6---3 331 332 Hexahedron 333 ---------- 334 Test 0: 335 Two hexes sharing a face 336 337 cell 9-----31------8-----42------13 cell 338 0 /| /| /| 1 339 32 | 15 30| 21 41| 340 / | / | / | 341 6-----29------7-----40------12 | 342 | | 18 | | 24 | | 343 | 36 | 35 | 44 344 |19 | |17 | |23 | 345 33 | 16 34 | 22 43 | 346 | 5-----27--|---4-----39--|---11 347 | / | / | / 348 | 28 14 | 26 20 | 38 349 |/ |/ |/ 350 2-----25------3-----37------10 351 352 should become two hexes separated by a zero-volume cell with 8 vertices 353 354 cell 2 355 cell 10-----41------9-----62------18----52------14 cell 356 0 /| /| /| /| 1 357 42 | 20 40| 32 56| 26 51| 358 / | / | / | / | 359 7-----39------8-----61------17--|-50------13 | 360 | | 23 | | | | 29 | | 361 | 46 | 45 | 58 | 54 362 |24 | |22 | |30 | |28 | 363 43 | 21 44 | 57 | 27 53 | 364 | 6-----37--|---5-----60--|---16----49--|---12 365 | / | / | / | / 366 | 38 19 | 36 31 | 55 25 | 48 367 |/ |/ |/ |/ 368 3-----35------4-----59------15----47------11 369 370 In parallel, 371 372 cell 2 373 cell 9-----31------8-----44------13 8----20------4 cell 374 0 /| /| /| /| /| 1 375 32 | 15 30| 22 38| 24 | 10 19| 376 / | / | / | / | / | 377 6-----29------7-----43------12 | 7----18------3 | 378 | | 18 | | | | | | 13 | | 379 | 36 | 35 | 40 | 26 | 22 380 |19 | |17 | |20 | |14 | |12 | 381 33 | 16 34 | 39 | 25 | 11 21 | 382 | 5-----27--|---4-----42--|---11 | 6----17--|---2 383 | / | / | / | / | / 384 | 28 14 | 26 21 | 37 |23 9 | 16 385 |/ |/ |/ |/ |/ 386 2-----25------3-----41------10 5----15------1 387 388 Test 1: 389 390 */ 391 392 typedef struct { 393 PetscInt debug; /* The debugging level */ 394 PetscInt dim; /* The topological mesh dimension */ 395 PetscBool cellSimplex; /* Use simplices or hexes */ 396 PetscBool testPartition; /* Use a fixed partitioning for testing */ 397 PetscInt testNum; /* The particular mesh to test */ 398 PetscInt cohesiveFields; /* The number of cohesive fields */ 399 } AppCtx; 400 401 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 402 { 403 PetscFunctionBegin; 404 options->debug = 0; 405 options->dim = 2; 406 options->cellSimplex = PETSC_TRUE; 407 options->testPartition = PETSC_TRUE; 408 options->testNum = 0; 409 options->cohesiveFields = 1; 410 411 PetscOptionsBegin(comm, "", "Meshing Problem Options", "DMPLEX"); 412 PetscCall(PetscOptionsBoundedInt("-debug", "The debugging level", "ex5.c", options->debug, &options->debug, NULL,0)); 413 PetscCall(PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex5.c", options->dim, &options->dim, NULL,1,3)); 414 PetscCall(PetscOptionsBool("-cell_simplex", "Use simplices if true, otherwise hexes", "ex5.c", options->cellSimplex, &options->cellSimplex, NULL)); 415 PetscCall(PetscOptionsBool("-test_partition", "Use a fixed partition for testing", "ex5.c", options->testPartition, &options->testPartition, NULL)); 416 PetscCall(PetscOptionsBoundedInt("-test_num", "The particular mesh to test", "ex5.c", options->testNum, &options->testNum, NULL,0)); 417 PetscCall(PetscOptionsBoundedInt("-cohesive_fields", "The number of cohesive fields", "ex5.c", options->cohesiveFields, &options->cohesiveFields, NULL, 0)); 418 PetscOptionsEnd(); 419 PetscFunctionReturn(0); 420 } 421 422 static PetscErrorCode CreateSimplex_2D(MPI_Comm comm, PetscInt testNum, DM *dm) 423 { 424 DM idm; 425 PetscInt p; 426 PetscMPIInt rank; 427 428 PetscFunctionBegin; 429 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 430 if (rank == 0) { 431 switch (testNum) { 432 case 0: 433 { 434 PetscInt numPoints[2] = {4, 2}; 435 PetscInt coneSize[6] = {3, 3, 0, 0, 0, 0}; 436 PetscInt cones[6] = {2, 3, 4, 5, 4, 3}; 437 PetscInt coneOrientations[6] = {0, 0, 0, 0, 0, 0}; 438 PetscScalar vertexCoords[8] = {-0.5, 0.5, 0.0, 0.0, 0.0, 1.0, 0.5, 0.5}; 439 PetscInt markerPoints[8] = {2, 1, 3, 1, 4, 1, 5, 1}; 440 PetscInt faultPoints[2] = {3, 4}; 441 442 PetscCall(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 443 for (p = 0; p < 4; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p*2], markerPoints[p*2+1])); 444 for (p = 0; p < 2; ++p) PetscCall(DMSetLabelValue(*dm, "fault", faultPoints[p], 1)); 445 PetscCall(DMSetLabelValue(*dm, "material", 0, 1)); 446 PetscCall(DMSetLabelValue(*dm, "material", 1, 2)); 447 } 448 break; 449 case 1: 450 { 451 PetscInt numPoints[2] = {6, 4}; 452 PetscInt coneSize[10] = {3, 3, 3, 3, 0, 0, 0, 0, 0, 0}; 453 PetscInt cones[12] = {4, 6, 5, 5, 6, 7, 8, 5, 9, 8, 4, 5}; 454 PetscInt coneOrientations[12] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 455 PetscScalar vertexCoords[12] = {-1.0, 0.0, 0.0, 1.0, 0.0, -1.0, 1.0, 0.0, -2.0, 1.0, -1.0, 2.0}; 456 PetscInt markerPoints[6] = {4, 1, 6, 1, 8, 1}; 457 PetscInt faultPoints[3] = {5, 6, 8}; 458 459 PetscCall(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 460 for (p = 0; p < 3; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p*2], markerPoints[p*2+1])); 461 for (p = 0; p < 3; ++p) PetscCall(DMSetLabelValue(*dm, "fault", faultPoints[p], 1)); 462 PetscCall(DMSetLabelValue(*dm, "material", 0, 1));PetscCall(DMSetLabelValue(*dm, "material", 3, 1)); 463 PetscCall(DMSetLabelValue(*dm, "material", 1, 2));PetscCall(DMSetLabelValue(*dm, "material", 2, 2)); 464 } 465 break; 466 case 2: 467 case 3: 468 case 4: 469 { 470 PetscInt numPoints[2] = {8, 6}; 471 PetscInt coneSize[14] = {3, 3, 3, 3, 3, 3, 0, 0, 0, 0, 0, 0, 0, 0}; 472 PetscInt cones[18] = {6, 7, 9, 9, 12, 11, 7, 12, 9, 473 7, 8, 10, 10, 13, 12, 7, 10, 12}; 474 PetscInt coneOrientations[18] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 475 PetscScalar vertexCoords[16] = {-1., -1., 0., -1., 1., -1., -1., 0., 1., 0., 476 -1., 1., 0., 1., 1., 1.,}; 477 PetscInt markerPoints[16] = {6, 1, 7, 1, 8, 1, 9, 1, 10, 1, 11, 1, 12, 1, 13, 1}; 478 PetscInt faultPoints[2] = {7, 12}; 479 480 PetscCall(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 481 PetscCall(DMSetLabelValue(*dm, "material", 0, 1)); 482 PetscCall(DMSetLabelValue(*dm, "material", 1, 1)); 483 PetscCall(DMSetLabelValue(*dm, "material", 2, 1)); 484 PetscCall(DMSetLabelValue(*dm, "material", 3, 2)); 485 PetscCall(DMSetLabelValue(*dm, "material", 4, 2)); 486 PetscCall(DMSetLabelValue(*dm, "material", 5, 2)); 487 for (p = 0; p < 8; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p*2], markerPoints[p*2+1])); 488 if (testNum == 2) for (p = 0; p < 2; ++p) PetscCall(DMSetLabelValue(*dm, "fault", faultPoints[p], 1)); 489 if (testNum == 3 || testNum == 4) for (p = 0; p < 2; ++p) PetscCall(DMSetLabelValue(*dm, "pfault", faultPoints[p], 1)); 490 } 491 break; 492 case 5: 493 { 494 PetscInt numPoints[2] = {7, 6}; 495 PetscInt coneSize[13] = {3, 3, 3, 3, 3, 3, 0, 0, 0, 0, 0, 0, 0}; 496 PetscInt cones[18] = {6, 7, 8, 8, 9, 6, 7, 10, 8, 9, 8, 12, 8, 10, 11, 11, 12, 8}; 497 PetscInt coneOrientations[18] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 498 PetscScalar vertexCoords[14] = {0.0, 2.0, -1.0, 1.0, 0.0, 0.0, 1.0, 1.0, -1.0, -1.0, 0.0, -2.0, 1.0, -1.0}; 499 PetscInt faultPoints[2] = {8, 11}; 500 PetscInt faultBdPoints[1] = {8}; 501 502 PetscCall(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 503 for (p = 0; p < 2; ++p) PetscCall(DMSetLabelValue(*dm, "fault", faultPoints[p], 1)); 504 for (p = 0; p < 1; ++p) PetscCall(DMSetLabelValue(*dm, "faultBd", faultBdPoints[p], 1)); 505 PetscCall(DMSetLabelValue(*dm, "material", 0, 0));PetscCall(DMSetLabelValue(*dm, "material", 2, 0));PetscCall(DMSetLabelValue(*dm, "material", 4, 0)); 506 PetscCall(DMSetLabelValue(*dm, "material", 1, 2));PetscCall(DMSetLabelValue(*dm, "material", 3, 2));PetscCall(DMSetLabelValue(*dm, "material", 5, 2)); 507 } 508 break; 509 default: 510 SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum); 511 } 512 } else { 513 PetscInt numPoints[3] = {0, 0, 0}; 514 515 PetscCall(DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL)); 516 if (testNum == 3 || testNum == 4) PetscCall(DMCreateLabel(*dm, "pfault")); 517 else PetscCall(DMCreateLabel(*dm, "fault")); 518 } 519 PetscCall(DMPlexInterpolate(*dm, &idm)); 520 PetscCall(DMViewFromOptions(idm, NULL, "-in_dm_view")); 521 PetscCall(DMDestroy(dm)); 522 *dm = idm; 523 PetscFunctionReturn(0); 524 } 525 526 static PetscErrorCode CreateSimplex_3D(MPI_Comm comm, AppCtx *user, DM dm) 527 { 528 PetscInt depth = 3, testNum = user->testNum, p; 529 PetscMPIInt rank; 530 531 PetscFunctionBegin; 532 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 533 if (rank == 0) { 534 switch (testNum) { 535 case 0: 536 { 537 PetscInt numPoints[4] = {5, 7, 9, 2}; 538 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}; 539 PetscInt cones[47] = {7, 8, 9, 10, 11, 10, 13, 12, 15, 17, 14, 16, 18, 15, 14, 19, 16, 17, 18, 19, 17, 21, 20, 18, 22, 21, 22, 19, 20, 2, 3, 2, 4, 2, 5, 3, 4, 4, 5, 5, 3, 3, 6, 4, 6, 5, 6}; 540 PetscInt coneOrientations[47] = {0, 0, 0, 0, 0, -2, 2, 2, 0, -1, -1, 0, -1, -1, 0, -1, -1, 0, 0, 0, 0, 0, -1, 0, 0, -1, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 541 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}; 542 PetscInt markerPoints[20] = {2, 1, 3, 1, 4, 1, 5, 1, 14, 1, 15, 1, 16, 1, 17, 1, 18, 1, 19, 1}; 543 PetscInt faultPoints[3] = {3, 4, 5}; 544 545 PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 546 for (p = 0; p < 10; ++p) { 547 PetscCall(DMSetLabelValue(dm, "marker", markerPoints[p*2], markerPoints[p*2+1])); 548 } 549 for (p = 0; p < 3; ++p) { 550 PetscCall(DMSetLabelValue(dm, "fault", faultPoints[p], 1)); 551 } 552 PetscCall(DMSetLabelValue(dm, "material", 0, 1)); 553 PetscCall(DMSetLabelValue(dm, "material", 1, 2)); 554 } 555 break; 556 case 1: 557 { 558 PetscInt numPoints[4] = {6, 13, 12, 4}; 559 PetscInt coneSize[35] = {4, 4, 4, 4, 0, 0, 0, 0, 0, 0, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2}; 560 PetscInt cones[78] = {10, 11, 12, 13, 10, 15, 16, 14, 17, 18, 14, 19, 20, 13, 19, 21, 22, 23, 24, 25, 26, 22, 24, 27, 25, 23, 26, 27, 28, 29, 23, 24, 30, 28, 22, 29, 30, 31, 32, 28, 29, 33, 31, 32, 33, 23, 26, 34, 33, 34, 27, 32, 6, 5, 5, 7, 7, 6, 6, 4, 4, 5, 7, 4, 7, 9, 9, 5, 6, 9, 9, 8, 8, 7, 5, 8, 4, 8}; 561 PetscInt coneOrientations[78] = { 0, 0, 0, 0, -2, 1, 0, 2, 0, 0, -3, 0, 0, -3, -1, 0, 0, 0, 0, 0, 0, -1, -1, 0, -1, -1, -1, -1, 0, 0, 0, 0, 0, -1, 0, -1, -1, 0, 0, 0, 0, 0, -1, -1, -1, 0, -1, 0, -1, -1, -1, -1, 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}; 562 PetscScalar vertexCoords[18] = {-1.0, 0.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 1.0, 0.0, 0.0, 0.0, -1.0, 1.0, 0.0, 0.0}; 563 PetscInt markerPoints[14] = {5, 1, 6, 1, 7, 1, 10, 1, 22, 1, 23, 1, 24, 1}; 564 PetscInt faultPoints[4] = {5, 6, 7, 8}; 565 566 PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 567 for (p = 0; p < 7; ++p) { 568 PetscCall(DMSetLabelValue(dm, "marker", markerPoints[p*2], markerPoints[p*2+1])); 569 } 570 for (p = 0; p < 4; ++p) { 571 PetscCall(DMSetLabelValue(dm, "fault", faultPoints[p], 1)); 572 } 573 PetscCall(DMSetLabelValue(dm, "material", 0, 1)); 574 PetscCall(DMSetLabelValue(dm, "material", 1, 2)); 575 } 576 break; 577 default: 578 SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum); 579 } 580 } else { 581 PetscInt numPoints[4] = {0, 0, 0, 0}; 582 583 PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, NULL, NULL, NULL, NULL)); 584 PetscCall(DMCreateLabel(dm, "fault")); 585 } 586 PetscFunctionReturn(0); 587 } 588 589 static PetscErrorCode CreateQuad_2D(MPI_Comm comm, PetscInt testNum, DM *dm) 590 { 591 DM idm; 592 PetscInt p; 593 PetscMPIInt rank; 594 595 PetscFunctionBegin; 596 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 597 if (rank == 0) { 598 switch (testNum) { 599 case 0: 600 case 2: 601 case 3: 602 { 603 PetscInt numPoints[2] = {6, 2}; 604 PetscInt coneSize[8] = {4, 4, 0, 0, 0, 0, 0, 0}; 605 PetscInt cones[8] = {2, 3, 4, 5, 3, 6, 7, 4}; 606 PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0}; 607 PetscScalar vertexCoords[12] = {-0.5, 0.0, 0.0, 0.0, 0.0, 1.0, -0.5, 1.0, 0.5, 0.0, 0.5, 1.0}; 608 PetscInt markerPoints[12] = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1, 7, 1}; 609 PetscInt faultPoints[2] = {3, 4}; 610 611 PetscCall(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 612 for (p = 0; p < 6; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p*2], markerPoints[p*2+1])); 613 if (testNum == 0) for (p = 0; p < 2; ++p) PetscCall(DMSetLabelValue(*dm, "fault", faultPoints[p], 1)); 614 if (testNum == 2 || testNum == 3) for (p = 0; p < 2; ++p) PetscCall(DMSetLabelValue(*dm, "pfault", faultPoints[p], 1)); 615 PetscCall(DMSetLabelValue(*dm, "material", 0, 1)); 616 PetscCall(DMSetLabelValue(*dm, "material", 1, 2)); 617 } 618 break; 619 case 1: 620 { 621 PetscInt numPoints[2] = {16, 9}; 622 PetscInt coneSize[25] = {4, 4, 4, 4, 4, 4, 4, 4, 4, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 623 PetscInt cones[36] = {9, 13, 14, 10, 624 10, 14, 15, 11, 625 11, 15, 16, 12, 626 13, 17, 18, 14, 627 14, 18, 19, 15, 628 15, 19, 20, 16, 629 17, 21, 22, 18, 630 18, 22, 23, 19, 631 19, 23, 24, 20}; 632 PetscInt coneOrientations[36] = {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}; 633 PetscScalar vertexCoords[32] = {-3.0, 3.0, -1.0, 3.0, 1.0, 3.0, 3.0, 3.0, -3.0, 1.0, -1.0, 1.0, 1.0, 1.0, 3.0, 1.0, 634 -3.0, -1.0, -1.0, -1.0, 1.0, -1.0, 3.0, -1.0, -3.0, -3.0, -1.0, -3.0, 1.0, -3.0, 3.0, -3.0}; 635 PetscInt faultPoints[4] = {11, 15, 19, 20}; 636 PetscInt fault2Points[2] = {17, 18}; 637 638 PetscCall(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 639 for (p = 0; p < 4; ++p) PetscCall(DMSetLabelValue(*dm, "fault", faultPoints[p], 1)); 640 for (p = 3; p < 4; ++p) PetscCall(DMSetLabelValue(*dm, "faultBd", faultPoints[p], 1)); 641 for (p = 0; p < 2; ++p) PetscCall(DMSetLabelValue(*dm, "fault2", fault2Points[p], 1)); 642 PetscCall(DMSetLabelValue(*dm, "material", 0, 1)); 643 PetscCall(DMSetLabelValue(*dm, "material", 1, 1)); 644 PetscCall(DMSetLabelValue(*dm, "material", 2, 1)); 645 PetscCall(DMSetLabelValue(*dm, "material", 3, 1)); 646 PetscCall(DMSetLabelValue(*dm, "material", 4, 1)); 647 PetscCall(DMSetLabelValue(*dm, "material", 5, 2)); 648 PetscCall(DMSetLabelValue(*dm, "material", 6, 2)); 649 PetscCall(DMSetLabelValue(*dm, "material", 7, 2)); 650 PetscCall(DMSetLabelValue(*dm, "material", 8, 2)); 651 } 652 break; 653 default: 654 SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum); 655 } 656 } else { 657 PetscInt numPoints[3] = {0, 0, 0}; 658 659 PetscCall(DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL)); 660 if (testNum == 2 || testNum == 3) PetscCall(DMCreateLabel(*dm, "pfault")); 661 else PetscCall(DMCreateLabel(*dm, "fault")); 662 } 663 PetscCall(DMPlexInterpolate(*dm, &idm)); 664 PetscCall(DMViewFromOptions(idm, NULL, "-in_dm_view")); 665 PetscCall(DMDestroy(dm)); 666 *dm = idm; 667 PetscFunctionReturn(0); 668 } 669 670 static PetscErrorCode CreateHex_3D(MPI_Comm comm, PetscInt testNum, DM *dm) 671 { 672 DM idm; 673 PetscInt depth = 3, p; 674 PetscMPIInt rank; 675 676 PetscFunctionBegin; 677 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 678 if (rank == 0) { 679 switch (testNum) { 680 case 0: 681 { 682 PetscInt numPoints[2] = {12, 2}; 683 PetscInt coneSize[14] = {8,8, 0,0,0,0,0,0,0,0,0,0,0,0}; 684 PetscInt cones[16] = {2,5,4,3,6,7,8,9, 3,4,11,10,7,12,13,8}; 685 PetscInt coneOrientations[16] = {0,0,0,0,0,0,0,0, 0,0, 0,0 ,0, 0, 0,0}; 686 PetscScalar vertexCoords[36] = {-0.5,0.0,0.0, 0.0,0.0,0.0, 0.0,1.0,0.0, -0.5,1.0,0.0, 687 -0.5,0.0,1.0, 0.0,0.0,1.0, 0.0,1.0,1.0, -0.5,1.0,1.0, 688 0.5,0.0,0.0, 0.5,1.0,0.0, 0.5,0.0,1.0, 0.5,1.0,1.0}; 689 PetscInt markerPoints[52] = {2,1,3,1,4,1,5,1,6,1,7,1,8,1,9,1}; 690 PetscInt faultPoints[4] = {3, 4, 7, 8}; 691 692 PetscCall(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 693 PetscCall(DMPlexInterpolate(*dm, &idm)); 694 for (p = 0; p < 8; ++p) PetscCall(DMSetLabelValue(idm, "marker", markerPoints[p*2], markerPoints[p*2+1])); 695 for (p = 0; p < 4; ++p) PetscCall(DMSetLabelValue(idm, "fault", faultPoints[p], 1)); 696 PetscCall(DMSetLabelValue(*dm, "material", 0, 1)); 697 PetscCall(DMSetLabelValue(*dm, "material", 1, 2)); 698 } 699 break; 700 case 1: 701 { 702 /* Cell Adjacency Graph: 703 0 -- { 8, 13, 21, 24} --> 1 704 0 -- {20, 21, 23, 24} --> 5 F 705 1 -- {10, 15, 21, 24} --> 2 706 1 -- {13, 14, 15, 24} --> 6 707 2 -- {21, 22, 24, 25} --> 4 F 708 3 -- {21, 24, 30, 35} --> 4 709 3 -- {21, 24, 28, 33} --> 5 710 */ 711 PetscInt numPoints[2] = {30, 7}; 712 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}; 713 PetscInt cones[56] = { 8, 21, 20, 7, 13, 12, 23, 24, 714 14, 15, 10, 9, 13, 8, 21, 24, 715 15, 16, 11, 10, 24, 21, 22, 25, 716 30, 29, 28, 21, 35, 24, 33, 34, 717 24, 21, 30, 35, 25, 36, 31, 22, 718 27, 20, 21, 28, 32, 33, 24, 23, 719 15, 24, 13, 14, 19, 18, 17, 26}; 720 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}; 721 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, 722 -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, 723 -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, 724 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, 725 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}; 726 PetscInt faultPoints[6] = {20, 21, 22, 23, 24, 25}; 727 728 PetscCall(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 729 PetscCall(DMPlexInterpolate(*dm, &idm)); 730 for (p = 0; p < 6; ++p) PetscCall(DMSetLabelValue(idm, "fault", faultPoints[p], 1)); 731 PetscCall(DMSetLabelValue(*dm, "material", 0, 1)); 732 PetscCall(DMSetLabelValue(*dm, "material", 1, 1)); 733 PetscCall(DMSetLabelValue(*dm, "material", 2, 1)); 734 PetscCall(DMSetLabelValue(*dm, "material", 3, 2)); 735 PetscCall(DMSetLabelValue(*dm, "material", 4, 2)); 736 PetscCall(DMSetLabelValue(*dm, "material", 5, 2)); 737 PetscCall(DMSetLabelValue(*dm, "material", 6, 2)); 738 } 739 break; 740 case 2: 741 { 742 /* Buried fault edge */ 743 PetscInt numPoints[2] = {18, 4}; 744 PetscInt coneSize[22] = {8,8,8,8, 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; 745 PetscInt cones[32] = { 4, 5, 8, 7, 13, 16, 17, 14, 746 5, 6, 9, 8, 14, 17, 18, 15, 747 7, 8, 11, 10, 16, 19, 20, 17, 748 8, 9, 12, 11, 17, 20, 21, 18}; 749 PetscInt coneOrientations[32] = {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}; 750 PetscScalar vertexCoords[54] = {-2.0, -2.0, 0.0, -2.0, 0.0, 0.0, -2.0, 2.0, 0.0, 0.0, -2.0, 0.0, 0.0, 0.0, 0.0, 0.0, 2.0, 0.0, 751 2.0, -2.0, 0.0, 2.0, 0.0, 0.0, 2.0, 2.0, 0.0, -2.0, -2.0, 2.0, -2.0, 0.0, 2.0, -2.0, 2.0, 2.0, 752 0.0, -2.0, 2.0, 0.0, 0.0, 2.0, 0.0, 2.0, 2.0, 2.0, -2.0, 2.0, 2.0, 0.0, 2.0, 2.0, 2.0, 2.0}; 753 PetscInt faultPoints[4] = {7, 8, 16, 17}; 754 755 PetscCall(DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords)); 756 PetscCall(DMPlexInterpolate(*dm, &idm)); 757 for (p = 0; p < 4; ++p) PetscCall(DMSetLabelValue(idm, "fault", faultPoints[p], 1)); 758 PetscCall(DMSetLabelValue(*dm, "material", 0, 1)); 759 PetscCall(DMSetLabelValue(*dm, "material", 1, 1)); 760 PetscCall(DMSetLabelValue(*dm, "material", 2, 2)); 761 PetscCall(DMSetLabelValue(*dm, "material", 3, 2)); 762 } 763 break; 764 default: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum); 765 } 766 } else { 767 PetscInt numPoints[4] = {0, 0, 0, 0}; 768 769 PetscCall(DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL)); 770 PetscCall(DMPlexInterpolate(*dm, &idm)); 771 PetscCall(DMCreateLabel(idm, "fault")); 772 } 773 PetscCall(DMViewFromOptions(idm, NULL, "-in_dm_view")); 774 PetscCall(DMDestroy(dm)); 775 *dm = idm; 776 PetscFunctionReturn(0); 777 } 778 779 static PetscErrorCode CreateFaultLabel(DM dm) 780 { 781 DMLabel label; 782 PetscInt dim, h, pStart, pEnd, pMax, p; 783 784 PetscFunctionBegin; 785 PetscCall(DMGetDimension(dm, &dim)); 786 PetscCall(DMCreateLabel(dm, "cohesive")); 787 PetscCall(DMGetLabel(dm, "cohesive", &label)); 788 for (h = 0; h <= dim; ++h) { 789 PetscCall(DMPlexGetSimplexOrBoxCells(dm, h, NULL, &pMax)); 790 PetscCall(DMPlexGetHeightStratum(dm, h, &pStart, &pEnd)); 791 for (p = pMax; p < pEnd; ++p) PetscCall(DMLabelSetValue(label, p, 1)); 792 } 793 PetscFunctionReturn(0); 794 } 795 796 static PetscErrorCode CreateDiscretization(DM dm, AppCtx *user) 797 { 798 PetscFE fe; 799 DMLabel fault; 800 PetscInt dim, Ncf = user->cohesiveFields, f; 801 802 PetscFunctionBegin; 803 PetscCall(DMGetDimension(dm, &dim)); 804 PetscCall(DMGetLabel(dm, "cohesive", &fault)); 805 PetscCall(DMLabelView(fault, PETSC_VIEWER_STDOUT_WORLD)); 806 807 PetscCall(PetscFECreateDefault(PETSC_COMM_SELF, dim, dim, user->cellSimplex, "displacement_", PETSC_DETERMINE, &fe)); 808 PetscCall(PetscFESetName(fe, "displacement")); 809 PetscCall(DMAddField(dm, NULL, (PetscObject) fe)); 810 PetscCall(PetscFEDestroy(&fe)); 811 812 if (Ncf > 0) { 813 PetscCall(PetscFECreateDefault(PETSC_COMM_SELF, dim-1, dim, user->cellSimplex, "faulttraction_", PETSC_DETERMINE, &fe)); 814 PetscCall(PetscFESetName(fe, "fault traction")); 815 PetscCall(DMAddField(dm, fault, (PetscObject) fe)); 816 PetscCall(PetscFEDestroy(&fe)); 817 } 818 for (f = 1; f < Ncf; ++f) { 819 char name[256], opt[256]; 820 821 PetscCall(PetscSNPrintf(name, 256, "fault field %" PetscInt_FMT, f)); 822 PetscCall(PetscSNPrintf(opt, 256, "faultfield_%" PetscInt_FMT "_", f)); 823 PetscCall(PetscFECreateDefault(PETSC_COMM_SELF, dim-1, dim, user->cellSimplex, opt, PETSC_DETERMINE, &fe)); 824 PetscCall(PetscFESetName(fe, name)); 825 PetscCall(DMAddField(dm, fault, (PetscObject) fe)); 826 PetscCall(PetscFEDestroy(&fe)); 827 } 828 829 PetscCall(DMCreateDS(dm)); 830 PetscFunctionReturn(0); 831 } 832 833 static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm) 834 { 835 PetscInt dim = user->dim; 836 PetscBool cellSimplex = user->cellSimplex, hasFault, hasFault2, hasParallelFault; 837 PetscMPIInt rank, size; 838 DMLabel matLabel; 839 840 PetscFunctionBegin; 841 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 842 PetscCallMPI(MPI_Comm_size(comm, &size)); 843 PetscCall(DMCreate(comm, dm)); 844 PetscCall(DMSetType(*dm, DMPLEX)); 845 PetscCall(DMSetDimension(*dm, dim)); 846 switch (dim) { 847 case 2: 848 if (cellSimplex) { 849 PetscCall(CreateSimplex_2D(comm, user->testNum, dm)); 850 } else { 851 PetscCall(CreateQuad_2D(comm, user->testNum, dm)); 852 } 853 break; 854 case 3: 855 if (cellSimplex) { 856 PetscCall(CreateSimplex_3D(comm, user, *dm)); 857 } else { 858 PetscCall(CreateHex_3D(comm, user->testNum, dm)); 859 } 860 break; 861 default: 862 SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "Cannot make hybrid meshes for dimension %" PetscInt_FMT, dim); 863 } 864 PetscCall(PetscObjectSetOptionsPrefix((PetscObject) *dm, "orig_")); 865 PetscCall(DMPlexDistributeSetDefault(*dm, PETSC_FALSE)); 866 PetscCall(DMSetFromOptions(*dm)); 867 PetscCall(DMGetLabel(*dm, "material", &matLabel)); 868 if (matLabel) { 869 PetscCall(DMPlexLabelComplete(*dm, matLabel)); 870 } 871 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 872 PetscCall(DMHasLabel(*dm, "fault", &hasFault)); 873 if (hasFault) { 874 DM dmHybrid = NULL, dmInterface = NULL; 875 DMLabel faultLabel, faultBdLabel, hybridLabel, splitLabel; 876 877 PetscCall(DMGetLabel(*dm, "fault", &faultLabel)); 878 PetscCall(DMGetLabel(*dm, "faultBd", &faultBdLabel)); 879 PetscCall(DMPlexCreateHybridMesh(*dm, faultLabel, faultBdLabel, 1, &hybridLabel, &splitLabel, &dmInterface, &dmHybrid)); 880 PetscCall(DMLabelView(hybridLabel, PETSC_VIEWER_STDOUT_WORLD)); 881 PetscCall(DMLabelDestroy(&hybridLabel)); 882 PetscCall(DMLabelView(splitLabel, PETSC_VIEWER_STDOUT_WORLD)); 883 PetscCall(DMLabelDestroy(&splitLabel)); 884 PetscCall(DMViewFromOptions(dmInterface, NULL, "-dm_interface_view")); 885 PetscCall(DMDestroy(&dmInterface)); 886 PetscCall(DMDestroy(dm)); 887 *dm = dmHybrid; 888 } 889 PetscCall(DMHasLabel(*dm, "fault2", &hasFault2)); 890 if (hasFault2) { 891 DM dmHybrid = NULL; 892 DMLabel faultLabel, faultBdLabel, hybridLabel; 893 894 PetscCall(PetscObjectSetOptionsPrefix((PetscObject) *dm, "faulted_")); 895 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view_pre")); 896 PetscCall(DMPlexDistributeSetDefault(*dm, PETSC_FALSE)); 897 PetscCall(DMSetFromOptions(*dm)); 898 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 899 PetscCall(DMGetLabel(*dm, "fault2", &faultLabel)); 900 PetscCall(DMGetLabel(*dm, "fault2Bd", &faultBdLabel)); 901 PetscCall(DMPlexCreateHybridMesh(*dm, faultLabel, faultBdLabel, 1, &hybridLabel, NULL, NULL, &dmHybrid)); 902 PetscCall(DMLabelView(hybridLabel, PETSC_VIEWER_STDOUT_WORLD)); 903 PetscCall(DMLabelDestroy(&hybridLabel)); 904 PetscCall(DMDestroy(dm)); 905 *dm = dmHybrid; 906 } 907 if (user->testPartition && size > 1) { 908 PetscPartitioner part; 909 PetscInt *sizes = NULL; 910 PetscInt *points = NULL; 911 912 if (rank == 0) { 913 if (dim == 2 && cellSimplex && size == 2) { 914 switch (user->testNum) { 915 case 0: { 916 PetscInt triSizes_p2[2] = {1, 2}; 917 PetscInt triPoints_p2[3] = {0, 1, 2}; 918 919 PetscCall(PetscMalloc2(2, &sizes, 3, &points)); 920 PetscCall(PetscArraycpy(sizes, triSizes_p2, 2)); 921 PetscCall(PetscArraycpy(points, triPoints_p2, 3));break;} 922 case 3: { 923 PetscInt triSizes_p2[2] = {3, 3}; 924 PetscInt triPoints_p2[6] = {0, 1, 2, 3, 4, 5}; 925 926 PetscCall(PetscMalloc2(2, &sizes, 6, &points)); 927 PetscCall(PetscArraycpy(sizes, triSizes_p2, 2)); 928 PetscCall(PetscArraycpy(points, triPoints_p2, 6));break;} 929 default: 930 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for triangular mesh on 2 procs", user->testNum); 931 } 932 } else if (dim == 2 && cellSimplex && size == 6) { 933 switch (user->testNum) { 934 case 4: { 935 PetscInt triSizes_p6[6] = {1, 1, 1, 1, 1, 1}; 936 PetscInt triPoints_p6[6] = {0, 1, 2, 3, 4, 5}; 937 938 PetscCall(PetscMalloc2(6, &sizes, 6, &points)); 939 PetscCall(PetscArraycpy(sizes, triSizes_p6, 6)); 940 PetscCall(PetscArraycpy(points, triPoints_p6, 6));break;} 941 default: 942 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for triangular mesh on 6 procs", user->testNum); 943 } 944 } else if (dim == 2 && !cellSimplex && size == 2) { 945 switch (user->testNum) { 946 case 0: { 947 PetscInt quadSizes_p2[2] = {1, 2}; 948 PetscInt quadPoints_p2[3] = {0, 1, 2}; 949 950 PetscCall(PetscMalloc2(2, &sizes, 3, &points)); 951 PetscCall(PetscArraycpy(sizes, quadSizes_p2, 2)); 952 PetscCall(PetscArraycpy(points, quadPoints_p2, 3));break;} 953 case 2: { 954 PetscInt quadSizes_p2[2] = {1, 1}; 955 PetscInt quadPoints_p2[2] = {0, 1}; 956 957 PetscCall(PetscMalloc2(2, &sizes, 2, &points)); 958 PetscCall(PetscArraycpy(sizes, quadSizes_p2, 2)); 959 PetscCall(PetscArraycpy(points, quadPoints_p2, 2));break;} 960 case 3: { 961 PetscInt quadSizes_p2[2] = {1, 1}; 962 PetscInt quadPoints_p2[2] = {1, 0}; 963 964 PetscCall(PetscMalloc2(2, &sizes, 2, &points)); 965 PetscCall(PetscArraycpy(sizes, quadSizes_p2, 2)); 966 PetscCall(PetscArraycpy(points, quadPoints_p2, 2));break;} 967 default: 968 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for quadrilateral mesh on 2 procs", user->testNum); 969 } 970 } else if (dim == 3 && cellSimplex && size == 2) { 971 switch (user->testNum) { 972 case 0: { 973 PetscInt tetSizes_p2[2] = {1, 2}; 974 PetscInt tetPoints_p2[3] = {0, 1, 2}; 975 976 PetscCall(PetscMalloc2(2, &sizes, 3, &points)); 977 PetscCall(PetscArraycpy(sizes, tetSizes_p2, 2)); 978 PetscCall(PetscArraycpy(points, tetPoints_p2, 3));break;} 979 default: 980 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for teterehedral mesh on 2 procs", user->testNum); 981 } 982 } else if (dim == 3 && !cellSimplex && size == 2) { 983 switch (user->testNum) { 984 case 0: { 985 PetscInt hexSizes_p2[2] = {1, 2}; 986 PetscInt hexPoints_p2[3] = {0, 1, 2}; 987 988 PetscCall(PetscMalloc2(2, &sizes, 3, &points)); 989 PetscCall(PetscArraycpy(sizes, hexSizes_p2, 2)); 990 PetscCall(PetscArraycpy(points, hexPoints_p2, 3));break;} 991 default: 992 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for hexahedral mesh on 2 procs", user->testNum); 993 } 994 } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test partition"); 995 } 996 PetscCall(DMPlexGetPartitioner(*dm, &part)); 997 PetscCall(PetscPartitionerSetType(part, PETSCPARTITIONERSHELL)); 998 PetscCall(PetscPartitionerShellSetPartition(part, size, sizes, points)); 999 PetscCall(PetscFree2(sizes, points)); 1000 } 1001 { 1002 DM pdm = NULL; 1003 1004 /* Distribute mesh over processes */ 1005 PetscCall(DMPlexDistribute(*dm, 0, NULL, &pdm)); 1006 if (pdm) { 1007 PetscCall(DMViewFromOptions(pdm, NULL, "-dm_view")); 1008 PetscCall(DMDestroy(dm)); 1009 *dm = pdm; 1010 } 1011 } 1012 PetscCall(DMHasLabel(*dm, "pfault", &hasParallelFault)); 1013 if (hasParallelFault) { 1014 DM dmHybrid = NULL, dmInterface; 1015 DMLabel faultLabel, faultBdLabel, hybridLabel; 1016 1017 PetscCall(DMGetLabel(*dm, "pfault", &faultLabel)); 1018 PetscCall(DMGetLabel(*dm, "pfaultBd", &faultBdLabel)); 1019 PetscCall(DMPlexCreateHybridMesh(*dm, faultLabel, faultBdLabel, 1, &hybridLabel, NULL, &dmInterface, &dmHybrid)); 1020 PetscCall(DMViewFromOptions(dmInterface, NULL, "-dm_fault_view")); 1021 { 1022 PetscViewer viewer; 1023 PetscMPIInt rank; 1024 1025 PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject) *dm), &rank)); 1026 PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &viewer)); 1027 PetscCall(PetscViewerASCIIPrintf(viewer, "Rank %d\n", rank)); 1028 PetscCall(DMLabelView(hybridLabel, viewer)); 1029 PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &viewer)); 1030 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 1031 } 1032 PetscCall(DMLabelDestroy(&hybridLabel)); 1033 PetscCall(DMDestroy(&dmInterface)); 1034 PetscCall(DMDestroy(dm)); 1035 *dm = dmHybrid; 1036 } 1037 PetscCall(PetscObjectSetName((PetscObject) *dm, "Hybrid Mesh")); 1038 PetscCall(CreateFaultLabel(*dm)); 1039 PetscCall(CreateDiscretization(*dm, user)); 1040 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view_pre")); 1041 PetscCall(DMPlexDistributeSetDefault(*dm, PETSC_FALSE)); 1042 PetscCall(DMSetFromOptions(*dm)); 1043 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 1044 PetscFunctionReturn(0); 1045 } 1046 1047 static PetscErrorCode TestMesh(DM dm, AppCtx *user) 1048 { 1049 PetscFunctionBegin; 1050 PetscCall(DMPlexCheckSymmetry(dm)); 1051 PetscCall(DMPlexCheckSkeleton(dm, 0)); 1052 PetscCall(DMPlexCheckFaces(dm, 0)); 1053 PetscFunctionReturn(0); 1054 } 1055 1056 static PetscErrorCode TestDiscretization(DM dm, AppCtx *user) 1057 { 1058 PetscSection s; 1059 1060 PetscFunctionBegin; 1061 PetscCall(DMGetSection(dm, &s)); 1062 PetscCall(PetscObjectViewFromOptions((PetscObject) s, NULL, "-local_section_view")); 1063 PetscFunctionReturn(0); 1064 } 1065 1066 static PetscErrorCode r(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx) 1067 { 1068 PetscInt d; 1069 for (d = 0; d < dim; ++d) u[d] = x[d]; 1070 return 0; 1071 } 1072 1073 static PetscErrorCode rp1(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx) 1074 { 1075 PetscInt d; 1076 for (d = 0; d < dim; ++d) u[d] = x[d] + (d > 0 ? 1.0 : 0.0); 1077 return 0; 1078 } 1079 1080 static PetscErrorCode phi(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx) 1081 { 1082 PetscInt d; 1083 u[0] = -x[1]; 1084 u[1] = x[0]; 1085 for (d = 2; d < dim; ++d) u[d] = x[d]; 1086 return 0; 1087 } 1088 1089 /* \lambda \cdot (\psi_u^- - \psi_u^+) */ 1090 static void f0_bd_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1091 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1092 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1093 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]) 1094 { 1095 const PetscInt Nc = uOff[1]-uOff[0]; 1096 PetscInt c; 1097 for (c = 0; c < Nc; ++c) { 1098 f0[c] = u[Nc*2+c] + x[Nc-c-1]; 1099 f0[Nc+c] = -(u[Nc*2+c] + x[Nc-c-1]); 1100 } 1101 } 1102 1103 /* (d - u^+ + u^-) \cdot \psi_\lambda */ 1104 static void f0_bd_l(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1105 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1106 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1107 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]) 1108 { 1109 const PetscInt Nc = uOff[2]-uOff[1]; 1110 PetscInt c; 1111 1112 for (c = 0; c < Nc; ++c) f0[c] = (c > 0 ? 1.0 : 0.0) + u[c] - u[Nc+c]; 1113 } 1114 1115 /* \psi_lambda \cdot (\psi_u^- - \psi_u^+) */ 1116 static void g0_bd_ul(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1117 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1118 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1119 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]) 1120 { 1121 const PetscInt Nc = uOff[1]-uOff[0]; 1122 PetscInt c; 1123 1124 for (c = 0; c < Nc; ++c) { 1125 g0[(0 +c)*Nc+c] = 1.0; 1126 g0[(Nc+c)*Nc+c] = -1.0; 1127 } 1128 } 1129 1130 /* (-\psi_u^+ + \psi_u^-) \cdot \psi_\lambda */ 1131 static void g0_bd_lu(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1132 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1133 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1134 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]) 1135 { 1136 const PetscInt Nc = uOff[2]-uOff[1]; 1137 PetscInt c; 1138 1139 for (c = 0; c < Nc; ++c) { 1140 g0[c*Nc*2+c] = 1.0; 1141 g0[c*Nc*2+Nc+c] = -1.0; 1142 } 1143 } 1144 1145 static PetscErrorCode TestAssembly(DM dm, AppCtx *user) 1146 { 1147 Mat J; 1148 Vec locX, locF; 1149 PetscDS probh; 1150 DMLabel fault, material; 1151 IS cohesiveCells; 1152 PetscWeakForm wf; 1153 PetscFormKey keys[3]; 1154 PetscErrorCode (*initialGuess[2])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx); 1155 PetscInt dim, Nf, cMax, cEnd, id; 1156 PetscMPIInt rank; 1157 1158 PetscFunctionBegin; 1159 PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject) dm), &rank)); 1160 PetscCall(DMGetDimension(dm, &dim)); 1161 PetscCall(DMPlexGetSimplexOrBoxCells(dm, 0, NULL, &cMax)); 1162 PetscCall(DMPlexGetHeightStratum(dm, 0, NULL, &cEnd)); 1163 PetscCall(ISCreateStride(PETSC_COMM_SELF, cEnd - cMax, cMax, 1, &cohesiveCells)); 1164 PetscCall(DMGetLabel(dm, "cohesive", &fault)); 1165 PetscCall(DMGetLocalVector(dm, &locX)); 1166 PetscCall(PetscObjectSetName((PetscObject) locX, "Local Solution")); 1167 PetscCall(DMGetLocalVector(dm, &locF)); 1168 PetscCall(PetscObjectSetName((PetscObject) locF, "Local Residual")); 1169 PetscCall(DMCreateMatrix(dm, &J)); 1170 PetscCall(PetscObjectSetName((PetscObject) J, "Jacobian")); 1171 1172 /* The initial guess has displacement shifted by one unit in each fault parallel direction across the fault */ 1173 PetscCall(DMGetLabel(dm, "material", &material)); 1174 id = 1; 1175 initialGuess[0] = r; 1176 initialGuess[1] = NULL; 1177 PetscCall(DMProjectFunctionLabelLocal(dm, 0.0, material, 1, &id, PETSC_DETERMINE, NULL, initialGuess, NULL, INSERT_VALUES, locX)); 1178 id = 2; 1179 initialGuess[0] = rp1; 1180 initialGuess[1] = NULL; 1181 PetscCall(DMProjectFunctionLabelLocal(dm, 0.0, material, 1, &id, PETSC_DETERMINE, NULL, initialGuess, NULL, INSERT_VALUES, locX)); 1182 id = 1; 1183 initialGuess[0] = NULL; 1184 initialGuess[1] = phi; 1185 PetscCall(DMProjectFunctionLabelLocal(dm, 0.0, fault, 1, &id, PETSC_DETERMINE, NULL, initialGuess, NULL, INSERT_VALUES, locX)); 1186 PetscCall(VecViewFromOptions(locX, NULL, "-local_solution_view")); 1187 1188 PetscCall(DMGetCellDS(dm, cMax, &probh)); 1189 PetscCall(PetscDSGetWeakForm(probh, &wf)); 1190 PetscCall(PetscDSGetNumFields(probh, &Nf)); 1191 PetscCall(PetscWeakFormSetIndexBdResidual(wf, material, 1, 0, 0, 0, f0_bd_u, 0, NULL)); 1192 PetscCall(PetscWeakFormSetIndexBdResidual(wf, material, 2, 0, 0, 0, f0_bd_u, 0, NULL)); 1193 PetscCall(PetscWeakFormSetIndexBdJacobian(wf, material, 1, 0, 1, 0, 0, g0_bd_ul, 0, NULL, 0, NULL, 0, NULL)); 1194 PetscCall(PetscWeakFormSetIndexBdJacobian(wf, material, 2, 0, 1, 0, 0, g0_bd_ul, 0, NULL, 0, NULL, 0, NULL)); 1195 if (Nf > 1) { 1196 PetscCall(PetscWeakFormSetIndexBdResidual(wf, fault, 1, 1, 0, 0, f0_bd_l, 0, NULL)); 1197 PetscCall(PetscWeakFormSetIndexBdJacobian(wf, fault, 1, 1, 0, 0, 0, g0_bd_lu, 0, NULL, 0, NULL, 0, NULL)); 1198 } 1199 if (!rank) PetscCall(PetscDSView(probh, NULL)); 1200 1201 keys[0].label = NULL; 1202 keys[0].value = 0; 1203 keys[0].field = 0; 1204 keys[0].part = 0; 1205 keys[1].label = material; 1206 keys[1].value = 2; 1207 keys[1].field = 0; 1208 keys[1].part = 0; 1209 keys[2].label = fault; 1210 keys[2].value = 1; 1211 keys[2].field = 1; 1212 keys[2].part = 0; 1213 PetscCall(VecSet(locF, 0.)); 1214 PetscCall(DMPlexComputeResidual_Hybrid_Internal(dm, keys, cohesiveCells, 0.0, locX, NULL, 0.0, locF, user)); 1215 PetscCall(VecViewFromOptions(locF, NULL, "-local_residual_view")); 1216 PetscCall(MatZeroEntries(J)); 1217 PetscCall(DMPlexComputeJacobian_Hybrid_Internal(dm, keys, cohesiveCells, 0.0, 0.0, locX, NULL, J, J, user)); 1218 PetscCall(MatViewFromOptions(J, NULL, "-local_jacobian_view")); 1219 1220 PetscCall(DMRestoreLocalVector(dm, &locX)); 1221 PetscCall(DMRestoreLocalVector(dm, &locF)); 1222 PetscCall(MatDestroy(&J)); 1223 PetscCall(ISDestroy(&cohesiveCells)); 1224 PetscFunctionReturn(0); 1225 } 1226 1227 int main(int argc, char **argv) 1228 { 1229 DM dm; 1230 AppCtx user; /* user-defined work context */ 1231 1232 PetscCall(PetscInitialize(&argc, &argv, NULL,help)); 1233 PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user)); 1234 PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm)); 1235 PetscCall(TestMesh(dm, &user)); 1236 PetscCall(TestDiscretization(dm, &user)); 1237 PetscCall(TestAssembly(dm, &user)); 1238 PetscCall(DMDestroy(&dm)); 1239 PetscCall(PetscFinalize()); 1240 return 0; 1241 } 1242 1243 /*TEST 1244 testset: 1245 args: -orig_dm_plex_check_all -dm_plex_check_all \ 1246 -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 -local_section_view \ 1247 -local_solution_view -local_residual_view -local_jacobian_view 1248 test: 1249 suffix: tri_0 1250 args: -dim 2 1251 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1252 test: 1253 suffix: tri_t1_0 1254 args: -dim 2 -test_num 1 1255 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1256 test: 1257 suffix: tri_t2_0 1258 args: -dim 2 -test_num 2 1259 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1260 test: 1261 suffix: tri_t5_0 1262 args: -dim 2 -test_num 5 1263 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1264 test: 1265 suffix: tet_0 1266 args: -dim 3 1267 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1268 test: 1269 suffix: tet_t1_0 1270 args: -dim 3 -test_num 1 1271 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1272 1273 testset: 1274 args: -orig_dm_plex_check_all -dm_plex_check_all \ 1275 -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 1276 test: 1277 suffix: tet_1 1278 nsize: 2 1279 args: -dim 3 1280 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1281 test: 1282 suffix: tri_1 1283 nsize: 2 1284 args: -dim 2 1285 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1286 test: 1287 suffix: tri_t3_0 1288 nsize: 2 1289 args: -dim 2 -test_num 3 1290 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1291 test: 1292 suffix: tri_t4_0 1293 nsize: 6 1294 args: -dim 2 -test_num 4 1295 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1296 1297 testset: 1298 args: -orig_dm_plex_check_all -dm_plex_check_all \ 1299 -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 1300 # 2D Quads 1301 test: 1302 suffix: quad_0 1303 args: -dim 2 -cell_simplex 0 1304 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1305 test: 1306 suffix: quad_1 1307 nsize: 2 1308 args: -dim 2 -cell_simplex 0 1309 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1310 test: 1311 suffix: quad_t1_0 1312 args: -dim 2 -cell_simplex 0 -test_num 1 -faulted_dm_plex_check_all 1313 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1314 test: 1315 suffix: quad_t2_0 1316 nsize: 2 1317 args: -dim 2 -cell_simplex 0 -test_num 2 1318 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1319 test: 1320 # TODO: The PetscSF is wrong here (connects to wrong side of split) 1321 suffix: quad_t3_0 1322 nsize: 2 1323 args: -dim 2 -cell_simplex 0 -test_num 3 1324 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1325 # 3D Hex 1326 test: 1327 suffix: hex_0 1328 args: -dim 3 -cell_simplex 0 1329 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1330 test: 1331 suffix: hex_1 1332 nsize: 2 1333 args: -dim 3 -cell_simplex 0 1334 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1335 test: 1336 suffix: hex_t1_0 1337 args: -dim 3 -cell_simplex 0 -test_num 1 1338 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1339 test: 1340 suffix: hex_t2_0 1341 args: -dim 3 -cell_simplex 0 -test_num 2 1342 filter: sed -e "s/_start//g" -e "s/f0_bd_u//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul//g" -e "s/g0_bd_lu//g" 1343 1344 TEST*/ 1345