Actual source code: ex69.c

  1: static char help[] = "Tests for creation of cohesive meshes by transforms\n\n";

  3: #include <petscdmplex.h>
  4: #include <petscsf.h>

  6: #include <petsc/private/dmpleximpl.h>

  8: PETSC_EXTERN char tri_2_cv[];
  9: char              tri_2_cv[] = "\
 10: 2 4 6 3 1\n\
 11: 0 2 1\n\
 12: 1 2 3\n\
 13: 4 1 5\n\
 14: 4 0 1\n\
 15: -1.0  0.0 0.0  1\n\
 16:  0.0  1.0 0.0 -1\n\
 17:  0.0 -1.0 0.0  1\n\
 18:  1.0  0.0 0.0 -1\n\
 19: -2.0  1.0 0.0  1\n\
 20: -1.0  2.0 0.0 -1";

 22: /* List of test meshes

 24: Test tri_0: triangle

 26:  4-10--5      8-16--7-14--4
 27:  |\  1 |      |\     \  1 |
 28:  | \   |      | \     \   |
 29:  6  8  9  ->  9 12  2  11 13
 30:  |   \ |      |   \     \ |
 31:  | 0  \|      | 0  \     \|
 32:  2--7--3      3-10--6-15--5

 34: Test tri_1: triangle, not tensor

 36:  4-10--5      8-10--7-16--4
 37:  |\  1 |      |\     \  1 |
 38:  | \   |      | \     \   |
 39:  6  8  9  -> 11 14  2  13 15
 40:  |   \ |      |   \     \ |
 41:  | 0  \|      | 0  \     \|
 42:  2--7--3      3-12--6--9--5

 44: Test tri_2: 4 triangles, non-oriented surface

 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:   becomes
 66:            8
 67:           / \
 68:          /   \
 69:         /     \
 70:       25   2   24
 71:       /         \
 72:      /           \
 73:    13-----18------9
 74: 28  |     5    26/ \
 75:    14----19----10   \
 76:      \         /|   |\
 77:       \       / |   | \
 78:       21  3  20 |   |  23
 79:         \   /   |   |   \
 80:          \ /    |   |    \
 81:           6  0 17 4 16 1  7
 82:            \    |   |    /
 83:             \   |   |   /
 84:             15  |   |  22
 85:               \ |   | /
 86:                \|   |/
 87:                12---11
 88:                  27

 90: Test tri_3: tri_2, in parallel

 92:            6
 93:           / \
 94:          /   \
 95:         /     \
 96:       12   1   11
 97:       /         \
 98:      /           \
 99:     5-----10------2
100:                    \
101:     5-----9-----3   2
102:      \         /|   |\
103:       \       / |   | \
104:       10  1  8  |   |  9
105:         \   /   |   |   \
106:          \ /    |   |    \
107:           2  0  7   7  0  4
108:            \    |   |    /
109:             \   |   |   /
110:              6  |   |  8
111:               \ |   | /
112:                \|   |/
113:                 4   3
114:   becomes
115:                  11
116:                 / \
117:                /   \
118:               /     \
119:             19   1   18
120:             /         \
121:            /           \
122:           8-----14------4
123:         22 \     3       |
124:             9------15    |\
125:                     \    | \
126:     9------14-----5  \  20 |
127:   20\    3     18/ \  \/   |
128:    10----15-----6   |  5   |
129:      \         /|   |  |   |\
130:       \       / |   |  |   | \
131:       17  1 16  |   |  |   |  17
132:         \   /   | 2 |  | 2 |   \
133:          \ /    |   |  |   |    \
134:           4  0  13 12  13  12 0 10
135:            \    |   |  |   |    /
136:             \   |   |  |   |   /
137:             11  |   |  |   |  16
138:               \ |   |  |   | /
139:                \|   |  |   |/
140:                 8---7  7---6
141:                  19      21

143: Test quad_0: quadrilateral

145:  5-10--6-11--7       5-12-10-20--9-14--6
146:  |     |     |       |     |     |     |
147: 12  0 13  1  14 --> 15  0 18  2 17  1  16
148:  |     |     |       |     |     |     |
149:  2--8--3--9--4       3-11--8-19--7-13--4

151: Test quad_1: quadrilateral, not tensor

153:  5-10--6-11--7       5-14-10-12--9-16--6
154:  |     |     |       |     |     |     |
155: 12  0 13  1  14 --> 17  0 20  2 19  1  18
156:  |     |     |       |     |     |     |
157:  2--8--3--9--4       3-13--8-11--7-15--4

159: Test quad_2: quadrilateral, 2 processes

161:  3--6--4  3--6--4       3--9--7-14--6   5-14--4--9--7
162:  |     |  |     |       |     |     |   |     |     |
163:  7  0  8  7  0  8  --> 10  0 12  1 11  12  1 11  0  10
164:  |     |  |     |       |     |     |   |     |     |
165:  1--5--2  1--5--2       2--8--5-13--4   3-13--2--8--6

167: Test quad_3: quadrilateral, 4 processes, non-oriented surface

169:  3--6--4  3--6--4      3--9--7-14--6   5-14--4--9--7
170:  |     |  |     |      |     |     |   |     |     |
171:  7  0  8  7  0  8     10  0  12 1  11 12  1 11  0  10
172:  |     |  |     |      |     |     |   |     |     |
173:  1--5--2  1--5--2      2--8--5-13--4   3-13--2--8--6
174:                    -->
175:  3--6--4  3--6--4      3--9--7-14--6   5-14--4--9--7
176:  |     |  |     |      |     |     |   |     |     |
177:  7  0  8  7  0  8     10  0  12 1  11 12  1 11  0  10
178:  |     |  |     |      |     |     |   |     |     |
179:  1--5--2  1--5--2      2--8--5-13--4   3-13--2--8--6

181: Test quad_4: embedded fault

183: 14-24-15-25-16-26--17
184:  |     |     |     |
185: 28  3 30  4 32  5  34
186:  |     |     |     |
187: 10-21-11-22-12-23--13
188:  |     |     |     |
189: 27  0 29  1 31  2  33
190:  |     |     |     |
191:  6-18--7-19--8-20--9

193: becomes

195:  13-26-14-27-15-28--16
196:   |     |     |     |
197:  30  3 32  4 39  5  40
198:   |     |     |     |
199:  12-25-17-36-19-38--21
200:         |     |     |
201:        41  6 42  7  43
202:         |     |     |
203:  12-25-17-35-18-37--20
204:   |     |     |     |
205:  29  0 31  1 33  2  34
206:   |     |     |     |
207:   8-22--9-23-10-24--11

209: Test quad_5: two faults

211: 14-24-15-25-16-26--17
212:  |     |     |     |
213: 28  3 30  4 32  5  34
214:  |     |     |     |
215: 10-21-11-22-12-23--13
216:  |     |     |     |
217: 27  0 29  1 31  2  33
218:  |     |     |     |
219:  6-18--7-19--8-20--9

221: becomes

223: 12-26-13-27-14-28--15
224:  |     |     |     |
225: 37  4 31  3 33  5  40
226:  |     |     |     |
227: 17-36-18-25-19-39--21
228:  |     |     |     |
229: 43  6  44   41  7  42
230:  |     |     |     |
231: 16-35-18-25-19-38--20
232:  |     |     |     |
233: 29  0 30  1 32  2  34
234:  |     |     |     |
235:  8-22--9-23-10-24--11

237: Test quad_6: T-junction

239: 14-24-15-25-16-26--17
240:  |     |     |     |
241: 28  3 30  4 32  5  34
242:  |     |     |     |
243: 10-21-11-22-12-23--13
244:  |     |     |     |
245: 27  0 29  1 31  2  33
246:  |     |     |     |
247:  6-18--7-19--8-20--9

249: becomes

251:  13-26-14-27-15-28--16
252:   |     |     |     |
253:  30  3 32  4 39  5  40
254:   |     |     |     |
255:  12-25-17-36-19-38--21
256:         |     |     |
257:        41  6 42  7  43
258:         |     |     |
259:  12-25-17-35-18-37--20
260:   |     |     |     |
261:  29  0 31  1 33  2  34
262:   |     |     |     |
263:   8-22--9-23-10-24--11

265: becomes

267:  14-28-15-41-21-44--20-29-16
268:   |     |     |     |     |
269:  31  3 33  5 43  8 42  4  40
270:   |     |     |     |     |
271:  13-27-17-37-23-46--23-39-19
272:         |     |     |     |
273:        47  6 48    48  7  49
274:         |     |     |     |
275:  13-27-17-36-22-45--22-38-18
276:   |     |     |     |     |
277:  30  0 32  1 34    34  2  35
278:   |     |     |     |     |
279:   9-24-10-25-11-----11-26-12

281: Test tet_0: Two tets sharing a face

283:  cell   5 _______    cell
284:  0    / | \      \      1
285:     19  |  16     20
286:     /  15   \      \
287:    2-17------4--22--6
288:     \   |   /      /
289:     18  |  14     21
290:       \ | /      /
291:         3-------

293: becomes

295:  cell  10 ___36____9______    cell
296:  0    / | \        |\      \     1
297:     29  |  27      | 26     31
298:     /  25   \     24  \      \
299:    3-28------8--35-----7--33--4
300:     \   |   /      |  /      /
301:     30  |  23      | 22     32
302:       \ | /        |/      /
303:         6----34----5------
304:          cell 2

306: Test tet_1: Two tets sharing a face in parallel

308:  cell   4          3______    cell
309:  0    / | \        |\      \     0
310:     14  |  11      | 11     12
311:     /  10   \     10  \      \
312:    1-12------3     |   2--14--4
313:     \   |   /      |  /      /
314:     13  |  9       | 9      13
315:       \ | /        |/      /
316:         2          1------

318: becomes
319:            cell 1              cell 1
320:  cell   8---28---7           7---28---6______    cell
321:  0    / | \      |\          |\       |\      \     0
322:     24  |  22    | 21        | 22     | 21     23
323:     /  20   \    |   \       |  \    19  \      \
324:    2-23------6---27---5     20  5---27---4--25--8
325:     \   |   /   19   /       |  /     |  /      /
326:     25  |  18    | 17        | 18     | 17     24
327:       \ | /      |/          |/       |/      /
328:         4---26---3           3---26---2------

330: Test hex_0: Two hexes sharing a face

332: cell  11-----31-----12-----32------13 cell
333: 0     /|            /|            /|     1
334:     36 |   22      37|   24      38|
335:     /  |          /  |          /  |
336:    8-----29------9-----30------10  |
337:    |   |     18  |   |     20  |   |
338:    |  42         |  43         |   44
339:    |14 |         |15 |         |16 |
340:   39   |  17    40   |   19   41   |
341:    |   5-----27--|---6-----28--|---7
342:    |  /          |  /          |  /
343:    | 33   21     | 34    23    | 35
344:    |/            |/            |/
345:    2-----25------3-----26------4

347: becomes

349:                          cell 2
350: cell   9-----38-----18-----62------17----42------10 cell
351: 0     /|            /|            /|            /|     1
352:     45 |   30      54|  32       53|   24      46|
353:     /  |          /  |          /  |          /  |
354:    7-----37-----16-----61------15--|-41------8   |
355:    |   |     28  |   |         |   |     22  |   |
356:    |  49         |  58         |   57        |   50
357:    |19 |         |26 |         |25 |         |20 |
358:   47   |  27    56   |        55   |   21   48   |
359:    |   5-----36--|--14-----60--|---13----40--|---6
360:    |  /          |  /          |  /          |  /
361:    | 43   29     | 52   31     | 51    23    | 44
362:    |/            |/            |/            |/
363:    3-----35-----12-----59------11----39------4

365: Test hex_1: Two hexes sharing a face, in parallel

367: cell   7-----18------8             7-----18------8 cell
368: 0     /|            /|            /|            /|    0
369:     21 |   14      22|           21|   14      22|
370:     /  |          /  |          /  |          /  |
371:    5-----17------6   |         5---|-17------6   |
372:    |   |     12  |   |         |   |     12  |   |
373:    |  25         |  26         |  25         |  26
374:    | 9 |         |10 |         | 9 |         |10 |
375:   23   |  11    24   |        23   |   11   24   |
376:    |   3-----16--|---4         |   3-----16--|---4
377:    |  /          |  /          |  /          |  /
378:    | 19   13     | 20          | 19    13    | 20
379:    |/            |/            |/            |/
380:    1-----15------2             1-----15------2

382: becomes
383:                         cell 1                      cell 1
384: cell   5-----28-----13-----44-----12             9-----44-----8-----28------13 cell
385: 0     /|            /|           /|             /|           /|            /|     0
386:     30 |   20      36|   22     35|            36|   22     35|   20      30|
387:     /  |          /  |         /  |           /  |         /  |          /  |
388:    4-----27-----11-----43-----10  |          7-----43-----6-----27------12  |
389:    |   |     18  |   |        |   |          |   |        |   |     18  |   |
390:    |  32         |  40        |   39         |  40        |   39        |   32
391:    |14 |         |16 |        | 15|          |15 |        |14 |         |16 |
392:   31   |  17    38   |        37  |         38   |       37   |   17   31   |
393:    |   3-----26--|---9-----42-|---8          |   5----42--|---4-----26--|---11
394:    |  /          |  /         |  /           |  /         |  /          |  /
395:    | 29   19     | 34    21   | 33           | 34    21   | 33    19    | 29
396:    |/            |/           |/             |/           |/            |/
397:    2-----25------7-----41-----6              3-----41-----2-----25------10

399: Test hex_2: hexahedra, 4 processes, non-oriented surface

401:           cell 0                  cell 0
402:        7-----18------8       7-----18------8
403:       /|            /|      /|            /|
404:     21 |   14      22|    21 |   14      22|
405:     /  |          /  |    /  |          /  |
406:    5-----17------6   |   5-----17------6   |
407:    |   |     12  |   |   |   |     12  |   |
408:    |  25         |  26   |  25         |   26
409:    |9  |         |10 |   |9  |         |10 |
410:   23   |  11    24   |  23   |  11    24   |
411:    |   3-----16--|---4   |   3-----16--|---4
412:    |  /          |  /    |  /          |  /
413:    | 19    13    | 20    | 19    13    | 20
414:    |/            |/      |/            |/
415:    1-----15------2       1-----15------2

417:        7-----18------8       7-----18------8
418:       /|            /|      /|            /|
419:     21 |   14      22|    21 |   14      22|
420:     /  |          /  |    /  |          /  |
421:    5-----17------6   |   5-----17------6   |
422:    |   |     12  |   |   |   |     12  |   |
423:    |  25         |  26   |  25         |  26
424:    |9  |         |10 |   |9  |         |10 |
425:   23   |  11    24   |  23   |   11   24   |
426:    |   3-----16--|---4   |   3-----16--|---4
427:    |  /          |  /    |  /          |  /
428:    | 19   13     | 20    | 19    13    | 20
429:    |/            |/      |/            |/
430:    1-----15------2       1-----15------2
431:       cell 0                cell 0

433: becomes

435:           cell 0         cell 1                cell 1        cell 0
436:        5-----28------13----44------12      9-----44------8-----28------13
437:       /|            /|            /|      /|            /|            /|
438:     30 |   20      36|   22      35|     36|   22     35 |   20      30|
439:     /  |          /  |          /  |    /  |          /  |          /  |
440:    4-----27------11----43------10  |   7-----43------6-----27------12  |
441:    |   |     18  |   |         |   |   |   |         |   |     18  |   |
442:    |  32         |  40         |  39   |  40         |  39         |   32
443:    |14 |         |16 |         |15 |   |15 |         |14 |         |16 |
444:   31   |  17    38   |         37  |   38  |        37   |  17    31   |
445:    |   3-----26--|---9-----42--|---8   |   5-----42--|---4-----26--|---11
446:    |  /          |  /          |  /    |  /          |  /          |  /
447:    | 29    19    | 34    21    |33     | 34    21    | 33    19    | 29
448:    |/            |/            |/      |/            |/            |/
449:    2-----25------7-----41------6       3-----41------2-----25------10

451:        5-----28------13----44------12      9-----44------8-----28------13
452:       /|            /|            /|      /|            /|            /|
453:     30 |   20      36|   22      35|     36|    22     35|   20      30|
454:     /  |          /  |          /  |    /  |          /  |          /  |
455:    4-----27------11----43------10  |   7-----43------6-----27------12  |
456:    |   |     18  |   |         |   |   |   |         |   |     18  |   |
457:    |  32         |  40         |   39  |   40        |  39         |   32
458:    |14 |         |16 |         |15 |   |15 |         |14 |         |16 |
459:   31   |  17    38   |         37  |   38  |        37   |  17    31   |
460:    |   3-----26--|---9-----42--|---8   |   5-----42--|---4-----26--|---11
461:    |  /          |  /          |  /    |  /          |  /          |  /
462:    | 29    19    | 34    21    |33     | 34    21    | 33    19    | 29
463:    |/            |/            |/      |/            |/            |/
464:    2-----25------7-----41------6       3-----41------2-----25------10
465:       cell 0         cell 1                cell 1        cell 0

467: Test hex_3: T-junction

469:       19-----52-----20-----53------21
470:       /|            /|            /|
471:     60 |   38      61|   41      62|
472:     /  |          /  |          /  |
473:   16-----50-----17-----51------18  |
474:    |   |     33  |   |     35  |   |
475:    |  70         |  72         |   74
476:    |25 |         |26 |         |27 |
477:   64   |  32    66   |  34    68   |
478:    |  13-----48--|--14-----49--|---15
479:    |  /|         |  /|         |  /|
480:    |57 |   37    | 58|   40    | 59|
481:    |/  |         |/  |         |/  |
482:   10-----46-----11-----47------12  |
483:    |   |     29  |   |     31  |   |
484:    |  69         |  71         |   73
485:    |22 |         |23 |         |24 |
486:   63   |  28    65   |   30   67   |
487:    |   7-----44--|---8-----45--|---9
488:    |  /          |  /          |  /
489:    | 54   36     | 55    39    | 56
490:    |/            |/            |/
491:    4-----42------5-----43------6
492:       cell 0         cell 1

494: becomes

496:       15----102-----28---112----___27-----73------16
497:       /|            /|         /   /             /|
498:     77 |   55     104|      ---  103    46      78|
499:     /  |          /  |     /     /             /  |
500:   13----101-----26---111--/----25-----72------14  |
501:    |   |     54  |   |  107   /           43  |   |
502:    |  81         |  108 / 51 /                |   82
503:    |40 |         |52 | /   105                |41 |
504:   79   |  53    106  |/   /            42    80   |
505:    |  21-----87--|--31---/-89------23-------/----/
506:    |  /|         |  /|  /         /|       /
507:    |91 |   47    |109|-- 49      93|  -----
508:    |/  |         |/ /|          /  | /
509:   17-----83-----29-----85------19----
510:    |   |         |   |         |   |
511:    |  120        |  121        |  122
512:    |   |         |26 |         |   |
513:  117   |        118  |        119  |
514:    |  22-----88--|--32-----90--|---24
515:    |  /|         |  /|         |  /|
516:    |92 |   48    |110|   50    | 94|
517:    |/  |         |/  |         |/  |
518:   18-----84-----30-----86------20  |
519:    |   |     37  |   |     39  |   |
520:    |  98         |  99         |   100
521:    |33 |         |34 |         |35 |
522:   95   |  36    96   |   38   97   |
523:    |  10-----70--|--11-----71--|---12
524:    |  /          |  /          |  /
525:    | 74   44     | 75    45    | 76
526:    |/            |/            |/
527:    7-----68------8-----69------9
528:       cell 0         cell 1

530: Test hex_4: Two non-intersecting faults

532:           cell 4         cell 5         cell 6        cell 7
533:       33-----96-----34-----97-----35-----98-----36-----99------37
534:       /|            /|            /|            /|            /|
535:     110|   66     111|   69     112|   72     113|   75     114|
536:     /  |          /  |          /  |          /  |          /  |
537:   28-----92-----29-----93-----30-----94-----31-----95------32  |
538:    |   |     57  |   |     59  |   |     61  |   |     63  |   |
539:    |  126        |  128        |  130        |  132        |  134
540:    |43 |         |44 |         |45 |         |46 |         |47 |
541:   116  |  56    118  |  58    120  |  60    122  |  62    124  |
542:    |  23-----88--|--24-----89--|--25-----90--|--26-----91--|---27
543:    |  /|         |  /|         |  /|         |  /|         |  /|
544:    |105|   65    |106|   68    |107|   71    |108|   74    |109|
545:    |/  |         |/  |         |/  |         |/  |         |/  |
546:   18-----84-----19-----95-----20-----86-----21-----87------22  |
547:    |   |     49  |   |     51  |   |     53  |   |     55  |   |
548:    |  125        |  127        |  129        |  131        |  133
549:    |38 |         |39 |         |40 |         |41 |         |42 |
550:   115  |  48    117  |  50    119  |  52    121  |  54    123  |
551:    |  13-----80--|--14-----81--|--15-----82--|--16-----83--|---17
552:    |  /          |  /          |  /          |  /          |  /
553:    |100    64    |101    67    |102    70    |103    73    |104
554:    |/            |/            |/            |/            |/
555:    8-----76------9-----77-----10-----78-----11-----79------12
556:       cell 0         cell 1        cell 2        cell 3

558: becomes

560:           cell 4         cell 5        cell 7        cell 10       cell 6
561:       27-----114----28-----115----29-----159----46-----170----45------116----30
562:       /|            /|            /|            /|            /|            /|
563:     123|   71     124|   73     125|   87     162|          161|    78    126|
564:     /  |          /  |          /  |          /  |          /  |          /  |
565:   23-----111----24-----112----25-----158----44-----169----43-----113-----26  |
566:    |   |     65  |   |    67   |   |    86   |   |         |   |     69  |   |
567:    |  134        |  135        |  137        |  166        |  165        |  140
568:    |56 |         |57 |         |58 |         |84 |         |83 |         |59 |
569:   127  |  64    128  |  66    130  |  85    164  |        163  |  68    133  |
570:    |  35-----143-|--37-----151-|--40-----109-|--42-----168-|--42-----110-|---22
571:    |  /|         |  /|         |  /|         |  /          |  /          |  /
572:    |145|   79    |147|   81    |153|   75    |160          |160    77    |122
573:    |/ 173        |/ 174        |/ 176        |/            |/            |/
574:   31-----141----33-----149----39-----107----41-----167----41-----108-----21
575: cell   |         |   |         |   | cell 9
576: 8  |  36-----144-|--38-----152-|--40-----109----42-----110-----22
577:   171 /|        172 /|        175 /|            /|            /|
578:    |146|   80    |148|   82    |153|    75    160|   77     122|
579:    |/  |         |/  |         |/  |          /  |          /  |
580:   32-----142----34-----150----39-----107----41-----108-----21  |
581:    |   |     50  |   |    52   |   |    61   |   |     63  |   |
582:    |  156        |  157        |  136        |  138        |  139
583:    |47 |         |48 |         |53 |         |54 |         |55 |
584:   154  |  49    155  |  51    129  |  60    131  |  62    132  |
585:    |  16-----103-|--17-----104-|--18-----105-|--19-----106-|---20
586:    |  /          |  /          |  /          |  /          |  /
587:    |117    70    |118    72    |119    74    |120    76    |121
588:    |/            |/            |/            |/            |/
589:   11-----99-----12-----100----13-----101----14-----102-----15
590:       cell 0         cell 1        cell 2        cell 3

592: */

594: typedef struct {
595:   PetscInt testNum;        // The mesh to test
596:   PetscInt cohesiveFields; // The number of fault fields
597: } AppCtx;

599: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
600: {
601:   PetscFunctionBegin;
602:   options->testNum        = 0;
603:   options->cohesiveFields = 1;

605:   PetscOptionsBegin(comm, "", "Cohesive Meshing Options", "DMPLEX");
606:   PetscCall(PetscOptionsBoundedInt("-test_num", "The particular mesh to test", __FILE__, options->testNum, &options->testNum, NULL, 0));
607:   PetscCall(PetscOptionsBoundedInt("-cohesive_fields", "The number of cohesive fields", __FILE__, options->cohesiveFields, &options->cohesiveFields, NULL, 0));
608:   PetscOptionsEnd();
609:   PetscFunctionReturn(PETSC_SUCCESS);
610: }

612: static PetscErrorCode CreateQuadMesh1(MPI_Comm comm, AppCtx *user, DM *dm)
613: {
614:   const PetscInt faces[2] = {1, 1};
615:   PetscReal      lower[2], upper[2];
616:   DMLabel        label;
617:   PetscMPIInt    rank;
618:   void          *get_tmp;
619:   PetscInt64    *cidx;
620:   PetscMPIInt    iflg;

622:   PetscFunctionBeginUser;
623:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
624:   // Create serial mesh
625:   lower[0] = (PetscReal)(rank % 2);
626:   lower[1] = (PetscReal)(rank / 2);
627:   upper[0] = (PetscReal)(rank % 2) + 1.;
628:   upper[1] = (PetscReal)(rank / 2) + 1.;
629:   PetscCall(DMPlexCreateBoxMesh(PETSC_COMM_SELF, 2, PETSC_FALSE, faces, lower, upper, NULL, PETSC_TRUE, 0, PETSC_TRUE, dm));
630:   PetscCall(PetscObjectSetName((PetscObject)*dm, "box"));
631:   // Flip edges to make fault non-oriented
632:   switch (rank) {
633:   case 2:
634:     PetscCall(DMPlexOrientPoint(*dm, 8, -1));
635:     break;
636:   case 3:
637:     PetscCall(DMPlexOrientPoint(*dm, 7, -1));
638:     break;
639:   default:
640:     break;
641:   }
642:   // Need this so that all procs create the cell types
643:   PetscCall(DMPlexGetCellTypeLabel(*dm, &label));
644:   // Replace comm in object (copied from PetscHeaderCreate/Destroy())
645:   PetscCall(PetscCommDestroy(&(*dm)->hdr.comm));
646:   PetscCall(PetscCommDuplicate(comm, &(*dm)->hdr.comm, &(*dm)->hdr.tag));
647:   PetscCallMPI(MPI_Comm_get_attr((*dm)->hdr.comm, Petsc_CreationIdx_keyval, &get_tmp, &iflg));
648:   PetscCheck(iflg, (*dm)->hdr.comm, PETSC_ERR_ARG_CORRUPT, "MPI_Comm does not have an object creation index");
649:   cidx            = (PetscInt64 *)get_tmp;
650:   (*dm)->hdr.cidx = (*cidx)++;
651:   // Create new pointSF
652:   {
653:     PetscSF      sf;
654:     PetscInt    *local  = NULL;
655:     PetscSFNode *remote = NULL;
656:     PetscInt     Nl;

658:     PetscCall(PetscSFCreate(comm, &sf));
659:     switch (rank) {
660:     case 0:
661:       Nl = 5;
662:       PetscCall(PetscMalloc1(Nl, &local));
663:       PetscCall(PetscMalloc1(Nl, &remote));
664:       local[0]        = 2;
665:       remote[0].index = 1;
666:       remote[0].rank  = 1;
667:       local[1]        = 3;
668:       remote[1].index = 1;
669:       remote[1].rank  = 2;
670:       local[2]        = 4;
671:       remote[2].index = 1;
672:       remote[2].rank  = 3;
673:       local[3]        = 6;
674:       remote[3].index = 5;
675:       remote[3].rank  = 2;
676:       local[4]        = 8;
677:       remote[4].index = 7;
678:       remote[4].rank  = 1;
679:       break;
680:     case 1:
681:       Nl = 3;
682:       PetscCall(PetscMalloc1(Nl, &local));
683:       PetscCall(PetscMalloc1(Nl, &remote));
684:       local[0]        = 3;
685:       remote[0].index = 1;
686:       remote[0].rank  = 3;
687:       local[1]        = 4;
688:       remote[1].index = 2;
689:       remote[1].rank  = 3;
690:       local[2]        = 6;
691:       remote[2].index = 5;
692:       remote[2].rank  = 3;
693:       break;
694:     case 2:
695:       Nl = 3;
696:       PetscCall(PetscMalloc1(Nl, &local));
697:       PetscCall(PetscMalloc1(Nl, &remote));
698:       local[0]        = 2;
699:       remote[0].index = 1;
700:       remote[0].rank  = 3;
701:       local[1]        = 4;
702:       remote[1].index = 3;
703:       remote[1].rank  = 3;
704:       local[2]        = 8;
705:       remote[2].index = 7;
706:       remote[2].rank  = 3;
707:       break;
708:     case 3:
709:       Nl = 0;
710:       break;
711:     default:
712:       SETERRQ(comm, PETSC_ERR_SUP, "This example only supports 4 ranks");
713:     }
714:     PetscCall(PetscSFSetGraph(sf, 9, Nl, local, PETSC_OWN_POINTER, remote, PETSC_OWN_POINTER));
715:     PetscCall(DMSetPointSF(*dm, sf));
716:     PetscCall(PetscSFDestroy(&sf));
717:   }
718:   // Create fault label
719:   PetscCall(DMCreateLabel(*dm, "fault"));
720:   PetscCall(DMGetLabel(*dm, "fault", &label));
721:   switch (rank) {
722:   case 0:
723:   case 2:
724:     PetscCall(DMLabelSetValue(label, 8, 1));
725:     PetscCall(DMLabelSetValue(label, 2, 0));
726:     PetscCall(DMLabelSetValue(label, 4, 0));
727:     break;
728:   case 1:
729:   case 3:
730:     PetscCall(DMLabelSetValue(label, 7, 1));
731:     PetscCall(DMLabelSetValue(label, 1, 0));
732:     PetscCall(DMLabelSetValue(label, 3, 0));
733:     break;
734:   default:
735:     break;
736:   }
737:   PetscCall(DMPlexOrientLabel(*dm, label));
738:   PetscCall(DMPlexLabelCohesiveComplete(*dm, label, NULL, 1, PETSC_FALSE, NULL));
739:   PetscCall(DMPlexDistributeSetDefault(*dm, PETSC_FALSE));
740:   PetscFunctionReturn(PETSC_SUCCESS);
741: }

743: static PetscErrorCode CreateHexMesh1(MPI_Comm comm, AppCtx *user, DM *dm)
744: {
745:   const PetscInt faces[3] = {1, 1, 1};
746:   PetscReal      lower[3], upper[3];
747:   DMLabel        label;
748:   PetscMPIInt    rank;
749:   void          *get_tmp;
750:   PetscInt64    *cidx;
751:   PetscMPIInt    iflg;

753:   PetscFunctionBeginUser;
754:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
755:   // Create serial mesh
756:   lower[0] = (PetscReal)(rank % 2);
757:   lower[1] = 0.;
758:   lower[2] = (PetscReal)(rank / 2);
759:   upper[0] = (PetscReal)(rank % 2) + 1.;
760:   upper[1] = 1.;
761:   upper[2] = (PetscReal)(rank / 2) + 1.;
762:   PetscCall(DMPlexCreateBoxMesh(PETSC_COMM_SELF, 3, PETSC_FALSE, faces, lower, upper, NULL, PETSC_TRUE, 0, PETSC_TRUE, dm));
763:   PetscCall(PetscObjectSetName((PetscObject)*dm, "box"));
764:   // Flip edges to make fault non-oriented
765:   switch (rank) {
766:   case 2:
767:     PetscCall(DMPlexOrientPoint(*dm, 10, -1));
768:     break;
769:   case 3:
770:     PetscCall(DMPlexOrientPoint(*dm, 9, -1));
771:     break;
772:   default:
773:     break;
774:   }
775:   // Need this so that all procs create the cell types
776:   PetscCall(DMPlexGetCellTypeLabel(*dm, &label));
777:   // Replace comm in object (copied from PetscHeaderCreate/Destroy())
778:   PetscCall(PetscCommDestroy(&(*dm)->hdr.comm));
779:   PetscCall(PetscCommDuplicate(comm, &(*dm)->hdr.comm, &(*dm)->hdr.tag));
780:   PetscCallMPI(MPI_Comm_get_attr((*dm)->hdr.comm, Petsc_CreationIdx_keyval, &get_tmp, &iflg));
781:   PetscCheck(iflg, (*dm)->hdr.comm, PETSC_ERR_ARG_CORRUPT, "MPI_Comm does not have an object creation index");
782:   cidx            = (PetscInt64 *)get_tmp;
783:   (*dm)->hdr.cidx = (*cidx)++;
784:   // Create new pointSF
785:   {
786:     PetscSF      sf;
787:     PetscInt    *local  = NULL;
788:     PetscSFNode *remote = NULL;
789:     PetscInt     Nl;

791:     PetscCall(PetscSFCreate(comm, &sf));
792:     switch (rank) {
793:     case 0:
794:       Nl = 15;
795:       PetscCall(PetscMalloc1(Nl, &local));
796:       PetscCall(PetscMalloc1(Nl, &remote));
797:       local[0]         = 2;
798:       remote[0].index  = 1;
799:       remote[0].rank   = 1;
800:       local[1]         = 4;
801:       remote[1].index  = 3;
802:       remote[1].rank   = 1;
803:       local[2]         = 5;
804:       remote[2].index  = 1;
805:       remote[2].rank   = 2;
806:       local[3]         = 6;
807:       remote[3].index  = 1;
808:       remote[3].rank   = 3;
809:       local[4]         = 7;
810:       remote[4].index  = 3;
811:       remote[4].rank   = 2;
812:       local[5]         = 8;
813:       remote[5].index  = 3;
814:       remote[5].rank   = 3;
815:       local[6]         = 17;
816:       remote[6].index  = 15;
817:       remote[6].rank   = 2;
818:       local[7]         = 18;
819:       remote[7].index  = 16;
820:       remote[7].rank   = 2;
821:       local[8]         = 20;
822:       remote[8].index  = 19;
823:       remote[8].rank   = 1;
824:       local[9]         = 21;
825:       remote[9].index  = 19;
826:       remote[9].rank   = 2;
827:       local[10]        = 22;
828:       remote[10].index = 19;
829:       remote[10].rank  = 3;
830:       local[11]        = 24;
831:       remote[11].index = 23;
832:       remote[11].rank  = 1;
833:       local[12]        = 26;
834:       remote[12].index = 25;
835:       remote[12].rank  = 1;
836:       local[13]        = 10;
837:       remote[13].index = 9;
838:       remote[13].rank  = 1;
839:       local[14]        = 14;
840:       remote[14].index = 13;
841:       remote[14].rank  = 2;
842:       break;
843:     case 1:
844:       Nl = 9;
845:       PetscCall(PetscMalloc1(Nl, &local));
846:       PetscCall(PetscMalloc1(Nl, &remote));
847:       local[0]        = 5;
848:       remote[0].index = 1;
849:       remote[0].rank  = 3;
850:       local[1]        = 6;
851:       remote[1].index = 2;
852:       remote[1].rank  = 3;
853:       local[2]        = 7;
854:       remote[2].index = 3;
855:       remote[2].rank  = 3;
856:       local[3]        = 8;
857:       remote[3].index = 4;
858:       remote[3].rank  = 3;
859:       local[4]        = 17;
860:       remote[4].index = 15;
861:       remote[4].rank  = 3;
862:       local[5]        = 18;
863:       remote[5].index = 16;
864:       remote[5].rank  = 3;
865:       local[6]        = 21;
866:       remote[6].index = 19;
867:       remote[6].rank  = 3;
868:       local[7]        = 22;
869:       remote[7].index = 20;
870:       remote[7].rank  = 3;
871:       local[8]        = 14;
872:       remote[8].index = 13;
873:       remote[8].rank  = 3;
874:       break;
875:     case 2:
876:       Nl = 9;
877:       PetscCall(PetscMalloc1(Nl, &local));
878:       PetscCall(PetscMalloc1(Nl, &remote));
879:       local[0]        = 2;
880:       remote[0].index = 1;
881:       remote[0].rank  = 3;
882:       local[1]        = 4;
883:       remote[1].index = 3;
884:       remote[1].rank  = 3;
885:       local[2]        = 6;
886:       remote[2].index = 5;
887:       remote[2].rank  = 3;
888:       local[3]        = 8;
889:       remote[3].index = 7;
890:       remote[3].rank  = 3;
891:       local[4]        = 20;
892:       remote[4].index = 19;
893:       remote[4].rank  = 3;
894:       local[5]        = 22;
895:       remote[5].index = 21;
896:       remote[5].rank  = 3;
897:       local[6]        = 24;
898:       remote[6].index = 23;
899:       remote[6].rank  = 3;
900:       local[7]        = 26;
901:       remote[7].index = 25;
902:       remote[7].rank  = 3;
903:       local[8]        = 10;
904:       remote[8].index = 9;
905:       remote[8].rank  = 3;
906:       break;
907:     case 3:
908:       Nl = 0;
909:       break;
910:     default:
911:       SETERRQ(comm, PETSC_ERR_SUP, "This example only supports 4 ranks");
912:     }
913:     PetscCall(PetscSFSetGraph(sf, 27, Nl, local, PETSC_OWN_POINTER, remote, PETSC_OWN_POINTER));
914:     PetscCall(DMSetPointSF(*dm, sf));
915:     PetscCall(PetscSFDestroy(&sf));
916:   }
917:   // Create fault label
918:   PetscCall(DMCreateLabel(*dm, "fault"));
919:   PetscCall(DMGetLabel(*dm, "fault", &label));
920:   switch (rank) {
921:   case 0:
922:   case 2:
923:     PetscCall(DMLabelSetValue(label, 10, 2));
924:     PetscCall(DMLabelSetValue(label, 20, 1));
925:     PetscCall(DMLabelSetValue(label, 22, 1));
926:     PetscCall(DMLabelSetValue(label, 24, 1));
927:     PetscCall(DMLabelSetValue(label, 26, 1));
928:     PetscCall(DMLabelSetValue(label, 2, 0));
929:     PetscCall(DMLabelSetValue(label, 4, 0));
930:     PetscCall(DMLabelSetValue(label, 6, 0));
931:     PetscCall(DMLabelSetValue(label, 8, 0));
932:     break;
933:   case 1:
934:   case 3:
935:     PetscCall(DMLabelSetValue(label, 9, 2));
936:     PetscCall(DMLabelSetValue(label, 19, 1));
937:     PetscCall(DMLabelSetValue(label, 21, 1));
938:     PetscCall(DMLabelSetValue(label, 23, 1));
939:     PetscCall(DMLabelSetValue(label, 25, 1));
940:     PetscCall(DMLabelSetValue(label, 1, 0));
941:     PetscCall(DMLabelSetValue(label, 3, 0));
942:     PetscCall(DMLabelSetValue(label, 5, 0));
943:     PetscCall(DMLabelSetValue(label, 7, 0));
944:     break;
945:   default:
946:     break;
947:   }
948:   PetscCall(DMPlexOrientLabel(*dm, label));
949:   PetscCall(DMPlexLabelCohesiveComplete(*dm, label, NULL, 1, PETSC_FALSE, NULL));
950:   PetscCall(DMPlexDistributeSetDefault(*dm, PETSC_FALSE));
951:   PetscFunctionReturn(PETSC_SUCCESS);
952: }

954: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
955: {
956:   PetscFunctionBegin;
957:   switch (user->testNum) {
958:   case 1:
959:     PetscCall(CreateQuadMesh1(comm, user, dm));
960:     break;
961:   case 2:
962:     PetscCall(CreateHexMesh1(comm, user, dm));
963:     break;
964:   default:
965:     PetscCall(DMCreate(comm, dm));
966:     PetscCall(DMSetType(*dm, DMPLEX));
967:     break;
968:   }
969:   PetscCall(DMSetFromOptions(*dm));
970:   {
971:     const char *prefix;

973:     // We cannot redistribute with cohesive cells in the SF
974:     PetscCall(DMPlexDistributeSetDefault(*dm, PETSC_FALSE));
975:     PetscCall(PetscObjectGetOptionsPrefix((PetscObject)*dm, &prefix));
976:     PetscCall(PetscObjectSetOptionsPrefix((PetscObject)*dm, "f0_"));
977:     PetscCall(DMSetFromOptions(*dm));
978:     PetscCall(PetscObjectSetOptionsPrefix((PetscObject)*dm, "f1_"));
979:     PetscCall(DMSetFromOptions(*dm));
980:     PetscCall(PetscObjectSetOptionsPrefix((PetscObject)*dm, prefix));
981:   }
982:   PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
983:   PetscFunctionReturn(PETSC_SUCCESS);
984: }

986: // Create a displacement field, and some number of vector fault fields
987: static PetscErrorCode CreateDiscretization(DM dm, AppCtx *user)
988: {
989:   PetscSection   s;
990:   DMLabel        fault, faultSpace;
991:   PetscFE        fe;
992:   DMPolytopeType ct, fct;
993:   PetscInt       dim, cStart, fStart, Ncf = user->cohesiveFields;

995:   PetscFunctionBegin;
996:   PetscCall(DMGetDimension(dm, &dim));
997:   PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, NULL));
998:   PetscCall(DMPlexGetCellType(dm, cStart, &ct));
999:   PetscCall(DMGetLabel(dm, "fault", &fault));
1000:   if (!fault) PetscFunctionReturn(PETSC_SUCCESS);
1001:   PetscCall(DMLabelView(fault, PETSC_VIEWER_STDOUT_WORLD));

1003:   PetscCall(PetscFECreateByCell(PETSC_COMM_SELF, dim, dim, ct, "displacement_", PETSC_DETERMINE, &fe));
1004:   PetscCall(PetscFESetName(fe, "displacement"));
1005:   PetscCall(DMAddField(dm, NULL, (PetscObject)fe));
1006:   PetscCall(PetscFEDestroy(&fe));

1008:   // Make label for fault space definition
1009:   PetscCall(DMCreateLabel(dm, "faultSpace"));
1010:   PetscCall(DMGetLabel(dm, "faultSpace", &faultSpace));
1011:   for (PetscInt d = 0; d <= dim; ++d) {
1012:     PetscInt pStart, pEnd, pMax;

1014:     PetscCall(DMPlexGetSimplexOrBoxCells(dm, d, NULL, &pMax));
1015:     PetscCall(DMPlexGetHeightStratum(dm, d, &pStart, &pEnd));
1016:     for (PetscInt p = pMax; p < pEnd; ++p) PetscCall(DMLabelSetValue(faultSpace, p, 1));
1017:   }
1018:   PetscCall(DMPlexGetHeightStratum(dm, 1, &fStart, NULL));
1019:   PetscCall(DMPlexGetCellType(dm, fStart, &fct));
1020:   if (Ncf > 0) {
1021:     PetscCall(PetscFECreateByCell(PETSC_COMM_SELF, dim - 1, dim, fct, "faulttraction_", PETSC_DETERMINE, &fe));
1022:     PetscCall(PetscFESetName(fe, "fault traction"));
1023:     PetscCall(DMAddField(dm, faultSpace, (PetscObject)fe));
1024:     PetscCall(PetscFEDestroy(&fe));
1025:   }
1026:   for (PetscInt f = 1; f < Ncf; ++f) {
1027:     char name[256], opt[256];

1029:     PetscCall(PetscSNPrintf(name, 256, "fault field %" PetscInt_FMT, f));
1030:     PetscCall(PetscSNPrintf(opt, 256, "faultfield_%" PetscInt_FMT "_", f));
1031:     PetscCall(PetscFECreateByCell(PETSC_COMM_SELF, dim - 1, dim, fct, opt, PETSC_DETERMINE, &fe));
1032:     PetscCall(PetscFESetName(fe, name));
1033:     PetscCall(DMAddField(dm, faultSpace, (PetscObject)fe));
1034:     PetscCall(PetscFEDestroy(&fe));
1035:   }
1036:   PetscCall(DMCreateDS(dm));

1038:   PetscCall(DMGetLocalSection(dm, &s));
1039:   PetscCall(PetscObjectViewFromOptions((PetscObject)s, NULL, "-local_section_view"));
1040:   PetscFunctionReturn(PETSC_SUCCESS);
1041: }

1043: // Label cells 1 for negative side, and 2 for positive side
1044: static PetscErrorCode CreateMaterialLabel(DM dm)
1045: {
1046:   DMLabel         fault, material;
1047:   IS              faceIS;
1048:   const PetscInt *faces;
1049:   PetscReal       fvol, fcentroid[3], fnormal[3];
1050:   PetscInt        dim, cStart, cEnd, Nf;

1052:   PetscFunctionBegin;
1053:   PetscCall(DMGetDimension(dm, &dim));
1054:   PetscCall(DMGetLabel(dm, "fault", &fault));
1055:   PetscCall(DMCreateLabel(dm, "material"));
1056:   PetscCall(DMGetLabel(dm, "material", &material));
1057:   for (PetscInt s = 1; s < 3; ++s) {
1058:     IS              pointIS;
1059:     const PetscInt *points;
1060:     PetscInt        n;

1062:     PetscCall(DMLabelGetStratumIS(fault, s > 1 ? 100 + dim : -(100 + dim), &pointIS));
1063:     if (!pointIS) continue;
1064:     PetscCall(ISGetLocalSize(pointIS, &n));
1065:     PetscCall(ISGetIndices(pointIS, &points));
1066:     for (PetscInt i = 0; i < n; ++i) {
1067:       PetscCall(DMLabelSetValue(material, points[i], s));
1068:     }
1069:     PetscCall(ISRestoreIndices(pointIS, &points));
1070:     PetscCall(ISDestroy(&pointIS));
1071:   }
1072:   // This simple algorithm will work for now (note that cohesive cells get added into this label)
1073:   PetscCall(DMLabelGetStratumIS(fault, dim - 1, &faceIS));
1074:   PetscCall(ISGetLocalSize(faceIS, &Nf));
1075:   PetscCheck(Nf > 0, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Fault label must contain at least one face");
1076:   PetscCall(ISGetIndices(faceIS, &faces));
1077:   for (PetscInt i = 0; i < Nf; ++i) {
1078:     const PetscInt face = faces[i];
1079:     DMPolytopeType ct;

1081:     PetscCall(DMPlexGetCellType(dm, face, &ct));
1082:     if (DMPolytopeTypeGetDim(ct) != dim - 1) continue;
1083:     PetscCall(DMPlexComputeCellGeometryFVM(dm, face, &fvol, fcentroid, fnormal));
1084:     break;
1085:   }
1086:   PetscCall(ISRestoreIndices(faceIS, &faces));
1087:   PetscCall(ISDestroy(&faceIS));
1088:   PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd));
1089:   for (PetscInt c = cStart; c < cEnd; ++c) {
1090:     PetscReal vol, centroid[3];
1091:     PetscInt  val;

1093:     PetscCall(DMLabelGetValue(fault, c, &val));
1094:     if (val >= 0) continue;
1095:     PetscCall(DMPlexComputeCellGeometryFVM(dm, c, &vol, centroid, NULL));
1096:     for (PetscInt e = 0; e < dim; ++e) centroid[e] -= fcentroid[e];
1097:     if (DMPlex_DotRealD_Internal(dim, centroid, fnormal) > 0) PetscCall(DMLabelSetValue(material, c, 1));
1098:     else PetscCall(DMLabelSetValue(material, c, 2));
1099:   }
1100:   PetscFunctionReturn(PETSC_SUCCESS);
1101: }

1103: // Label cohesive cells and endcap faces 1
1104: static PetscErrorCode CreateFaultLabel(DM dm)
1105: {
1106:   DMLabel  fault;
1107:   PetscInt cMax, cEnd;

1109:   PetscFunctionBegin;
1110:   PetscCall(DMCreateLabel(dm, "faultCells"));
1111:   PetscCall(DMGetLabel(dm, "faultCells", &fault));
1112:   PetscCall(DMPlexGetSimplexOrBoxCells(dm, 0, &cEnd, &cMax));
1113:   for (PetscInt c = cMax; c < cEnd; ++c) {
1114:     const PetscInt *cone;

1116:     PetscCall(DMLabelSetValue(fault, c, 1));
1117:     PetscCall(DMPlexGetCone(dm, c, &cone));
1118:     PetscCall(DMLabelSetValue(fault, cone[0], 1));
1119:     PetscCall(DMLabelSetValue(fault, cone[1], 1));
1120:   }
1121:   PetscFunctionReturn(PETSC_SUCCESS);
1122: }

1124: static PetscErrorCode r(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, PetscCtx ctx)
1125: {
1126:   PetscInt d;
1127:   for (d = 0; d < dim; ++d) u[d] = x[d];
1128:   return PETSC_SUCCESS;
1129: }

1131: static PetscErrorCode rp1(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, PetscCtx ctx)
1132: {
1133:   PetscInt d;
1134:   for (d = 0; d < dim; ++d) u[d] = x[d] + (d > 0 ? 1.0 : 0.0);
1135:   return PETSC_SUCCESS;
1136: }

1138: static PetscErrorCode phi(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, PetscCtx ctx)
1139: {
1140:   PetscInt d;
1141:   u[0] = -x[1];
1142:   u[1] = x[0];
1143:   for (d = 2; d < dim; ++d) u[d] = x[d];
1144:   return PETSC_SUCCESS;
1145: }

1147: static void add_fields(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f[])
1148: {
1149:   PetscInt       d;
1150:   const PetscInt offN = 0;
1151:   const PetscInt offP = dim;
1152:   for (d = 0; d < dim; ++d) f[d] = u[offN + d] + u[offP + d];
1153: }

1155: static void normal_field(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f[])
1156: {
1157:   PetscInt d;
1158:   for (d = 0; d < dim; ++d) f[d] = n[d];
1159: }

1161: /* \lambda \cdot (\psi_u^- - \psi_u^+) */
1162: static void f0_bd_u_neg(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
1163: {
1164:   const PetscInt Nc = dim + 1;
1165:   for (PetscInt c = 0; c < Nc; ++c) f0[c] = -u[uOff[1] + c];
1166: }

1168: static void f0_bd_u_pos(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
1169: {
1170:   const PetscInt Nc = dim + 1;
1171:   for (PetscInt c = 0; c < Nc; ++c) f0[c] = u[uOff[1] + c];
1172: }

1174: /* (d - u^+ + u^-) \cdot \psi_\lambda */
1175: static void f0_bd_l(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
1176: {
1177:   const PetscInt Nc = uOff[2] - uOff[1];

1179:   for (PetscInt c = 0; c < Nc; ++c) f0[c] = (c > 0 ? 1.0 : 0.0) + u[c] - u[Nc + c];
1180: }

1182: /* \psi_lambda \cdot (\psi_u^- - \psi_u^+) */
1183: static void g0_bd_ul_neg(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[])
1184: {
1185:   const PetscInt Nc = dim + 1;
1186:   for (PetscInt c = 0; c < Nc; ++c) g0[c * Nc + c] = -1.0;
1187: }

1189: static void g0_bd_ul_pos(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[])
1190: {
1191:   const PetscInt Nc = dim + 1;
1192:   for (PetscInt c = 0; c < Nc; ++c) g0[c * Nc + c] = 1.0;
1193: }

1195: /* (-\psi_u^+ + \psi_u^-) \cdot \psi_\lambda */
1196: static void g0_bd_lu(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[])
1197: {
1198:   const PetscInt Nc = uOff[2] - uOff[1];

1200:   for (PetscInt c = 0; c < Nc; ++c) {
1201:     g0[c * Nc + c]           = -1.0;
1202:     g0[Nc * Nc + c * Nc + c] = 1.0;
1203:   }
1204: }

1206: static PetscErrorCode TestAssembly(DM dm, AppCtx *user)
1207: {
1208:   Mat           J;
1209:   Vec           locX, locF, locW;
1210:   PetscDS       probh;
1211:   DMLabel       fault, material;
1212:   DM            dmFault;
1213:   IS            cohesiveCells;
1214:   PetscFE       fe;
1215:   PetscWeakForm wf;
1216:   PetscFormKey  keys[3];
1217:   PetscErrorCode (*initialGuess[2])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], PetscCtx ctx);
1218:   DMPolytopeType fct;
1219:   PetscInt       dim, fStart, Nf, cMax, cEnd, id;
1220:   PetscMPIInt    rank, size;

1222:   PetscFunctionBegin;
1223:   PetscCall(DMGetNumFields(dm, &Nf));
1224:   if (Nf <= 0) PetscFunctionReturn(PETSC_SUCCESS);
1225:   PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
1226:   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)dm), &size));
1227:   PetscCall(DMGetDimension(dm, &dim));
1228:   PetscCall(DMPlexGetSimplexOrBoxCells(dm, 0, NULL, &cMax));
1229:   PetscCall(DMPlexGetHeightStratum(dm, 0, NULL, &cEnd));
1230:   if (size > 1) {
1231:     PetscSF         sf;
1232:     const PetscInt *leaves;
1233:     PetscInt       *points;
1234:     PetscInt        Nl, l, Ncoh = 0;

1236:     PetscCall(DMGetPointSF(dm, &sf));
1237:     PetscCall(PetscSFGetGraph(sf, NULL, &Nl, &leaves, NULL));
1238:     for (PetscInt c = cMax; c < cEnd; ++c) {
1239:       PetscCall(PetscFindInt(c, Nl, leaves, &l));
1240:       if (l < 0) ++Ncoh;
1241:     }
1242:     PetscCall(PetscMalloc1(Ncoh, &points));
1243:     Ncoh = 0;
1244:     for (PetscInt c = cMax; c < cEnd; ++c) {
1245:       PetscCall(PetscFindInt(c, Nl, leaves, &l));
1246:       if (l < 0) points[Ncoh++] = c;
1247:     }
1248:     PetscCall(ISCreateGeneral(PETSC_COMM_SELF, Ncoh, points, PETSC_OWN_POINTER, &cohesiveCells));
1249:   } else {
1250:     PetscCall(ISCreateStride(PETSC_COMM_SELF, cEnd - cMax, cMax, 1, &cohesiveCells));
1251:   }
1252:   PetscCall(CreateFaultLabel(dm));
1253:   PetscCall(DMGetLabel(dm, "faultCells", &fault));
1254:   PetscCall(DMGetLocalVector(dm, &locX));
1255:   PetscCall(PetscObjectSetName((PetscObject)locX, "Local Solution"));
1256:   PetscCall(DMGetLocalVector(dm, &locF));
1257:   PetscCall(PetscObjectSetName((PetscObject)locF, "Local Residual"));
1258:   PetscCall(DMCreateMatrix(dm, &J));
1259:   PetscCall(PetscObjectSetName((PetscObject)J, "Jacobian"));

1261:   /* The initial guess has displacement shifted by one unit in each fault parallel direction across the fault */
1262:   PetscCall(CreateMaterialLabel(dm));
1263:   PetscCall(DMGetLabel(dm, "material", &material));
1264:   id              = 1;
1265:   initialGuess[0] = r;
1266:   initialGuess[1] = NULL;
1267:   PetscCall(DMProjectFunctionLabelLocal(dm, 0.0, material, 1, &id, PETSC_DETERMINE, NULL, initialGuess, NULL, INSERT_VALUES, locX));
1268:   id              = 2;
1269:   initialGuess[0] = rp1;
1270:   initialGuess[1] = NULL;
1271:   PetscCall(DMProjectFunctionLabelLocal(dm, 0.0, material, 1, &id, PETSC_DETERMINE, NULL, initialGuess, NULL, INSERT_VALUES, locX));
1272:   id              = 1;
1273:   initialGuess[0] = NULL;
1274:   initialGuess[1] = phi;
1275:   PetscCall(DMProjectFunctionLabelLocal(dm, 0.0, fault, 1, &id, PETSC_DETERMINE, NULL, initialGuess, NULL, INSERT_VALUES, locX));
1276:   PetscCall(PetscObjectViewSynchronizedFromOptions((PetscObject)locX, (PetscObject)dm, "-local_solution_view"));

1278:   // Test projection to fault mesh
1279:   if (cMax < cEnd) {
1280:     PetscCall(DMPlexCreateCohesiveSubmesh(dm, PETSC_FALSE, NULL, 0, &dmFault));
1281:     PetscCall(PetscObjectSetName((PetscObject)dmFault, "Fault Mesh"));
1282:     PetscCall(DMViewFromOptions(dmFault, NULL, "-fault_view"));
1283:     PetscCall(DMPlexOrient(dmFault));
1284:     PetscCall(DMPlexGetHeightStratum(dm, 1, &fStart, NULL));
1285:     PetscCall(DMPlexGetCellType(dm, fStart, &fct));
1286:     //PetscCall(PetscFECreateByCell(PETSC_COMM_SELF, dim - 1, dim, fct, "fault_field_", PETSC_DETERMINE, &fe));
1287:     PetscCall(PetscFECreateDefault(PETSC_COMM_SELF, dim - 1, dim, PETSC_TRUE, "fault_field_", PETSC_DETERMINE, &fe));
1288:     PetscCall(PetscFESetName(fe, "fault_field"));
1289:     PetscCall(DMAddField(dmFault, NULL, (PetscObject)fe));
1290:     PetscCall(PetscFEDestroy(&fe));
1291:     PetscCall(DMCreateDS(dmFault));
1292:     PetscCall(DMGetLocalVector(dmFault, &locW));
1293:     PetscCall(DMViewFromOptions(dmFault, NULL, "-cohesive_view"));
1294:     void (*faultFuncs[1])(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f[]);

1296:     DMLabel  depthLabel;
1297:     PetscInt depth;
1298:     PetscCall(DMPlexGetDepthLabel(dmFault, &depthLabel));
1299:     PetscCall(DMPlexGetDepth(dmFault, &depth));
1300:     id = depth - 1;
1301:     /* w = r + rp1 */
1302:     faultFuncs[0] = add_fields;
1303:     PetscCall(DMProjectBdFieldLabelLocal(dmFault, 0.0, depthLabel, 1, &id, PETSC_DETERMINE, NULL, locX, faultFuncs, INSERT_VALUES, locW));
1304:     PetscCall(PetscObjectViewSynchronizedFromOptions((PetscObject)locW, (PetscObject)dm, "-local_projection_view"));

1306:     /* w = fault_normal */
1307:     faultFuncs[0] = normal_field;
1308:     PetscCall(DMProjectBdFieldLabelLocal(dmFault, 0.0, depthLabel, 1, &id, PETSC_DETERMINE, NULL, locX, faultFuncs, INSERT_VALUES, locW));
1309:     PetscCall(PetscObjectViewSynchronizedFromOptions((PetscObject)locW, (PetscObject)dm, "-local_projection_view"));
1310:     PetscCall(DMRestoreLocalVector(dmFault, &locW));
1311:     PetscCall(DMDestroy(&dmFault));
1312:   }

1314:   PetscCall(DMGetCellDS(dm, cMax, &probh, NULL));
1315:   PetscCall(PetscDSGetWeakForm(probh, &wf));
1316:   PetscCall(PetscDSGetNumFields(probh, &Nf));
1317:   PetscCall(PetscWeakFormSetIndexBdResidual(wf, material, 1, 0, 0, 0, f0_bd_u_neg, 0, NULL));
1318:   PetscCall(PetscWeakFormSetIndexBdResidual(wf, material, 2, 0, 0, 0, f0_bd_u_pos, 0, NULL));
1319:   PetscCall(PetscWeakFormSetIndexBdJacobian(wf, material, 1, 0, 1, 0, 0, g0_bd_ul_neg, 0, NULL, 0, NULL, 0, NULL));
1320:   PetscCall(PetscWeakFormSetIndexBdJacobian(wf, material, 2, 0, 1, 0, 0, g0_bd_ul_pos, 0, NULL, 0, NULL, 0, NULL));
1321:   if (Nf > 1) {
1322:     PetscCall(PetscWeakFormSetIndexBdResidual(wf, fault, 1, 1, 0, 0, f0_bd_l, 0, NULL));
1323:     PetscCall(PetscWeakFormSetIndexBdJacobian(wf, fault, 1, 1, 0, 0, 0, g0_bd_lu, 0, NULL, 0, NULL, 0, NULL));
1324:   }
1325:   if (rank == 0) PetscCall(PetscDSView(probh, NULL));

1327:   keys[0].label = material;
1328:   keys[0].value = 1;
1329:   keys[0].field = 0;
1330:   keys[0].part  = 0;
1331:   keys[1].label = material;
1332:   keys[1].value = 2;
1333:   keys[1].field = 0;
1334:   keys[1].part  = 0;
1335:   keys[2].label = fault;
1336:   keys[2].value = 1;
1337:   keys[2].field = 1;
1338:   keys[2].part  = 0;
1339:   PetscCall(VecSet(locF, 0.));
1340:   PetscCall(DMPlexComputeResidualHybridByKey(dm, keys, cohesiveCells, 0.0, locX, NULL, 0.0, locF, user));
1341:   PetscCall(PetscObjectViewSynchronizedFromOptions((PetscObject)locF, (PetscObject)dm, "-local_residual_view"));
1342:   PetscCall(MatZeroEntries(J));
1343:   PetscCall(DMPlexComputeJacobianHybridByKey(dm, keys, cohesiveCells, 0.0, 0.0, locX, NULL, J, J, user));
1344:   PetscCall(MatAssemblyBegin(J, MAT_FINAL_ASSEMBLY));
1345:   PetscCall(MatAssemblyEnd(J, MAT_FINAL_ASSEMBLY));
1346:   PetscCall(MatViewFromOptions(J, NULL, "-local_jacobian_view"));

1348:   PetscCall(DMRestoreLocalVector(dm, &locX));
1349:   PetscCall(DMRestoreLocalVector(dm, &locF));
1350:   PetscCall(MatDestroy(&J));
1351:   PetscCall(ISDestroy(&cohesiveCells));

1353:   if (cMax < cEnd) {
1354:     PetscDS         ds;
1355:     PetscFE         fe;
1356:     PetscQuadrature quad;
1357:     IS             *perm;
1358:     const PetscInt *cone;
1359:     PetscInt        Na, a;

1361:     PetscCall(DMPlexGetCone(dm, cMax, &cone));
1362:     PetscCall(DMGetCellDS(dm, cMax, &ds, NULL));
1363:     PetscCall(PetscDSGetDiscretization(ds, 0, (PetscObject *)&fe));
1364:     PetscCall(PetscFEGetQuadrature(fe, &quad));
1365:     PetscCall(PetscQuadratureComputePermutations(quad, &Na, &perm));
1366:     for (a = 0; a < Na; ++a) PetscCall(ISDestroy(&perm[a]));
1367:     PetscCall(PetscFree(perm));
1368:   }
1369:   PetscFunctionReturn(PETSC_SUCCESS);
1370: }

1372: int main(int argc, char **argv)
1373: {
1374:   DM     dm;
1375:   AppCtx user;

1377:   PetscFunctionBeginUser;
1378:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
1379:   PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user));
1380:   PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm));
1381:   PetscCall(CreateDiscretization(dm, &user));
1382:   PetscCall(TestAssembly(dm, &user));
1383:   PetscCall(DMDestroy(&dm));
1384:   PetscCall(PetscFinalize());
1385:   return 0;
1386: }

1388: /*TEST

1390:   testset:
1391:     requires: triangle
1392:     args: -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1393:             -dm_plex_transform_active fault -dm_plex_save_transform -dm_plex_check_transform \
1394:           -dm_view ::ascii_info_detail -coarse_dm_view ::ascii_info_detail \
1395:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1396:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1397:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1399:     test:
1400:       suffix: tri_0
1401:       args: -dm_plex_box_faces 1,1 -dm_plex_cohesive_label_fault 8
1402:     test:
1403:       suffix: tri_1
1404:       args: -dm_plex_box_faces 1,1 -dm_plex_cohesive_label_fault 8 \
1405:               -dm_plex_transform_extrude_use_tensor 0
1406:     test:
1407:       suffix: tri_2
1408:       args: -dm_plex_file_contents dat:tri_2_cv -dm_plex_cohesive_label_fault 11,15
1409:     test:
1410:       suffix: tri_2_perm
1411:       args: -dm_plex_file_contents dat:tri_2_cv -dm_plex_cohesive_label_fault 11,15 \
1412:             -dm_reorder_section -dm_reorder_section_type cohesive
1413:     # Note that the mesh is not parallel when the cohesive label is oriented
1414:     test:
1415:       suffix: tri_3
1416:       nsize: 2
1417:       args: -dm_plex_file_contents dat:tri_2_cv -dm_plex_cohesive_label_fault 11,15 \
1418:               -petscpartitioner_type shell -petscpartitioner_shell_sizes 2,2 \
1419:               -petscpartitioner_shell_points 0,3,1,2

1421:   testset:
1422:     requires: triangle
1423:     args: -dm_plex_option_phases coh_,ref_ \
1424:             -coh_dm_refine 1 -coh_dm_plex_transform_type cohesive_extrude \
1425:               -coh_dm_plex_transform_active fault \
1426:             -ref_dm_refine 1 -ref_dm_plex_transform_type refine_regular \
1427:           -dm_view ::ascii_info_detail \
1428:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1429:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1430:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1432:     test:
1433:       suffix: tri_0_ref
1434:       args: -dm_plex_box_faces 1,1 -dm_plex_cohesive_label_fault 8
1435:     test:
1436:       suffix: tri_2_ref
1437:       args: -dm_plex_file_contents dat:tri_2_cv -dm_plex_cohesive_label_fault 11,15
1438:     test:
1439:       suffix: tri_2_ref_perm
1440:       args: -dm_plex_file_contents dat:tri_2_cv -dm_plex_cohesive_label_fault 11,15 \
1441:             -dm_reorder_section -dm_reorder_section_type cohesive

1443:   testset:
1444:     args: -dm_plex_simplex 0 -dm_plex_box_faces 2,1 \
1445:           -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1446:             -dm_plex_transform_active fault -dm_plex_cohesive_label_fault 13 \
1447:             -dm_plex_save_transform -dm_plex_check_transform \
1448:           -dm_view ::ascii_info_detail -coarse_dm_view ::ascii_info_detail \
1449:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1450:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1451:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1453:     test:
1454:       suffix: quad_0
1455:     test:
1456:       suffix: quad_1
1457:       args: -dm_plex_transform_extrude_use_tensor 0
1458:     test:
1459:       suffix: quad_2
1460:       nsize: 2
1461:       args: -petscpartitioner_type simple

1463:   test:
1464:     suffix: quad_3
1465:     nsize: 4
1466:     args: -test_num 1 \
1467:           -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1468:             -dm_plex_transform_active fault -dm_plex_save_transform -dm_plex_check_transform \
1469:           -dm_view ::ascii_info_detail -coarse_dm_view ::ascii_info_detail \
1470:           -orientation_view -orientation_view_synchronized \
1471:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1472:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1473:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1475:   test:
1476:     suffix: quad_4
1477:     args: -dm_plex_simplex 0 -dm_plex_box_faces 3,2 \
1478:           -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1479:             -dm_plex_transform_active fault -dm_plex_cohesive_label_fault 22,23 \
1480:             -dm_plex_save_transform -dm_plex_check_transform \
1481:           -dm_view ::ascii_info_detail -coarse_dm_view ::ascii_info_detail \
1482:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1483:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1484:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1486:   test:
1487:     suffix: quad_5
1488:     args: -dm_plex_simplex 0 -dm_plex_box_faces 3,2 \
1489:             -dm_plex_cohesive_label_fault0 21 \
1490:             -dm_plex_cohesive_label_fault1 23 \
1491:           -f0_dm_refine 1 -f0_dm_plex_transform_type cohesive_extrude \
1492:             -f0_dm_plex_transform_active fault0  -f0_coarse_dm_view ::ascii_info_detail \
1493:           -f1_dm_refine 1 -f1_dm_plex_transform_type cohesive_extrude \
1494:             -f1_dm_plex_transform_active fault1  -f1_coarse_dm_view ::ascii_info_detail \
1495:           -dm_plex_save_transform -dm_plex_check_transform \
1496:           -dm_view ::ascii_info_detail
1497:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1499:   test:
1500:     suffix: quad_6
1501:     args: -dm_plex_simplex 0 -dm_plex_box_faces 3,2 \
1502:             -dm_plex_cohesive_label_fault0 22,23 \
1503:             -dm_plex_cohesive_label_fault1 32 \
1504:           -f0_dm_refine 1 -f0_dm_plex_transform_type cohesive_extrude \
1505:             -f0_dm_plex_transform_active fault0  -f0_coarse_dm_view ::ascii_info_detail \
1506:           -f1_dm_refine 1 -f1_dm_plex_transform_type cohesive_extrude \
1507:             -f1_dm_plex_transform_active fault1  -f1_coarse_dm_view ::ascii_info_detail \
1508:           -dm_plex_save_transform -dm_plex_check_transform \
1509:           -dm_view ::ascii_info_detail
1510:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1512:   test:
1513:     suffix: quad_6w
1514:     args: -dm_plex_simplex 0 -dm_plex_box_faces 3,2 \
1515:             -dm_plex_cohesive_label_fault0 22,23 \
1516:             -dm_plex_cohesive_label_fault1 32 \
1517:           -f0_dm_refine 1 -f0_dm_plex_transform_type cohesive_extrude \
1518:             -f0_dm_plex_transform_active fault0  -f0_coarse_dm_view ::ascii_info_detail \
1519:             -f0_dm_plex_transform_cohesive_width 0.05 \
1520:           -f1_dm_refine 1 -f1_dm_plex_transform_type cohesive_extrude \
1521:             -f1_dm_plex_transform_active fault1  -f1_coarse_dm_view ::ascii_info_detail \
1522:             -f1_dm_plex_transform_cohesive_width 0.05 \
1523:           -dm_plex_save_transform -dm_plex_check_transform \
1524:           -dm_view ::ascii_info_detail
1525:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1527:   testset:
1528:     args: -dm_plex_simplex 0 -dm_plex_box_faces 2,1 -dm_plex_cohesive_label_fault 13 \
1529:           -dm_plex_option_phases coh_,ref_ \
1530:             -coh_dm_refine 1 -coh_dm_plex_transform_type cohesive_extrude \
1531:               -coh_dm_plex_transform_active fault \
1532:             -ref_dm_refine 1 -ref_dm_plex_transform_type refine_regular \
1533:           -dm_view ::ascii_info_detail \
1534:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1535:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1536:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1538:     test:
1539:       suffix: quad_0_ref

1541:   testset:
1542:     args: -dm_plex_dim 3 -dm_plex_shape doublet \
1543:           -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1544:             -dm_plex_transform_active fault -dm_plex_cohesive_label_fault 7 \
1545:             -dm_plex_save_transform -dm_plex_check_transform \
1546:           -dm_view ::ascii_info_detail -coarse_dm_view ::ascii_info_detail \
1547:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1548:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1549:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1551:     test:
1552:       suffix: tet_0
1553:     test:
1554:       suffix: tet_1
1555:       nsize: 2
1556:       args: -petscpartitioner_type simple

1558:   testset:
1559:     args: -dm_plex_dim 3 -dm_plex_shape doublet -dm_plex_cohesive_label_fault 7 \
1560:           -dm_plex_option_phases coh_,ref_ \
1561:             -coh_dm_refine 1 -coh_dm_plex_transform_type cohesive_extrude \
1562:               -coh_dm_plex_transform_active fault \
1563:             -ref_dm_refine 1 -ref_dm_plex_transform_type refine_regular \
1564:           -dm_view ::ascii_info_detail \
1565:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1566:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1567:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1569:     test:
1570:       suffix: tet_0_ref

1572:   testset:
1573:     args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_plex_box_faces 2,1,1 -dm_plex_box_upper 2,1,1 \
1574:           -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1575:             -dm_plex_transform_active fault -dm_plex_cohesive_label_fault 15 \
1576:             -dm_plex_save_transform -dm_plex_check_transform \
1577:           -dm_view ::ascii_info_detail -coarse_dm_view ::ascii_info_detail \
1578:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1579:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1580:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1582:     test:
1583:       suffix: hex_0
1584:     test:
1585:       suffix: hex_1
1586:       nsize: 2
1587:       args: -petscpartitioner_type simple

1589:   test:
1590:     suffix: hex_2
1591:     nsize: 4
1592:     args: -test_num 2 \
1593:           -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1594:             -dm_plex_transform_active fault -dm_plex_save_transform -dm_plex_check_transform \
1595:           -dm_view ::ascii_info_detail -coarse_dm_view ::ascii_info_detail \
1596:           -orientation_view -orientation_view_synchronized \
1597:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1598:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1599:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1601:   test:
1602:     suffix: hex_3
1603:     args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_plex_box_faces 2,1,2 -dm_plex_box_upper 2.,1.,2. \
1604:             -dm_plex_cohesive_label_fault0 37,40 \
1605:             -dm_plex_cohesive_label_fault1 26 \
1606:           -f0_dm_refine 1 -f0_dm_plex_transform_type cohesive_extrude \
1607:             -f0_dm_plex_transform_active fault0  -f0_coarse_dm_view ::ascii_info_detail \
1608:           -f1_dm_refine 1 -f1_dm_plex_transform_type cohesive_extrude \
1609:             -f1_dm_plex_transform_active fault1  -f1_coarse_dm_view ::ascii_info_detail \
1610:           -dm_plex_save_transform -dm_plex_check_transform \
1611:           -dm_view ::ascii_info_detail
1612:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1614:   test:
1615:     suffix: hex_4
1616:     args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_plex_box_faces 4,1,2 -dm_plex_box_upper 4.,1.,2. \
1617:             -dm_plex_cohesive_label_fault0 65,68 \
1618:             -dm_plex_cohesive_label_fault1 46 \
1619:           -f0_dm_refine 1 -f0_dm_plex_transform_type cohesive_extrude \
1620:             -f0_dm_plex_transform_active fault0  -f0_coarse_dm_view ::ascii_info_detail \
1621:           -f1_dm_refine 1 -f1_dm_plex_transform_type cohesive_extrude \
1622:             -f1_dm_plex_transform_active fault1  -f1_coarse_dm_view ::ascii_info_detail \
1623:           -dm_plex_save_transform -dm_plex_check_transform \
1624:           -dm_view ::ascii_info_detail
1625:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1627:   testset:
1628:     args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_plex_box_faces 2,1,1 -dm_plex_box_upper 2,1,1 -dm_plex_cohesive_label_fault 15 \
1629:           -dm_plex_option_phases coh_,ref_ \
1630:             -coh_dm_refine 1 -coh_dm_plex_transform_type cohesive_extrude \
1631:               -coh_dm_plex_transform_active fault \
1632:             -ref_dm_refine 1 -ref_dm_plex_transform_type refine_regular \
1633:           -dm_view ::ascii_info_detail \
1634:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1635:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1636:     filter: sed -e "s/_start//g" -e "s/f0_bd_u_neg//g" -e "s/f0_bd_u_pos//g" -e "s/f0_bd_l//g" -e "s/g0_bd_ul_neg//g" -e "s/g0_bd_ul_pos//g" -e "s/g0_bd_lu//g" -e "s~_ZL.*~~g"

1638:     test:
1639:       suffix: hex_0_ref

1641: TEST*/