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) PetscCall(DMPlexLabelComplete(*dm, matLabel)); 869 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 870 PetscCall(DMHasLabel(*dm, "fault", &hasFault)); 871 if (hasFault) { 872 DM dmHybrid = NULL, dmInterface = NULL; 873 DMLabel faultLabel, faultBdLabel, hybridLabel, splitLabel; 874 875 PetscCall(DMGetLabel(*dm, "fault", &faultLabel)); 876 PetscCall(DMGetLabel(*dm, "faultBd", &faultBdLabel)); 877 PetscCall(DMPlexCreateHybridMesh(*dm, faultLabel, faultBdLabel, 1, &hybridLabel, &splitLabel, &dmInterface, &dmHybrid)); 878 PetscCall(DMLabelView(hybridLabel, PETSC_VIEWER_STDOUT_WORLD)); 879 PetscCall(DMLabelDestroy(&hybridLabel)); 880 PetscCall(DMLabelView(splitLabel, PETSC_VIEWER_STDOUT_WORLD)); 881 PetscCall(DMLabelDestroy(&splitLabel)); 882 PetscCall(DMViewFromOptions(dmInterface, NULL, "-dm_interface_view")); 883 PetscCall(DMDestroy(&dmInterface)); 884 PetscCall(DMDestroy(dm)); 885 *dm = dmHybrid; 886 } 887 PetscCall(DMHasLabel(*dm, "fault2", &hasFault2)); 888 if (hasFault2) { 889 DM dmHybrid = NULL; 890 DMLabel faultLabel, faultBdLabel, hybridLabel; 891 892 PetscCall(PetscObjectSetOptionsPrefix((PetscObject) *dm, "faulted_")); 893 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view_pre")); 894 PetscCall(DMPlexDistributeSetDefault(*dm, PETSC_FALSE)); 895 PetscCall(DMSetFromOptions(*dm)); 896 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 897 PetscCall(DMGetLabel(*dm, "fault2", &faultLabel)); 898 PetscCall(DMGetLabel(*dm, "fault2Bd", &faultBdLabel)); 899 PetscCall(DMPlexCreateHybridMesh(*dm, faultLabel, faultBdLabel, 1, &hybridLabel, NULL, NULL, &dmHybrid)); 900 PetscCall(DMLabelView(hybridLabel, PETSC_VIEWER_STDOUT_WORLD)); 901 PetscCall(DMLabelDestroy(&hybridLabel)); 902 PetscCall(DMDestroy(dm)); 903 *dm = dmHybrid; 904 } 905 if (user->testPartition && size > 1) { 906 PetscPartitioner part; 907 PetscInt *sizes = NULL; 908 PetscInt *points = NULL; 909 910 if (rank == 0) { 911 if (dim == 2 && cellSimplex && size == 2) { 912 switch (user->testNum) { 913 case 0: { 914 PetscInt triSizes_p2[2] = {1, 2}; 915 PetscInt triPoints_p2[3] = {0, 1, 2}; 916 917 PetscCall(PetscMalloc2(2, &sizes, 3, &points)); 918 PetscCall(PetscArraycpy(sizes, triSizes_p2, 2)); 919 PetscCall(PetscArraycpy(points, triPoints_p2, 3));break;} 920 case 3: { 921 PetscInt triSizes_p2[2] = {3, 3}; 922 PetscInt triPoints_p2[6] = {0, 1, 2, 3, 4, 5}; 923 924 PetscCall(PetscMalloc2(2, &sizes, 6, &points)); 925 PetscCall(PetscArraycpy(sizes, triSizes_p2, 2)); 926 PetscCall(PetscArraycpy(points, triPoints_p2, 6));break;} 927 default: 928 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for triangular mesh on 2 procs", user->testNum); 929 } 930 } else if (dim == 2 && cellSimplex && size == 6) { 931 switch (user->testNum) { 932 case 4: { 933 PetscInt triSizes_p6[6] = {1, 1, 1, 1, 1, 1}; 934 PetscInt triPoints_p6[6] = {0, 1, 2, 3, 4, 5}; 935 936 PetscCall(PetscMalloc2(6, &sizes, 6, &points)); 937 PetscCall(PetscArraycpy(sizes, triSizes_p6, 6)); 938 PetscCall(PetscArraycpy(points, triPoints_p6, 6));break;} 939 default: 940 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for triangular mesh on 6 procs", user->testNum); 941 } 942 } else if (dim == 2 && !cellSimplex && size == 2) { 943 switch (user->testNum) { 944 case 0: { 945 PetscInt quadSizes_p2[2] = {1, 2}; 946 PetscInt quadPoints_p2[3] = {0, 1, 2}; 947 948 PetscCall(PetscMalloc2(2, &sizes, 3, &points)); 949 PetscCall(PetscArraycpy(sizes, quadSizes_p2, 2)); 950 PetscCall(PetscArraycpy(points, quadPoints_p2, 3));break;} 951 case 2: { 952 PetscInt quadSizes_p2[2] = {1, 1}; 953 PetscInt quadPoints_p2[2] = {0, 1}; 954 955 PetscCall(PetscMalloc2(2, &sizes, 2, &points)); 956 PetscCall(PetscArraycpy(sizes, quadSizes_p2, 2)); 957 PetscCall(PetscArraycpy(points, quadPoints_p2, 2));break;} 958 case 3: { 959 PetscInt quadSizes_p2[2] = {1, 1}; 960 PetscInt quadPoints_p2[2] = {1, 0}; 961 962 PetscCall(PetscMalloc2(2, &sizes, 2, &points)); 963 PetscCall(PetscArraycpy(sizes, quadSizes_p2, 2)); 964 PetscCall(PetscArraycpy(points, quadPoints_p2, 2));break;} 965 default: 966 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for quadrilateral mesh on 2 procs", user->testNum); 967 } 968 } else if (dim == 3 && cellSimplex && size == 2) { 969 switch (user->testNum) { 970 case 0: { 971 PetscInt tetSizes_p2[2] = {1, 2}; 972 PetscInt tetPoints_p2[3] = {0, 1, 2}; 973 974 PetscCall(PetscMalloc2(2, &sizes, 3, &points)); 975 PetscCall(PetscArraycpy(sizes, tetSizes_p2, 2)); 976 PetscCall(PetscArraycpy(points, tetPoints_p2, 3));break;} 977 default: 978 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for teterehedral mesh on 2 procs", user->testNum); 979 } 980 } else if (dim == 3 && !cellSimplex && size == 2) { 981 switch (user->testNum) { 982 case 0: { 983 PetscInt hexSizes_p2[2] = {1, 2}; 984 PetscInt hexPoints_p2[3] = {0, 1, 2}; 985 986 PetscCall(PetscMalloc2(2, &sizes, 3, &points)); 987 PetscCall(PetscArraycpy(sizes, hexSizes_p2, 2)); 988 PetscCall(PetscArraycpy(points, hexPoints_p2, 3));break;} 989 default: 990 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for hexahedral mesh on 2 procs", user->testNum); 991 } 992 } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test partition"); 993 } 994 PetscCall(DMPlexGetPartitioner(*dm, &part)); 995 PetscCall(PetscPartitionerSetType(part, PETSCPARTITIONERSHELL)); 996 PetscCall(PetscPartitionerShellSetPartition(part, size, sizes, points)); 997 PetscCall(PetscFree2(sizes, points)); 998 } 999 { 1000 DM pdm = NULL; 1001 1002 /* Distribute mesh over processes */ 1003 PetscCall(DMPlexDistribute(*dm, 0, NULL, &pdm)); 1004 if (pdm) { 1005 PetscCall(DMViewFromOptions(pdm, NULL, "-dm_view")); 1006 PetscCall(DMDestroy(dm)); 1007 *dm = pdm; 1008 } 1009 } 1010 PetscCall(DMHasLabel(*dm, "pfault", &hasParallelFault)); 1011 if (hasParallelFault) { 1012 DM dmHybrid = NULL, dmInterface; 1013 DMLabel faultLabel, faultBdLabel, hybridLabel; 1014 1015 PetscCall(DMGetLabel(*dm, "pfault", &faultLabel)); 1016 PetscCall(DMGetLabel(*dm, "pfaultBd", &faultBdLabel)); 1017 PetscCall(DMPlexCreateHybridMesh(*dm, faultLabel, faultBdLabel, 1, &hybridLabel, NULL, &dmInterface, &dmHybrid)); 1018 PetscCall(DMViewFromOptions(dmInterface, NULL, "-dm_fault_view")); 1019 { 1020 PetscViewer viewer; 1021 PetscMPIInt rank; 1022 1023 PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject) *dm), &rank)); 1024 PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &viewer)); 1025 PetscCall(PetscViewerASCIIPrintf(viewer, "Rank %d\n", rank)); 1026 PetscCall(DMLabelView(hybridLabel, viewer)); 1027 PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &viewer)); 1028 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 1029 } 1030 PetscCall(DMLabelDestroy(&hybridLabel)); 1031 PetscCall(DMDestroy(&dmInterface)); 1032 PetscCall(DMDestroy(dm)); 1033 *dm = dmHybrid; 1034 } 1035 PetscCall(PetscObjectSetName((PetscObject) *dm, "Hybrid Mesh")); 1036 PetscCall(CreateFaultLabel(*dm)); 1037 PetscCall(CreateDiscretization(*dm, user)); 1038 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view_pre")); 1039 PetscCall(DMPlexDistributeSetDefault(*dm, PETSC_FALSE)); 1040 PetscCall(DMSetFromOptions(*dm)); 1041 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 1042 PetscFunctionReturn(0); 1043 } 1044 1045 static PetscErrorCode TestMesh(DM dm, AppCtx *user) 1046 { 1047 PetscFunctionBegin; 1048 PetscCall(DMPlexCheckSymmetry(dm)); 1049 PetscCall(DMPlexCheckSkeleton(dm, 0)); 1050 PetscCall(DMPlexCheckFaces(dm, 0)); 1051 PetscFunctionReturn(0); 1052 } 1053 1054 static PetscErrorCode TestDiscretization(DM dm, AppCtx *user) 1055 { 1056 PetscSection s; 1057 1058 PetscFunctionBegin; 1059 PetscCall(DMGetSection(dm, &s)); 1060 PetscCall(PetscObjectViewFromOptions((PetscObject) s, NULL, "-local_section_view")); 1061 PetscFunctionReturn(0); 1062 } 1063 1064 static PetscErrorCode r(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx) 1065 { 1066 PetscInt d; 1067 for (d = 0; d < dim; ++d) u[d] = x[d]; 1068 return 0; 1069 } 1070 1071 static PetscErrorCode rp1(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx) 1072 { 1073 PetscInt d; 1074 for (d = 0; d < dim; ++d) u[d] = x[d] + (d > 0 ? 1.0 : 0.0); 1075 return 0; 1076 } 1077 1078 static PetscErrorCode phi(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx) 1079 { 1080 PetscInt d; 1081 u[0] = -x[1]; 1082 u[1] = x[0]; 1083 for (d = 2; d < dim; ++d) u[d] = x[d]; 1084 return 0; 1085 } 1086 1087 /* \lambda \cdot (\psi_u^- - \psi_u^+) */ 1088 static void f0_bd_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1089 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1090 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1091 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]) 1092 { 1093 const PetscInt Nc = uOff[1]-uOff[0]; 1094 PetscInt c; 1095 for (c = 0; c < Nc; ++c) { 1096 f0[c] = u[Nc*2+c] + x[Nc-c-1]; 1097 f0[Nc+c] = -(u[Nc*2+c] + x[Nc-c-1]); 1098 } 1099 } 1100 1101 /* (d - u^+ + u^-) \cdot \psi_\lambda */ 1102 static void f0_bd_l(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1103 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1104 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1105 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]) 1106 { 1107 const PetscInt Nc = uOff[2]-uOff[1]; 1108 PetscInt c; 1109 1110 for (c = 0; c < Nc; ++c) f0[c] = (c > 0 ? 1.0 : 0.0) + u[c] - u[Nc+c]; 1111 } 1112 1113 /* \psi_lambda \cdot (\psi_u^- - \psi_u^+) */ 1114 static void g0_bd_ul(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1115 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1116 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1117 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]) 1118 { 1119 const PetscInt Nc = uOff[1]-uOff[0]; 1120 PetscInt c; 1121 1122 for (c = 0; c < Nc; ++c) { 1123 g0[(0 +c)*Nc+c] = 1.0; 1124 g0[(Nc+c)*Nc+c] = -1.0; 1125 } 1126 } 1127 1128 /* (-\psi_u^+ + \psi_u^-) \cdot \psi_\lambda */ 1129 static void g0_bd_lu(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1130 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1131 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1132 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]) 1133 { 1134 const PetscInt Nc = uOff[2]-uOff[1]; 1135 PetscInt c; 1136 1137 for (c = 0; c < Nc; ++c) { 1138 g0[c*Nc*2+c] = 1.0; 1139 g0[c*Nc*2+Nc+c] = -1.0; 1140 } 1141 } 1142 1143 static PetscErrorCode TestAssembly(DM dm, AppCtx *user) 1144 { 1145 Mat J; 1146 Vec locX, locF; 1147 PetscDS probh; 1148 DMLabel fault, material; 1149 IS cohesiveCells; 1150 PetscWeakForm wf; 1151 PetscFormKey keys[3]; 1152 PetscErrorCode (*initialGuess[2])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx); 1153 PetscInt dim, Nf, cMax, cEnd, id; 1154 PetscMPIInt rank; 1155 1156 PetscFunctionBegin; 1157 PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject) dm), &rank)); 1158 PetscCall(DMGetDimension(dm, &dim)); 1159 PetscCall(DMPlexGetSimplexOrBoxCells(dm, 0, NULL, &cMax)); 1160 PetscCall(DMPlexGetHeightStratum(dm, 0, NULL, &cEnd)); 1161 PetscCall(ISCreateStride(PETSC_COMM_SELF, cEnd - cMax, cMax, 1, &cohesiveCells)); 1162 PetscCall(DMGetLabel(dm, "cohesive", &fault)); 1163 PetscCall(DMGetLocalVector(dm, &locX)); 1164 PetscCall(PetscObjectSetName((PetscObject) locX, "Local Solution")); 1165 PetscCall(DMGetLocalVector(dm, &locF)); 1166 PetscCall(PetscObjectSetName((PetscObject) locF, "Local Residual")); 1167 PetscCall(DMCreateMatrix(dm, &J)); 1168 PetscCall(PetscObjectSetName((PetscObject) J, "Jacobian")); 1169 1170 /* The initial guess has displacement shifted by one unit in each fault parallel direction across the fault */ 1171 PetscCall(DMGetLabel(dm, "material", &material)); 1172 id = 1; 1173 initialGuess[0] = r; 1174 initialGuess[1] = NULL; 1175 PetscCall(DMProjectFunctionLabelLocal(dm, 0.0, material, 1, &id, PETSC_DETERMINE, NULL, initialGuess, NULL, INSERT_VALUES, locX)); 1176 id = 2; 1177 initialGuess[0] = rp1; 1178 initialGuess[1] = NULL; 1179 PetscCall(DMProjectFunctionLabelLocal(dm, 0.0, material, 1, &id, PETSC_DETERMINE, NULL, initialGuess, NULL, INSERT_VALUES, locX)); 1180 id = 1; 1181 initialGuess[0] = NULL; 1182 initialGuess[1] = phi; 1183 PetscCall(DMProjectFunctionLabelLocal(dm, 0.0, fault, 1, &id, PETSC_DETERMINE, NULL, initialGuess, NULL, INSERT_VALUES, locX)); 1184 PetscCall(VecViewFromOptions(locX, NULL, "-local_solution_view")); 1185 1186 PetscCall(DMGetCellDS(dm, cMax, &probh)); 1187 PetscCall(PetscDSGetWeakForm(probh, &wf)); 1188 PetscCall(PetscDSGetNumFields(probh, &Nf)); 1189 PetscCall(PetscWeakFormSetIndexBdResidual(wf, material, 1, 0, 0, 0, f0_bd_u, 0, NULL)); 1190 PetscCall(PetscWeakFormSetIndexBdResidual(wf, material, 2, 0, 0, 0, f0_bd_u, 0, NULL)); 1191 PetscCall(PetscWeakFormSetIndexBdJacobian(wf, material, 1, 0, 1, 0, 0, g0_bd_ul, 0, NULL, 0, NULL, 0, NULL)); 1192 PetscCall(PetscWeakFormSetIndexBdJacobian(wf, material, 2, 0, 1, 0, 0, g0_bd_ul, 0, NULL, 0, NULL, 0, NULL)); 1193 if (Nf > 1) { 1194 PetscCall(PetscWeakFormSetIndexBdResidual(wf, fault, 1, 1, 0, 0, f0_bd_l, 0, NULL)); 1195 PetscCall(PetscWeakFormSetIndexBdJacobian(wf, fault, 1, 1, 0, 0, 0, g0_bd_lu, 0, NULL, 0, NULL, 0, NULL)); 1196 } 1197 if (rank == 0) PetscCall(PetscDSView(probh, NULL)); 1198 1199 keys[0].label = NULL; 1200 keys[0].value = 0; 1201 keys[0].field = 0; 1202 keys[0].part = 0; 1203 keys[1].label = material; 1204 keys[1].value = 2; 1205 keys[1].field = 0; 1206 keys[1].part = 0; 1207 keys[2].label = fault; 1208 keys[2].value = 1; 1209 keys[2].field = 1; 1210 keys[2].part = 0; 1211 PetscCall(VecSet(locF, 0.)); 1212 PetscCall(DMPlexComputeResidual_Hybrid_Internal(dm, keys, cohesiveCells, 0.0, locX, NULL, 0.0, locF, user)); 1213 PetscCall(VecViewFromOptions(locF, NULL, "-local_residual_view")); 1214 PetscCall(MatZeroEntries(J)); 1215 PetscCall(DMPlexComputeJacobian_Hybrid_Internal(dm, keys, cohesiveCells, 0.0, 0.0, locX, NULL, J, J, user)); 1216 PetscCall(MatAssemblyBegin(J, MAT_FINAL_ASSEMBLY)); 1217 PetscCall(MatAssemblyEnd(J, MAT_FINAL_ASSEMBLY)); 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 PetscFunctionBeginUser; 1233 PetscCall(PetscInitialize(&argc, &argv, NULL,help)); 1234 PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user)); 1235 PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm)); 1236 PetscCall(TestMesh(dm, &user)); 1237 PetscCall(TestDiscretization(dm, &user)); 1238 PetscCall(TestAssembly(dm, &user)); 1239 PetscCall(DMDestroy(&dm)); 1240 PetscCall(PetscFinalize()); 1241 return 0; 1242 } 1243 1244 /*TEST 1245 testset: 1246 args: -orig_dm_plex_check_all -dm_plex_check_all \ 1247 -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 -local_section_view \ 1248 -local_solution_view -local_residual_view -local_jacobian_view 1249 test: 1250 suffix: tri_0 1251 args: -dim 2 1252 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" 1253 test: 1254 suffix: tri_t1_0 1255 args: -dim 2 -test_num 1 1256 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" 1257 test: 1258 suffix: tri_t2_0 1259 args: -dim 2 -test_num 2 1260 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" 1261 test: 1262 suffix: tri_t5_0 1263 args: -dim 2 -test_num 5 1264 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" 1265 test: 1266 suffix: tet_0 1267 args: -dim 3 1268 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" 1269 test: 1270 suffix: tet_t1_0 1271 args: -dim 3 -test_num 1 1272 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" 1273 1274 testset: 1275 args: -orig_dm_plex_check_all -dm_plex_check_all \ 1276 -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 1277 test: 1278 suffix: tet_1 1279 nsize: 2 1280 args: -dim 3 1281 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" 1282 test: 1283 suffix: tri_1 1284 nsize: 2 1285 args: -dim 2 1286 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" 1287 test: 1288 suffix: tri_t3_0 1289 nsize: 2 1290 args: -dim 2 -test_num 3 1291 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" 1292 test: 1293 suffix: tri_t4_0 1294 nsize: 6 1295 args: -dim 2 -test_num 4 1296 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" 1297 1298 testset: 1299 args: -orig_dm_plex_check_all -dm_plex_check_all \ 1300 -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 1301 # 2D Quads 1302 test: 1303 suffix: quad_0 1304 args: -dim 2 -cell_simplex 0 1305 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" 1306 test: 1307 suffix: quad_1 1308 nsize: 2 1309 args: -dim 2 -cell_simplex 0 1310 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" 1311 test: 1312 suffix: quad_t1_0 1313 args: -dim 2 -cell_simplex 0 -test_num 1 -faulted_dm_plex_check_all 1314 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" 1315 test: 1316 suffix: quad_t2_0 1317 nsize: 2 1318 args: -dim 2 -cell_simplex 0 -test_num 2 1319 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" 1320 test: 1321 # TODO: The PetscSF is wrong here (connects to wrong side of split) 1322 suffix: quad_t3_0 1323 nsize: 2 1324 args: -dim 2 -cell_simplex 0 -test_num 3 1325 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" 1326 # 3D Hex 1327 test: 1328 suffix: hex_0 1329 args: -dim 3 -cell_simplex 0 1330 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" 1331 test: 1332 suffix: hex_1 1333 nsize: 2 1334 args: -dim 3 -cell_simplex 0 1335 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" 1336 test: 1337 suffix: hex_t1_0 1338 args: -dim 3 -cell_simplex 0 -test_num 1 1339 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" 1340 test: 1341 suffix: hex_t2_0 1342 args: -dim 3 -cell_simplex 0 -test_num 2 1343 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" 1344 1345 TEST*/ 1346