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: PETSC_EXTERN char tri_3x3_cv[];
 23: char              tri_3x3_cv[] = "\
 24: 2 18 16 3 0\n\
 25: 0 1 5\n\
 26: 0 5 4\n\
 27: 1 2 6\n\
 28: 1 6 5\n\
 29: 2 3 7\n\
 30: 2 7 6\n\
 31: 4 5 9\n\
 32: 4 9 8\n\
 33: 5 6 10\n\
 34: 5 10 9\n\
 35: 6 7 11\n\
 36: 6 11 10\n\
 37: 8 9 13\n\
 38: 8 13 12\n\
 39: 9 10 14\n\
 40: 9 14 13\n\
 41: 10 11 15\n\
 42: 10 15 14\n\
 43: 0 0 0\n\
 44: 1 0 0\n\
 45: 2 0 0\n\
 46: 3 0 0\n\
 47: 0 1 0\n\
 48: 1 1 0\n\
 49: 2 1 0\n\
 50: 3 1 0\n\
 51: 0 2 0\n\
 52: 1 2 0\n\
 53: 2 2 0\n\
 54: 3 2 0\n\
 55: 0 3 0\n\
 56: 1 3 0\n\
 57: 2 3 0\n\
 58: 3 3 0";

 60: /* List of test meshes

 62: Test tri_0: triangle

 64:  4-10--5      8-16--7-14--4
 65:  |\  1 |      |\     \  1 |
 66:  | \   |      | \     \   |
 67:  6  8  9  ->  9 12  2  11 13
 68:  |   \ |      |   \     \ |
 69:  | 0  \|      | 0  \     \|
 70:  2--7--3      3-10--6-15--5

 72: Test tri_1: triangle, not tensor

 74:  4-10--5      8-10--7-16--4
 75:  |\  1 |      |\     \  1 |
 76:  | \   |      | \     \   |
 77:  6  8  9  -> 11 14  2  13 15
 78:  |   \ |      |   \     \ |
 79:  | 0  \|      | 0  \     \|
 80:  2--7--3      3-12--6--9--5

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

 84:            9
 85:           / \
 86:          /   \
 87:        17  2  16
 88:        /       \
 89:       /         \
 90:      8-----15----5
 91:       \         /|\
 92:        \       / | \
 93:        18  3  12 |  14
 94:          \   /   |   \
 95:           \ /    |    \
 96:            4  0 11  1  7
 97:             \    |    /
 98:              \   |   /
 99:              10  |  13
100:                \ | /
101:                 \|/
102:                  6
103:   becomes
104:            8
105:           / \
106:          /   \
107:         /     \
108:       25   2   24
109:       /         \
110:      /           \
111:    13-----18------9
112: 28  |     5    26/ \
113:    14----19----10   \
114:      \         /|   |\
115:       \       / |   | \
116:       21  3  20 |   |  23
117:         \   /   |   |   \
118:          \ /    |   |    \
119:           6  0 17 4 16 1  7
120:            \    |   |    /
121:             \   |   |   /
122:             15  |   |  22
123:               \ |   | /
124:                \|   |/
125:                12---11
126:                  27

128: Test tri_3: tri_2, in parallel

130:            6
131:           / \
132:          /   \
133:         /     \
134:       12   1   11
135:       /         \
136:      /           \
137:     5-----10------2
138:                    \
139:     5-----9-----3   2
140:      \         /|   |\
141:       \       / |   | \
142:       10  1  8  |   |  9
143:         \   /   |   |   \
144:          \ /    |   |    \
145:           2  0  7   7  0  4
146:            \    |   |    /
147:             \   |   |   /
148:              6  |   |  8
149:               \ |   | /
150:                \|   |/
151:                 4   3
152:   becomes
153:                  11
154:                 / \
155:                /   \
156:               /     \
157:             19   1   18
158:             /         \
159:            /           \
160:           8-----14------4
161:         22 \     3       |
162:             9------15    |\
163:                     \    | \
164:     9------14-----5  \  20 |
165:   20\    3     18/ \  \/   |
166:    10----15-----6   |  5   |
167:      \         /|   |  |   |\
168:       \       / |   |  |   | \
169:       17  1 16  |   |  |   |  17
170:         \   /   | 2 |  | 2 |   \
171:          \ /    |   |  |   |    \
172:           4  0  13 12  13  12 0 10
173:            \    |   |  |   |    /
174:             \   |   |  |   |   /
175:             11  |   |  |   |  16
176:               \ |   |  |   | /
177:                \|   |  |   |/
178:                 8---7  7---6
179:                  19      21

181: Test tri_4 to tri_8: 3x3 triangles on 3 processes

183: The fault is the horizontal line y = 1. The fault is oriented after distribution, so each process holds
184: only part of it. The tests use different partitions, some with overlap.

186:  +-----+-----+-----+
187:  |    /|    /|    /|
188:  |  /  |  /  |  /  |
189:  |/    |/    |/    |
190:  +-----+-----+-----+
191:  |    /|    /|    /|
192:  |  /  |  /  |  /  |
193:  |/    |/    |/    |
194:  +=====+=====+=====+
195:  |    /|    /|    /|
196:  |  /  |  /  |  /  |
197:  |/    |/    |/    |
198:  +-----+-----+-----+

200: Test quad_0: quadrilateral

202:  5-10--6-11--7       5-12-10-20--9-14--6
203:  |     |     |       |     |     |     |
204: 12  0 13  1  14 --> 15  0 18  2 17  1  16
205:  |     |     |       |     |     |     |
206:  2--8--3--9--4       3-11--8-19--7-13--4

208: Test quad_1: quadrilateral, not tensor

210:  5-10--6-11--7       5-14-10-12--9-16--6
211:  |     |     |       |     |     |     |
212: 12  0 13  1  14 --> 17  0 20  2 19  1  18
213:  |     |     |       |     |     |     |
214:  2--8--3--9--4       3-13--8-11--7-15--4

216: Test quad_2: quadrilateral, 2 processes

218:  3--6--4  3--6--4       3--9--7-14--6   5-14--4--9--7
219:  |     |  |     |       |     |     |   |     |     |
220:  7  0  8  7  0  8  --> 10  0 12  1 11  12  1 11  0  10
221:  |     |  |     |       |     |     |   |     |     |
222:  1--5--2  1--5--2       2--8--5-13--4   3-13--2--8--6

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

226:  3--6--4  3--6--4      3--9--7-14--6   5-14--4--9--7
227:  |     |  |     |      |     |     |   |     |     |
228:  7  0  8  7  0  8     10  0  12 1  11 12  1 11  0  10
229:  |     |  |     |      |     |     |   |     |     |
230:  1--5--2  1--5--2      2--8--5-13--4   3-13--2--8--6
231:                    -->
232:  3--6--4  3--6--4      3--9--7-14--6   5-14--4--9--7
233:  |     |  |     |      |     |     |   |     |     |
234:  7  0  8  7  0  8     10  0  12 1  11 12  1 11  0  10
235:  |     |  |     |      |     |     |   |     |     |
236:  1--5--2  1--5--2      2--8--5-13--4   3-13--2--8--6

238: Test quad_4: embedded fault

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

250: becomes

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

266: Test quad_5: two faults

268: 14-24-15-25-16-26--17
269:  |     |     |     |
270: 28  3 30  4 32  5  34
271:  |     |     |     |
272: 10-21-11-22-12-23--13
273:  |     |     |     |
274: 27  0 29  1 31  2  33
275:  |     |     |     |
276:  6-18--7-19--8-20--9

278: becomes

280: 12-26-13-27-14-28--15
281:  |     |     |     |
282: 37  4 31  3 33  5  40
283:  |     |     |     |
284: 17-36-18-25-19-39--21
285:  |     |     |     |
286: 43  6  44   41  7  42
287:  |     |     |     |
288: 16-35-18-25-19-38--20
289:  |     |     |     |
290: 29  0 30  1 32  2  34
291:  |     |     |     |
292:  8-22--9-23-10-24--11

294: Test quad_6: T-junction

296: 14-24-15-25-16-26--17
297:  |     |     |     |
298: 28  3 30  4 32  5  34
299:  |     |     |     |
300: 10-21-11-22-12-23--13
301:  |     |     |     |
302: 27  0 29  1 31  2  33
303:  |     |     |     |
304:  6-18--7-19--8-20--9

306: becomes

308:  13-26-14-27-15-28--16
309:   |     |     |     |
310:  30  3 32  4 39  5  40
311:   |     |     |     |
312:  12-25-17-36-19-38--21
313:         |     |     |
314:        41  6 42  7  43
315:         |     |     |
316:  12-25-17-35-18-37--20
317:   |     |     |     |
318:  29  0 31  1 33  2  34
319:   |     |     |     |
320:   8-22--9-23-10-24--11

322: becomes

324:  14-28-15-41-21-44--20-29-16
325:   |     |     |     |     |
326:  31  3 33  5 43  8 42  4  40
327:   |     |     |     |     |
328:  13-27-17-37-23-46--23-39-19
329:         |     |     |     |
330:        47  6 48    48  7  49
331:         |     |     |     |
332:  13-27-17-36-22-45--22-38-18
333:   |     |     |     |     |
334:  30  0 32  1 34    34  2  35
335:   |     |     |     |     |
336:   9-24-10-25-11-----11-26-12

338: Test tet_0: Two tets sharing a face

340:  cell   5 _______    cell
341:  0    / | \      \      1
342:     19  |  16     20
343:     /  15   \      \
344:    2-17------4--22--6
345:     \   |   /      /
346:     18  |  14     21
347:       \ | /      /
348:         3-------

350: becomes

352:  cell  10 ___36____9______    cell
353:  0    / | \        |\      \     1
354:     29  |  27      | 26     31
355:     /  25   \     24  \      \
356:    3-28------8--35-----7--33--4
357:     \   |   /      |  /      /
358:     30  |  23      | 22     32
359:       \ | /        |/      /
360:         6----34----5------
361:          cell 2

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

365:  cell   4          3______    cell
366:  0    / | \        |\      \     0
367:     14  |  11      | 11     12
368:     /  10   \     10  \      \
369:    1-12------3     |   2--14--4
370:     \   |   /      |  /      /
371:     13  |  9       | 9      13
372:       \ | /        |/      /
373:         2          1------

375: becomes
376:            cell 1              cell 1
377:  cell   8---28---7           7---28---6______    cell
378:  0    / | \      |\          |\       |\      \     0
379:     24  |  22    | 21        | 22     | 21     23
380:     /  20   \    |   \       |  \    19  \      \
381:    2-23------6---27---5     20  5---27---4--25--8
382:     \   |   /   19   /       |  /     |  /      /
383:     25  |  18    | 17        | 18     | 17     24
384:       \ | /      |/          |/       |/      /
385:         4---26---3           3---26---2------

387: Test hex_0: Two hexes sharing a face

389: cell  11-----31-----12-----32------13 cell
390: 0     /|            /|            /|     1
391:     36 |   22      37|   24      38|
392:     /  |          /  |          /  |
393:    8-----29------9-----30------10  |
394:    |   |     18  |   |     20  |   |
395:    |  42         |  43         |   44
396:    |14 |         |15 |         |16 |
397:   39   |  17    40   |   19   41   |
398:    |   5-----27--|---6-----28--|---7
399:    |  /          |  /          |  /
400:    | 33   21     | 34    23    | 35
401:    |/            |/            |/
402:    2-----25------3-----26------4

404: becomes

406:                          cell 2
407: cell   9-----38-----18-----62------17----42------10 cell
408: 0     /|            /|            /|            /|     1
409:     45 |   30      54|  32       53|   24      46|
410:     /  |          /  |          /  |          /  |
411:    7-----37-----16-----61------15--|-41------8   |
412:    |   |     28  |   |         |   |     22  |   |
413:    |  49         |  58         |   57        |   50
414:    |19 |         |26 |         |25 |         |20 |
415:   47   |  27    56   |        55   |   21   48   |
416:    |   5-----36--|--14-----60--|---13----40--|---6
417:    |  /          |  /          |  /          |  /
418:    | 43   29     | 52   31     | 51    23    | 44
419:    |/            |/            |/            |/
420:    3-----35-----12-----59------11----39------4

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

424: cell   7-----18------8             7-----18------8 cell
425: 0     /|            /|            /|            /|    0
426:     21 |   14      22|           21|   14      22|
427:     /  |          /  |          /  |          /  |
428:    5-----17------6   |         5---|-17------6   |
429:    |   |     12  |   |         |   |     12  |   |
430:    |  25         |  26         |  25         |  26
431:    | 9 |         |10 |         | 9 |         |10 |
432:   23   |  11    24   |        23   |   11   24   |
433:    |   3-----16--|---4         |   3-----16--|---4
434:    |  /          |  /          |  /          |  /
435:    | 19   13     | 20          | 19    13    | 20
436:    |/            |/            |/            |/
437:    1-----15------2             1-----15------2

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

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

458:           cell 0                  cell 0
459:        7-----18------8       7-----18------8
460:       /|            /|      /|            /|
461:     21 |   14      22|    21 |   14      22|
462:     /  |          /  |    /  |          /  |
463:    5-----17------6   |   5-----17------6   |
464:    |   |     12  |   |   |   |     12  |   |
465:    |  25         |  26   |  25         |   26
466:    |9  |         |10 |   |9  |         |10 |
467:   23   |  11    24   |  23   |  11    24   |
468:    |   3-----16--|---4   |   3-----16--|---4
469:    |  /          |  /    |  /          |  /
470:    | 19    13    | 20    | 19    13    | 20
471:    |/            |/      |/            |/
472:    1-----15------2       1-----15------2

474:        7-----18------8       7-----18------8
475:       /|            /|      /|            /|
476:     21 |   14      22|    21 |   14      22|
477:     /  |          /  |    /  |          /  |
478:    5-----17------6   |   5-----17------6   |
479:    |   |     12  |   |   |   |     12  |   |
480:    |  25         |  26   |  25         |  26
481:    |9  |         |10 |   |9  |         |10 |
482:   23   |  11    24   |  23   |   11   24   |
483:    |   3-----16--|---4   |   3-----16--|---4
484:    |  /          |  /    |  /          |  /
485:    | 19   13     | 20    | 19    13    | 20
486:    |/            |/      |/            |/
487:    1-----15------2       1-----15------2
488:       cell 0                cell 0

490: becomes

492:           cell 0         cell 1                cell 1        cell 0
493:        5-----28------13----44------12      9-----44------8-----28------13
494:       /|            /|            /|      /|            /|            /|
495:     30 |   20      36|   22      35|     36|   22     35 |   20      30|
496:     /  |          /  |          /  |    /  |          /  |          /  |
497:    4-----27------11----43------10  |   7-----43------6-----27------12  |
498:    |   |     18  |   |         |   |   |   |         |   |     18  |   |
499:    |  32         |  40         |  39   |  40         |  39         |   32
500:    |14 |         |16 |         |15 |   |15 |         |14 |         |16 |
501:   31   |  17    38   |         37  |   38  |        37   |  17    31   |
502:    |   3-----26--|---9-----42--|---8   |   5-----42--|---4-----26--|---11
503:    |  /          |  /          |  /    |  /          |  /          |  /
504:    | 29    19    | 34    21    |33     | 34    21    | 33    19    | 29
505:    |/            |/            |/      |/            |/            |/
506:    2-----25------7-----41------6       3-----41------2-----25------10

508:        5-----28------13----44------12      9-----44------8-----28------13
509:       /|            /|            /|      /|            /|            /|
510:     30 |   20      36|   22      35|     36|    22     35|   20      30|
511:     /  |          /  |          /  |    /  |          /  |          /  |
512:    4-----27------11----43------10  |   7-----43------6-----27------12  |
513:    |   |     18  |   |         |   |   |   |         |   |     18  |   |
514:    |  32         |  40         |   39  |   40        |  39         |   32
515:    |14 |         |16 |         |15 |   |15 |         |14 |         |16 |
516:   31   |  17    38   |         37  |   38  |        37   |  17    31   |
517:    |   3-----26--|---9-----42--|---8   |   5-----42--|---4-----26--|---11
518:    |  /          |  /          |  /    |  /          |  /          |  /
519:    | 29    19    | 34    21    |33     | 34    21    | 33    19    | 29
520:    |/            |/            |/      |/            |/            |/
521:    2-----25------7-----41------6       3-----41------2-----25------10
522:       cell 0         cell 1                cell 1        cell 0

524: Test hex_3: T-junction

526:       19-----52-----20-----53------21
527:       /|            /|            /|
528:     60 |   38      61|   41      62|
529:     /  |          /  |          /  |
530:   16-----50-----17-----51------18  |
531:    |   |     33  |   |     35  |   |
532:    |  70         |  72         |   74
533:    |25 |         |26 |         |27 |
534:   64   |  32    66   |  34    68   |
535:    |  13-----48--|--14-----49--|---15
536:    |  /|         |  /|         |  /|
537:    |57 |   37    | 58|   40    | 59|
538:    |/  |         |/  |         |/  |
539:   10-----46-----11-----47------12  |
540:    |   |     29  |   |     31  |   |
541:    |  69         |  71         |   73
542:    |22 |         |23 |         |24 |
543:   63   |  28    65   |   30   67   |
544:    |   7-----44--|---8-----45--|---9
545:    |  /          |  /          |  /
546:    | 54   36     | 55    39    | 56
547:    |/            |/            |/
548:    4-----42------5-----43------6
549:       cell 0         cell 1

551: becomes

553:       15----102-----28---112----___27-----73------16
554:       /|            /|         /   /             /|
555:     77 |   55     104|      ---  103    46      78|
556:     /  |          /  |     /     /             /  |
557:   13----101-----26---111--/----25-----72------14  |
558:    |   |     54  |   |  107   /           43  |   |
559:    |  81         |  108 / 51 /                |   82
560:    |40 |         |52 | /   105                |41 |
561:   79   |  53    106  |/   /            42    80   |
562:    |  21-----87--|--31---/-89------23-------/----/
563:    |  /|         |  /|  /         /|       /
564:    |91 |   47    |109|-- 49      93|  -----
565:    |/  |         |/ /|          /  | /
566:   17-----83-----29-----85------19----
567:    |   |         |   |         |   |
568:    |  120        |  121        |  122
569:    |   |         |26 |         |   |
570:  117   |        118  |        119  |
571:    |  22-----88--|--32-----90--|---24
572:    |  /|         |  /|         |  /|
573:    |92 |   48    |110|   50    | 94|
574:    |/  |         |/  |         |/  |
575:   18-----84-----30-----86------20  |
576:    |   |     37  |   |     39  |   |
577:    |  98         |  99         |   100
578:    |33 |         |34 |         |35 |
579:   95   |  36    96   |   38   97   |
580:    |  10-----70--|--11-----71--|---12
581:    |  /          |  /          |  /
582:    | 74   44     | 75    45    | 76
583:    |/            |/            |/
584:    7-----68------8-----69------9
585:       cell 0         cell 1

587: Test hex_4: Two non-intersecting faults

589:           cell 4         cell 5         cell 6        cell 7
590:       33-----96-----34-----97-----35-----98-----36-----99------37
591:       /|            /|            /|            /|            /|
592:     110|   66     111|   69     112|   72     113|   75     114|
593:     /  |          /  |          /  |          /  |          /  |
594:   28-----92-----29-----93-----30-----94-----31-----95------32  |
595:    |   |     57  |   |     59  |   |     61  |   |     63  |   |
596:    |  126        |  128        |  130        |  132        |  134
597:    |43 |         |44 |         |45 |         |46 |         |47 |
598:   116  |  56    118  |  58    120  |  60    122  |  62    124  |
599:    |  23-----88--|--24-----89--|--25-----90--|--26-----91--|---27
600:    |  /|         |  /|         |  /|         |  /|         |  /|
601:    |105|   65    |106|   68    |107|   71    |108|   74    |109|
602:    |/  |         |/  |         |/  |         |/  |         |/  |
603:   18-----84-----19-----95-----20-----86-----21-----87------22  |
604:    |   |     49  |   |     51  |   |     53  |   |     55  |   |
605:    |  125        |  127        |  129        |  131        |  133
606:    |38 |         |39 |         |40 |         |41 |         |42 |
607:   115  |  48    117  |  50    119  |  52    121  |  54    123  |
608:    |  13-----80--|--14-----81--|--15-----82--|--16-----83--|---17
609:    |  /          |  /          |  /          |  /          |  /
610:    |100    64    |101    67    |102    70    |103    73    |104
611:    |/            |/            |/            |/            |/
612:    8-----76------9-----77-----10-----78-----11-----79------12
613:       cell 0         cell 1        cell 2        cell 3

615: becomes

617:           cell 4         cell 5        cell 7        cell 10       cell 6
618:       27-----114----28-----115----29-----159----46-----170----45------116----30
619:       /|            /|            /|            /|            /|            /|
620:     123|   71     124|   73     125|   87     162|          161|    78    126|
621:     /  |          /  |          /  |          /  |          /  |          /  |
622:   23-----111----24-----112----25-----158----44-----169----43-----113-----26  |
623:    |   |     65  |   |    67   |   |    86   |   |         |   |     69  |   |
624:    |  134        |  135        |  137        |  166        |  165        |  140
625:    |56 |         |57 |         |58 |         |84 |         |83 |         |59 |
626:   127  |  64    128  |  66    130  |  85    164  |        163  |  68    133  |
627:    |  35-----143-|--37-----151-|--40-----109-|--42-----168-|--42-----110-|---22
628:    |  /|         |  /|         |  /|         |  /          |  /          |  /
629:    |145|   79    |147|   81    |153|   75    |160          |160    77    |122
630:    |/ 173        |/ 174        |/ 176        |/            |/            |/
631:   31-----141----33-----149----39-----107----41-----167----41-----108-----21
632: cell   |         |   |         |   | cell 9
633: 8  |  36-----144-|--38-----152-|--40-----109----42-----110-----22
634:   171 /|        172 /|        175 /|            /|            /|
635:    |146|   80    |148|   82    |153|    75    160|   77     122|
636:    |/  |         |/  |         |/  |          /  |          /  |
637:   32-----142----34-----150----39-----107----41-----108-----21  |
638:    |   |     50  |   |    52   |   |    61   |   |     63  |   |
639:    |  156        |  157        |  136        |  138        |  139
640:    |47 |         |48 |         |53 |         |54 |         |55 |
641:   154  |  49    155  |  51    129  |  60    131  |  62    132  |
642:    |  16-----103-|--17-----104-|--18-----105-|--19-----106-|---20
643:    |  /          |  /          |  /          |  /          |  /
644:    |117    70    |118    72    |119    74    |120    76    |121
645:    |/            |/            |/            |/            |/
646:   11-----99-----12-----100----13-----101----14-----102-----15
647:       cell 0         cell 1        cell 2        cell 3

649: */

651: typedef struct {
652:   PetscInt testNum;        // The mesh to test
653:   PetscInt cohesiveFields; // The number of fault fields
654: } AppCtx;

656: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
657: {
658:   PetscFunctionBegin;
659:   options->testNum        = 0;
660:   options->cohesiveFields = 1;

662:   PetscOptionsBegin(comm, "", "Cohesive Meshing Options", "DMPLEX");
663:   PetscCall(PetscOptionsBoundedInt("-test_num", "The particular mesh to test", __FILE__, options->testNum, &options->testNum, NULL, 0));
664:   PetscCall(PetscOptionsBoundedInt("-cohesive_fields", "The number of cohesive fields", __FILE__, options->cohesiveFields, &options->cohesiveFields, NULL, 0));
665:   PetscOptionsEnd();
666:   PetscFunctionReturn(PETSC_SUCCESS);
667: }

669: static PetscErrorCode CreateQuadMesh1(MPI_Comm comm, AppCtx *user, DM *dm)
670: {
671:   const PetscInt faces[2] = {1, 1};
672:   PetscReal      lower[2], upper[2];
673:   DMLabel        label;
674:   PetscMPIInt    rank;
675:   void          *get_tmp;
676:   PetscInt64    *cidx;
677:   PetscMPIInt    iflg;

679:   PetscFunctionBeginUser;
680:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
681:   // Create serial mesh
682:   lower[0] = (PetscReal)(rank % 2);
683:   lower[1] = (PetscReal)(rank / 2);
684:   upper[0] = (PetscReal)(rank % 2) + 1.;
685:   upper[1] = (PetscReal)(rank / 2) + 1.;
686:   PetscCall(DMPlexCreateBoxMesh(PETSC_COMM_SELF, 2, PETSC_FALSE, faces, lower, upper, NULL, PETSC_TRUE, 0, PETSC_TRUE, dm));
687:   PetscCall(PetscObjectSetName((PetscObject)*dm, "box"));
688:   // Flip edges to make fault non-oriented
689:   switch (rank) {
690:   case 2:
691:     PetscCall(DMPlexOrientPoint(*dm, 8, -1));
692:     break;
693:   case 3:
694:     PetscCall(DMPlexOrientPoint(*dm, 7, -1));
695:     break;
696:   default:
697:     break;
698:   }
699:   // Need this so that all procs create the cell types
700:   PetscCall(DMPlexGetCellTypeLabel(*dm, &label));
701:   // Replace comm in object (copied from PetscHeaderCreate/Destroy())
702:   PetscCall(PetscCommDestroy(&(*dm)->hdr.comm));
703:   PetscCall(PetscCommDuplicate(comm, &(*dm)->hdr.comm, &(*dm)->hdr.tag));
704:   PetscCallMPI(MPI_Comm_get_attr((*dm)->hdr.comm, Petsc_CreationIdx_keyval, &get_tmp, &iflg));
705:   PetscCheck(iflg, (*dm)->hdr.comm, PETSC_ERR_ARG_CORRUPT, "MPI_Comm does not have an object creation index");
706:   cidx            = (PetscInt64 *)get_tmp;
707:   (*dm)->hdr.cidx = (*cidx)++;
708:   // Create new pointSF
709:   {
710:     PetscSF      sf;
711:     PetscInt    *local  = NULL;
712:     PetscSFNode *remote = NULL;
713:     PetscInt     Nl;

715:     PetscCall(PetscSFCreate(comm, &sf));
716:     switch (rank) {
717:     case 0:
718:       Nl = 5;
719:       PetscCall(PetscMalloc1(Nl, &local));
720:       PetscCall(PetscMalloc1(Nl, &remote));
721:       local[0]        = 2;
722:       remote[0].index = 1;
723:       remote[0].rank  = 1;
724:       local[1]        = 3;
725:       remote[1].index = 1;
726:       remote[1].rank  = 2;
727:       local[2]        = 4;
728:       remote[2].index = 1;
729:       remote[2].rank  = 3;
730:       local[3]        = 6;
731:       remote[3].index = 5;
732:       remote[3].rank  = 2;
733:       local[4]        = 8;
734:       remote[4].index = 7;
735:       remote[4].rank  = 1;
736:       break;
737:     case 1:
738:       Nl = 3;
739:       PetscCall(PetscMalloc1(Nl, &local));
740:       PetscCall(PetscMalloc1(Nl, &remote));
741:       local[0]        = 3;
742:       remote[0].index = 1;
743:       remote[0].rank  = 3;
744:       local[1]        = 4;
745:       remote[1].index = 2;
746:       remote[1].rank  = 3;
747:       local[2]        = 6;
748:       remote[2].index = 5;
749:       remote[2].rank  = 3;
750:       break;
751:     case 2:
752:       Nl = 3;
753:       PetscCall(PetscMalloc1(Nl, &local));
754:       PetscCall(PetscMalloc1(Nl, &remote));
755:       local[0]        = 2;
756:       remote[0].index = 1;
757:       remote[0].rank  = 3;
758:       local[1]        = 4;
759:       remote[1].index = 3;
760:       remote[1].rank  = 3;
761:       local[2]        = 8;
762:       remote[2].index = 7;
763:       remote[2].rank  = 3;
764:       break;
765:     case 3:
766:       Nl = 0;
767:       break;
768:     default:
769:       SETERRQ(comm, PETSC_ERR_SUP, "This example only supports 4 ranks");
770:     }
771:     PetscCall(PetscSFSetGraph(sf, 9, Nl, local, PETSC_OWN_POINTER, remote, PETSC_OWN_POINTER));
772:     PetscCall(DMSetPointSF(*dm, sf));
773:     PetscCall(PetscSFDestroy(&sf));
774:   }
775:   // Create fault label
776:   PetscCall(DMCreateLabel(*dm, "fault"));
777:   PetscCall(DMGetLabel(*dm, "fault", &label));
778:   switch (rank) {
779:   case 0:
780:   case 2:
781:     PetscCall(DMLabelSetValue(label, 8, 1));
782:     PetscCall(DMLabelSetValue(label, 2, 0));
783:     PetscCall(DMLabelSetValue(label, 4, 0));
784:     break;
785:   case 1:
786:   case 3:
787:     PetscCall(DMLabelSetValue(label, 7, 1));
788:     PetscCall(DMLabelSetValue(label, 1, 0));
789:     PetscCall(DMLabelSetValue(label, 3, 0));
790:     break;
791:   default:
792:     break;
793:   }
794:   PetscCall(DMPlexOrientLabel(*dm, label));
795:   PetscCall(DMPlexLabelCohesiveComplete(*dm, label, NULL, 1, PETSC_FALSE, NULL));
796:   PetscCall(DMPlexDistributeSetDefault(*dm, PETSC_FALSE));
797:   PetscFunctionReturn(PETSC_SUCCESS);
798: }

800: static PetscErrorCode CreateHexMesh1(MPI_Comm comm, AppCtx *user, DM *dm)
801: {
802:   const PetscInt faces[3] = {1, 1, 1};
803:   PetscReal      lower[3], upper[3];
804:   DMLabel        label;
805:   PetscMPIInt    rank;
806:   void          *get_tmp;
807:   PetscInt64    *cidx;
808:   PetscMPIInt    iflg;

810:   PetscFunctionBeginUser;
811:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
812:   // Create serial mesh
813:   lower[0] = (PetscReal)(rank % 2);
814:   lower[1] = 0.;
815:   lower[2] = (PetscReal)(rank / 2);
816:   upper[0] = (PetscReal)(rank % 2) + 1.;
817:   upper[1] = 1.;
818:   upper[2] = (PetscReal)(rank / 2) + 1.;
819:   PetscCall(DMPlexCreateBoxMesh(PETSC_COMM_SELF, 3, PETSC_FALSE, faces, lower, upper, NULL, PETSC_TRUE, 0, PETSC_TRUE, dm));
820:   PetscCall(PetscObjectSetName((PetscObject)*dm, "box"));
821:   // Flip edges to make fault non-oriented
822:   switch (rank) {
823:   case 2:
824:     PetscCall(DMPlexOrientPoint(*dm, 10, -1));
825:     break;
826:   case 3:
827:     PetscCall(DMPlexOrientPoint(*dm, 9, -1));
828:     break;
829:   default:
830:     break;
831:   }
832:   // Need this so that all procs create the cell types
833:   PetscCall(DMPlexGetCellTypeLabel(*dm, &label));
834:   // Replace comm in object (copied from PetscHeaderCreate/Destroy())
835:   PetscCall(PetscCommDestroy(&(*dm)->hdr.comm));
836:   PetscCall(PetscCommDuplicate(comm, &(*dm)->hdr.comm, &(*dm)->hdr.tag));
837:   PetscCallMPI(MPI_Comm_get_attr((*dm)->hdr.comm, Petsc_CreationIdx_keyval, &get_tmp, &iflg));
838:   PetscCheck(iflg, (*dm)->hdr.comm, PETSC_ERR_ARG_CORRUPT, "MPI_Comm does not have an object creation index");
839:   cidx            = (PetscInt64 *)get_tmp;
840:   (*dm)->hdr.cidx = (*cidx)++;
841:   // Create new pointSF
842:   {
843:     PetscSF      sf;
844:     PetscInt    *local  = NULL;
845:     PetscSFNode *remote = NULL;
846:     PetscInt     Nl;

848:     PetscCall(PetscSFCreate(comm, &sf));
849:     switch (rank) {
850:     case 0:
851:       Nl = 15;
852:       PetscCall(PetscMalloc1(Nl, &local));
853:       PetscCall(PetscMalloc1(Nl, &remote));
854:       local[0]         = 2;
855:       remote[0].index  = 1;
856:       remote[0].rank   = 1;
857:       local[1]         = 4;
858:       remote[1].index  = 3;
859:       remote[1].rank   = 1;
860:       local[2]         = 5;
861:       remote[2].index  = 1;
862:       remote[2].rank   = 2;
863:       local[3]         = 6;
864:       remote[3].index  = 1;
865:       remote[3].rank   = 3;
866:       local[4]         = 7;
867:       remote[4].index  = 3;
868:       remote[4].rank   = 2;
869:       local[5]         = 8;
870:       remote[5].index  = 3;
871:       remote[5].rank   = 3;
872:       local[6]         = 17;
873:       remote[6].index  = 15;
874:       remote[6].rank   = 2;
875:       local[7]         = 18;
876:       remote[7].index  = 16;
877:       remote[7].rank   = 2;
878:       local[8]         = 20;
879:       remote[8].index  = 19;
880:       remote[8].rank   = 1;
881:       local[9]         = 21;
882:       remote[9].index  = 19;
883:       remote[9].rank   = 2;
884:       local[10]        = 22;
885:       remote[10].index = 19;
886:       remote[10].rank  = 3;
887:       local[11]        = 24;
888:       remote[11].index = 23;
889:       remote[11].rank  = 1;
890:       local[12]        = 26;
891:       remote[12].index = 25;
892:       remote[12].rank  = 1;
893:       local[13]        = 10;
894:       remote[13].index = 9;
895:       remote[13].rank  = 1;
896:       local[14]        = 14;
897:       remote[14].index = 13;
898:       remote[14].rank  = 2;
899:       break;
900:     case 1:
901:       Nl = 9;
902:       PetscCall(PetscMalloc1(Nl, &local));
903:       PetscCall(PetscMalloc1(Nl, &remote));
904:       local[0]        = 5;
905:       remote[0].index = 1;
906:       remote[0].rank  = 3;
907:       local[1]        = 6;
908:       remote[1].index = 2;
909:       remote[1].rank  = 3;
910:       local[2]        = 7;
911:       remote[2].index = 3;
912:       remote[2].rank  = 3;
913:       local[3]        = 8;
914:       remote[3].index = 4;
915:       remote[3].rank  = 3;
916:       local[4]        = 17;
917:       remote[4].index = 15;
918:       remote[4].rank  = 3;
919:       local[5]        = 18;
920:       remote[5].index = 16;
921:       remote[5].rank  = 3;
922:       local[6]        = 21;
923:       remote[6].index = 19;
924:       remote[6].rank  = 3;
925:       local[7]        = 22;
926:       remote[7].index = 20;
927:       remote[7].rank  = 3;
928:       local[8]        = 14;
929:       remote[8].index = 13;
930:       remote[8].rank  = 3;
931:       break;
932:     case 2:
933:       Nl = 9;
934:       PetscCall(PetscMalloc1(Nl, &local));
935:       PetscCall(PetscMalloc1(Nl, &remote));
936:       local[0]        = 2;
937:       remote[0].index = 1;
938:       remote[0].rank  = 3;
939:       local[1]        = 4;
940:       remote[1].index = 3;
941:       remote[1].rank  = 3;
942:       local[2]        = 6;
943:       remote[2].index = 5;
944:       remote[2].rank  = 3;
945:       local[3]        = 8;
946:       remote[3].index = 7;
947:       remote[3].rank  = 3;
948:       local[4]        = 20;
949:       remote[4].index = 19;
950:       remote[4].rank  = 3;
951:       local[5]        = 22;
952:       remote[5].index = 21;
953:       remote[5].rank  = 3;
954:       local[6]        = 24;
955:       remote[6].index = 23;
956:       remote[6].rank  = 3;
957:       local[7]        = 26;
958:       remote[7].index = 25;
959:       remote[7].rank  = 3;
960:       local[8]        = 10;
961:       remote[8].index = 9;
962:       remote[8].rank  = 3;
963:       break;
964:     case 3:
965:       Nl = 0;
966:       break;
967:     default:
968:       SETERRQ(comm, PETSC_ERR_SUP, "This example only supports 4 ranks");
969:     }
970:     PetscCall(PetscSFSetGraph(sf, 27, Nl, local, PETSC_OWN_POINTER, remote, PETSC_OWN_POINTER));
971:     PetscCall(DMSetPointSF(*dm, sf));
972:     PetscCall(PetscSFDestroy(&sf));
973:   }
974:   // Create fault label
975:   PetscCall(DMCreateLabel(*dm, "fault"));
976:   PetscCall(DMGetLabel(*dm, "fault", &label));
977:   switch (rank) {
978:   case 0:
979:   case 2:
980:     PetscCall(DMLabelSetValue(label, 10, 2));
981:     PetscCall(DMLabelSetValue(label, 20, 1));
982:     PetscCall(DMLabelSetValue(label, 22, 1));
983:     PetscCall(DMLabelSetValue(label, 24, 1));
984:     PetscCall(DMLabelSetValue(label, 26, 1));
985:     PetscCall(DMLabelSetValue(label, 2, 0));
986:     PetscCall(DMLabelSetValue(label, 4, 0));
987:     PetscCall(DMLabelSetValue(label, 6, 0));
988:     PetscCall(DMLabelSetValue(label, 8, 0));
989:     break;
990:   case 1:
991:   case 3:
992:     PetscCall(DMLabelSetValue(label, 9, 2));
993:     PetscCall(DMLabelSetValue(label, 19, 1));
994:     PetscCall(DMLabelSetValue(label, 21, 1));
995:     PetscCall(DMLabelSetValue(label, 23, 1));
996:     PetscCall(DMLabelSetValue(label, 25, 1));
997:     PetscCall(DMLabelSetValue(label, 1, 0));
998:     PetscCall(DMLabelSetValue(label, 3, 0));
999:     PetscCall(DMLabelSetValue(label, 5, 0));
1000:     PetscCall(DMLabelSetValue(label, 7, 0));
1001:     break;
1002:   default:
1003:     break;
1004:   }
1005:   PetscCall(DMPlexOrientLabel(*dm, label));
1006:   PetscCall(DMPlexLabelCohesiveComplete(*dm, label, NULL, 1, PETSC_FALSE, NULL));
1007:   PetscCall(DMPlexDistributeSetDefault(*dm, PETSC_FALSE));
1008:   PetscFunctionReturn(PETSC_SUCCESS);
1009: }

1011: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
1012: {
1013:   PetscFunctionBegin;
1014:   switch (user->testNum) {
1015:   case 1:
1016:     PetscCall(CreateQuadMesh1(comm, user, dm));
1017:     break;
1018:   case 2:
1019:     PetscCall(CreateHexMesh1(comm, user, dm));
1020:     break;
1021:   default:
1022:     PetscCall(DMCreate(comm, dm));
1023:     PetscCall(DMSetType(*dm, DMPLEX));
1024:     break;
1025:   }
1026:   PetscCall(DMSetFromOptions(*dm));
1027:   {
1028:     const char *prefix;

1030:     // We cannot redistribute with cohesive cells in the SF
1031:     PetscCall(DMPlexDistributeSetDefault(*dm, PETSC_FALSE));
1032:     PetscCall(PetscObjectGetOptionsPrefix((PetscObject)*dm, &prefix));
1033:     PetscCall(PetscObjectSetOptionsPrefix((PetscObject)*dm, "f0_"));
1034:     PetscCall(DMSetFromOptions(*dm));
1035:     PetscCall(PetscObjectSetOptionsPrefix((PetscObject)*dm, "f1_"));
1036:     PetscCall(DMSetFromOptions(*dm));
1037:     PetscCall(PetscObjectSetOptionsPrefix((PetscObject)*dm, prefix));
1038:   }
1039:   PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
1040:   PetscFunctionReturn(PETSC_SUCCESS);
1041: }

1043: // Create a displacement field, and some number of vector fault fields
1044: static PetscErrorCode CreateDiscretization(DM dm, AppCtx *user)
1045: {
1046:   PetscSection   s;
1047:   DMLabel        fault, faultSpace;
1048:   PetscFE        fe;
1049:   DMPolytopeType ct, fct;
1050:   PetscInt       dim, cStart, fStart, Ncf = user->cohesiveFields;

1052:   PetscFunctionBegin;
1053:   PetscCall(DMGetDimension(dm, &dim));
1054:   PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, NULL));
1055:   PetscCall(DMPlexGetCellType(dm, cStart, &ct));
1056:   PetscCall(DMGetLabel(dm, "fault", &fault));
1057:   if (!fault) PetscFunctionReturn(PETSC_SUCCESS);
1058:   PetscCall(DMLabelView(fault, PETSC_VIEWER_STDOUT_WORLD));

1060:   PetscCall(PetscFECreateByCell(PETSC_COMM_SELF, dim, dim, ct, "displacement_", PETSC_DETERMINE, &fe));
1061:   PetscCall(PetscFESetName(fe, "displacement"));
1062:   PetscCall(DMAddField(dm, NULL, (PetscObject)fe));
1063:   PetscCall(PetscFEDestroy(&fe));

1065:   // Make label for fault space definition
1066:   PetscCall(DMCreateLabel(dm, "faultSpace"));
1067:   PetscCall(DMGetLabel(dm, "faultSpace", &faultSpace));
1068:   for (PetscInt d = 0; d <= dim; ++d) {
1069:     PetscInt pStart, pEnd, pMax;

1071:     PetscCall(DMPlexGetSimplexOrBoxCells(dm, d, NULL, &pMax));
1072:     PetscCall(DMPlexGetHeightStratum(dm, d, &pStart, &pEnd));
1073:     for (PetscInt p = pMax; p < pEnd; ++p) PetscCall(DMLabelSetValue(faultSpace, p, 1));
1074:   }
1075:   PetscCall(DMPlexGetHeightStratum(dm, 1, &fStart, NULL));
1076:   PetscCall(DMPlexGetCellType(dm, fStart, &fct));
1077:   if (Ncf > 0) {
1078:     PetscCall(PetscFECreateByCell(PETSC_COMM_SELF, dim - 1, dim, fct, "faulttraction_", PETSC_DETERMINE, &fe));
1079:     PetscCall(PetscFESetName(fe, "fault traction"));
1080:     PetscCall(DMAddField(dm, faultSpace, (PetscObject)fe));
1081:     PetscCall(PetscFEDestroy(&fe));
1082:   }
1083:   for (PetscInt f = 1; f < Ncf; ++f) {
1084:     char name[256], opt[256];

1086:     PetscCall(PetscSNPrintf(name, 256, "fault field %" PetscInt_FMT, f));
1087:     PetscCall(PetscSNPrintf(opt, 256, "faultfield_%" PetscInt_FMT "_", f));
1088:     PetscCall(PetscFECreateByCell(PETSC_COMM_SELF, dim - 1, dim, fct, opt, PETSC_DETERMINE, &fe));
1089:     PetscCall(PetscFESetName(fe, name));
1090:     PetscCall(DMAddField(dm, faultSpace, (PetscObject)fe));
1091:     PetscCall(PetscFEDestroy(&fe));
1092:   }
1093:   PetscCall(DMCreateDS(dm));

1095:   PetscCall(DMGetLocalSection(dm, &s));
1096:   PetscCall(PetscObjectViewFromOptions((PetscObject)s, NULL, "-local_section_view"));
1097:   PetscFunctionReturn(PETSC_SUCCESS);
1098: }

1100: // Label cells 1 for negative side, and 2 for positive side
1101: static PetscErrorCode CreateMaterialLabel(DM dm)
1102: {
1103:   DMLabel         fault, material;
1104:   IS              faceIS;
1105:   const PetscInt *faces;
1106:   PetscReal       fvol, fcentroid[3] = {0., 0., 0.}, fnormal[3] = {0., 0., 0.};
1107:   PetscInt        dim, cStart, cEnd, Nf = 0;
1108:   PetscMPIInt     rank, size, root;

1110:   PetscFunctionBegin;
1111:   PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
1112:   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)dm), &size));
1113:   PetscCall(DMGetDimension(dm, &dim));
1114:   PetscCall(DMGetLabel(dm, "fault", &fault));
1115:   PetscCall(DMCreateLabel(dm, "material"));
1116:   PetscCall(DMGetLabel(dm, "material", &material));
1117:   for (PetscInt s = 1; s < 3; ++s) {
1118:     IS              pointIS;
1119:     const PetscInt *points;
1120:     PetscInt        n;

1122:     PetscCall(DMLabelGetStratumIS(fault, s > 1 ? 100 + dim : -(100 + dim), &pointIS));
1123:     if (!pointIS) continue;
1124:     PetscCall(ISGetLocalSize(pointIS, &n));
1125:     PetscCall(ISGetIndices(pointIS, &points));
1126:     for (PetscInt i = 0; i < n; ++i) {
1127:       PetscCall(DMLabelSetValue(material, points[i], s));
1128:     }
1129:     PetscCall(ISRestoreIndices(pointIS, &points));
1130:     PetscCall(ISDestroy(&pointIS));
1131:   }
1132:   // This simple algorithm will work for now (note that cohesive cells get added into this label)
1133:   // A process can hold only vertices of the fault, so all processes use a face from the lowest process that has one
1134:   PetscCall(DMLabelGetStratumIS(fault, dim - 1, &faceIS));
1135:   if (faceIS) PetscCall(ISGetLocalSize(faceIS, &Nf));
1136:   root = Nf ? rank : size;
1137:   PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, &root, 1, MPI_INT, MPI_MIN, PetscObjectComm((PetscObject)dm)));
1138:   PetscCheck(root < size, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Fault label must contain at least one face");
1139:   if (rank == root) {
1140:     PetscCall(ISGetIndices(faceIS, &faces));
1141:     for (PetscInt i = 0; i < Nf; ++i) {
1142:       const PetscInt face = faces[i];
1143:       DMPolytopeType ct;

1145:       PetscCall(DMPlexGetCellType(dm, face, &ct));
1146:       if (DMPolytopeTypeGetDim(ct) != dim - 1) continue;
1147:       PetscCall(DMPlexComputeCellGeometryFVM(dm, face, &fvol, fcentroid, fnormal));
1148:       break;
1149:     }
1150:     PetscCall(ISRestoreIndices(faceIS, &faces));
1151:   }
1152:   PetscCall(ISDestroy(&faceIS));
1153:   PetscCallMPI(MPI_Bcast(fcentroid, 3, MPIU_REAL, root, PetscObjectComm((PetscObject)dm)));
1154:   PetscCallMPI(MPI_Bcast(fnormal, 3, MPIU_REAL, root, PetscObjectComm((PetscObject)dm)));
1155:   PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd));
1156:   for (PetscInt c = cStart; c < cEnd; ++c) {
1157:     PetscReal vol, centroid[3];
1158:     PetscInt  val;

1160:     PetscCall(DMLabelGetValue(fault, c, &val));
1161:     if (val >= 0) continue;
1162:     PetscCall(DMPlexComputeCellGeometryFVM(dm, c, &vol, centroid, NULL));
1163:     for (PetscInt e = 0; e < dim; ++e) centroid[e] -= fcentroid[e];
1164:     if (DMPlex_DotRealD_Internal(dim, centroid, fnormal) > 0) PetscCall(DMLabelSetValue(material, c, 1));
1165:     else PetscCall(DMLabelSetValue(material, c, 2));
1166:   }
1167:   PetscFunctionReturn(PETSC_SUCCESS);
1168: }

1170: // Label cohesive cells and endcap faces 1
1171: static PetscErrorCode CreateFaultLabel(DM dm)
1172: {
1173:   DMLabel  fault;
1174:   PetscInt cMax, cEnd;

1176:   PetscFunctionBegin;
1177:   PetscCall(DMCreateLabel(dm, "faultCells"));
1178:   PetscCall(DMGetLabel(dm, "faultCells", &fault));
1179:   PetscCall(DMPlexGetSimplexOrBoxCells(dm, 0, &cEnd, &cMax));
1180:   for (PetscInt c = cMax; c < cEnd; ++c) {
1181:     const PetscInt *cone;

1183:     PetscCall(DMLabelSetValue(fault, c, 1));
1184:     PetscCall(DMPlexGetCone(dm, c, &cone));
1185:     PetscCall(DMLabelSetValue(fault, cone[0], 1));
1186:     PetscCall(DMLabelSetValue(fault, cone[1], 1));
1187:   }
1188:   PetscFunctionReturn(PETSC_SUCCESS);
1189: }

1191: static PetscErrorCode r(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, PetscCtx ctx)
1192: {
1193:   PetscInt d;
1194:   for (d = 0; d < dim; ++d) u[d] = x[d];
1195:   return PETSC_SUCCESS;
1196: }

1198: static PetscErrorCode rp1(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, PetscCtx ctx)
1199: {
1200:   PetscInt d;
1201:   for (d = 0; d < dim; ++d) u[d] = x[d] + (d > 0 ? 1.0 : 0.0);
1202:   return PETSC_SUCCESS;
1203: }

1205: static PetscErrorCode phi(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, PetscCtx ctx)
1206: {
1207:   PetscInt d;
1208:   u[0] = -x[1];
1209:   u[1] = x[0];
1210:   for (d = 2; d < dim; ++d) u[d] = x[d];
1211:   return PETSC_SUCCESS;
1212: }

1214: 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[])
1215: {
1216:   PetscInt       d;
1217:   const PetscInt offN = 0;
1218:   const PetscInt offP = dim;
1219:   for (d = 0; d < dim; ++d) f[d] = u[offN + d] + u[offP + d];
1220: }

1222: 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[])
1223: {
1224:   PetscInt d;
1225:   for (d = 0; d < dim; ++d) f[d] = n[d];
1226: }

1228: /* \lambda \cdot (\psi_u^- - \psi_u^+) */
1229: 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[])
1230: {
1231:   const PetscInt Nc = dim + 1;
1232:   for (PetscInt c = 0; c < Nc; ++c) f0[c] = -u[uOff[1] + c];
1233: }

1235: 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[])
1236: {
1237:   const PetscInt Nc = dim + 1;
1238:   for (PetscInt c = 0; c < Nc; ++c) f0[c] = u[uOff[1] + c];
1239: }

1241: /* (d - u^+ + u^-) \cdot \psi_\lambda */
1242: 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[])
1243: {
1244:   const PetscInt Nc = uOff[2] - uOff[1];

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

1249: /* \psi_lambda \cdot (\psi_u^- - \psi_u^+) */
1250: 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[])
1251: {
1252:   const PetscInt Nc = dim + 1;
1253:   for (PetscInt c = 0; c < Nc; ++c) g0[c * Nc + c] = -1.0;
1254: }

1256: 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[])
1257: {
1258:   const PetscInt Nc = dim + 1;
1259:   for (PetscInt c = 0; c < Nc; ++c) g0[c * Nc + c] = 1.0;
1260: }

1262: /* (-\psi_u^+ + \psi_u^-) \cdot \psi_\lambda */
1263: 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[])
1264: {
1265:   const PetscInt Nc = uOff[2] - uOff[1];

1267:   for (PetscInt c = 0; c < Nc; ++c) {
1268:     g0[c * Nc + c]           = -1.0;
1269:     g0[Nc * Nc + c * Nc + c] = 1.0;
1270:   }
1271: }

1273: static PetscErrorCode TestAssembly(DM dm, AppCtx *user)
1274: {
1275:   Mat           J;
1276:   Vec           locX, locF, locW;
1277:   PetscDS       probh;
1278:   DMLabel       fault, material;
1279:   DM            dmFault;
1280:   IS            cohesiveCells;
1281:   PetscFE       fe;
1282:   PetscWeakForm wf;
1283:   PetscFormKey  keys[3];
1284:   PetscErrorCode (*initialGuess[2])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], PetscCtx ctx);
1285:   DMPolytopeType fct;
1286:   PetscInt       dim, fStart, Nf, cMax, cEnd, id;
1287:   PetscBool      hasCohesive;
1288:   PetscMPIInt    rank, size;

1290:   PetscFunctionBegin;
1291:   PetscCall(DMGetNumFields(dm, &Nf));
1292:   if (Nf <= 0) PetscFunctionReturn(PETSC_SUCCESS);
1293:   PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
1294:   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)dm), &size));
1295:   PetscCall(DMGetDimension(dm, &dim));
1296:   PetscCall(DMPlexGetSimplexOrBoxCells(dm, 0, NULL, &cMax));
1297:   PetscCall(DMPlexGetHeightStratum(dm, 0, NULL, &cEnd));
1298:   if (size > 1) {
1299:     PetscSF         sf;
1300:     const PetscInt *leaves;
1301:     PetscInt       *points;
1302:     PetscInt        Nl, l, Ncoh = 0;

1304:     PetscCall(DMGetPointSF(dm, &sf));
1305:     PetscCall(PetscSFGetGraph(sf, NULL, &Nl, &leaves, NULL));
1306:     for (PetscInt c = cMax; c < cEnd; ++c) {
1307:       PetscCall(PetscFindInt(c, Nl, leaves, &l));
1308:       if (l < 0) ++Ncoh;
1309:     }
1310:     PetscCall(PetscMalloc1(Ncoh, &points));
1311:     Ncoh = 0;
1312:     for (PetscInt c = cMax; c < cEnd; ++c) {
1313:       PetscCall(PetscFindInt(c, Nl, leaves, &l));
1314:       if (l < 0) points[Ncoh++] = c;
1315:     }
1316:     PetscCall(ISCreateGeneral(PETSC_COMM_SELF, Ncoh, points, PETSC_OWN_POINTER, &cohesiveCells));
1317:   } else {
1318:     PetscCall(ISCreateStride(PETSC_COMM_SELF, cEnd - cMax, cMax, 1, &cohesiveCells));
1319:   }
1320:   PetscCall(CreateFaultLabel(dm));
1321:   PetscCall(DMGetLabel(dm, "faultCells", &fault));
1322:   PetscCall(DMGetLocalVector(dm, &locX));
1323:   PetscCall(PetscObjectSetName((PetscObject)locX, "Local Solution"));
1324:   PetscCall(DMGetLocalVector(dm, &locF));
1325:   PetscCall(PetscObjectSetName((PetscObject)locF, "Local Residual"));
1326:   PetscCall(DMCreateMatrix(dm, &J));
1327:   PetscCall(PetscObjectSetName((PetscObject)J, "Jacobian"));

1329:   /* The initial guess has displacement shifted by one unit in each fault parallel direction across the fault */
1330:   PetscCall(CreateMaterialLabel(dm));
1331:   PetscCall(DMGetLabel(dm, "material", &material));
1332:   id              = 1;
1333:   initialGuess[0] = r;
1334:   initialGuess[1] = NULL;
1335:   PetscCall(DMProjectFunctionLabelLocal(dm, 0.0, material, 1, &id, PETSC_DETERMINE, NULL, initialGuess, NULL, INSERT_VALUES, locX));
1336:   id              = 2;
1337:   initialGuess[0] = rp1;
1338:   initialGuess[1] = NULL;
1339:   PetscCall(DMProjectFunctionLabelLocal(dm, 0.0, material, 1, &id, PETSC_DETERMINE, NULL, initialGuess, NULL, INSERT_VALUES, locX));
1340:   id              = 1;
1341:   initialGuess[0] = NULL;
1342:   initialGuess[1] = phi;
1343:   PetscCall(DMProjectFunctionLabelLocal(dm, 0.0, fault, 1, &id, PETSC_DETERMINE, NULL, initialGuess, NULL, INSERT_VALUES, locX));
1344:   PetscCall(PetscObjectViewSynchronizedFromOptions((PetscObject)locX, (PetscObject)dm, "-local_solution_view"));

1346:   // Test projection to fault mesh, which is collective although a process can have no cohesive cells
1347:   hasCohesive = cMax < cEnd ? PETSC_TRUE : PETSC_FALSE;
1348:   PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, &hasCohesive, 1, MPI_C_BOOL, MPI_LOR, PetscObjectComm((PetscObject)dm)));
1349:   if (hasCohesive) {
1350:     PetscCall(DMPlexCreateCohesiveSubmesh(dm, PETSC_FALSE, NULL, 0, &dmFault));
1351:     PetscCall(PetscObjectSetName((PetscObject)dmFault, "Fault Mesh"));
1352:     PetscCall(DMViewFromOptions(dmFault, NULL, "-fault_view"));
1353:     PetscCall(DMPlexOrient(dmFault));
1354:     PetscCall(DMPlexGetHeightStratum(dm, 1, &fStart, NULL));
1355:     PetscCall(DMPlexGetCellType(dm, fStart, &fct));
1356:     //PetscCall(PetscFECreateByCell(PETSC_COMM_SELF, dim - 1, dim, fct, "fault_field_", PETSC_DETERMINE, &fe));
1357:     PetscCall(PetscFECreateDefault(PETSC_COMM_SELF, dim - 1, dim, PETSC_TRUE, "fault_field_", PETSC_DETERMINE, &fe));
1358:     PetscCall(PetscFESetName(fe, "fault_field"));
1359:     PetscCall(DMAddField(dmFault, NULL, (PetscObject)fe));
1360:     PetscCall(PetscFEDestroy(&fe));
1361:     PetscCall(DMCreateDS(dmFault));
1362:     PetscCall(DMGetLocalVector(dmFault, &locW));
1363:     PetscCall(DMViewFromOptions(dmFault, NULL, "-cohesive_view"));
1364:     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[]);

1366:     DMLabel  depthLabel;
1367:     PetscInt depth;
1368:     PetscCall(DMPlexGetDepthLabel(dmFault, &depthLabel));
1369:     PetscCall(DMPlexGetDepth(dmFault, &depth));
1370:     id = depth - 1;
1371:     /* w = r + rp1 */
1372:     faultFuncs[0] = add_fields;
1373:     PetscCall(DMProjectBdFieldLabelLocal(dmFault, 0.0, depthLabel, 1, &id, PETSC_DETERMINE, NULL, locX, faultFuncs, INSERT_VALUES, locW));
1374:     PetscCall(PetscObjectViewSynchronizedFromOptions((PetscObject)locW, (PetscObject)dm, "-local_projection_view"));

1376:     /* w = fault_normal */
1377:     faultFuncs[0] = normal_field;
1378:     PetscCall(DMProjectBdFieldLabelLocal(dmFault, 0.0, depthLabel, 1, &id, PETSC_DETERMINE, NULL, locX, faultFuncs, INSERT_VALUES, locW));
1379:     PetscCall(PetscObjectViewSynchronizedFromOptions((PetscObject)locW, (PetscObject)dm, "-local_projection_view"));
1380:     PetscCall(DMRestoreLocalVector(dmFault, &locW));
1381:     PetscCall(DMDestroy(&dmFault));
1382:   }

1384:   PetscCall(DMGetCellDS(dm, cMax, &probh, NULL));
1385:   PetscCall(PetscDSGetWeakForm(probh, &wf));
1386:   PetscCall(PetscDSGetNumFields(probh, &Nf));
1387:   PetscCall(PetscWeakFormSetIndexBdResidual(wf, material, 1, 0, 0, 0, f0_bd_u_neg, 0, NULL));
1388:   PetscCall(PetscWeakFormSetIndexBdResidual(wf, material, 2, 0, 0, 0, f0_bd_u_pos, 0, NULL));
1389:   PetscCall(PetscWeakFormSetIndexBdJacobian(wf, material, 1, 0, 1, 0, 0, g0_bd_ul_neg, 0, NULL, 0, NULL, 0, NULL));
1390:   PetscCall(PetscWeakFormSetIndexBdJacobian(wf, material, 2, 0, 1, 0, 0, g0_bd_ul_pos, 0, NULL, 0, NULL, 0, NULL));
1391:   if (Nf > 1) {
1392:     PetscCall(PetscWeakFormSetIndexBdResidual(wf, fault, 1, 1, 0, 0, f0_bd_l, 0, NULL));
1393:     PetscCall(PetscWeakFormSetIndexBdJacobian(wf, fault, 1, 1, 0, 0, 0, g0_bd_lu, 0, NULL, 0, NULL, 0, NULL));
1394:   }
1395:   if (rank == 0) PetscCall(PetscDSView(probh, NULL));

1397:   keys[0].label = material;
1398:   keys[0].value = 1;
1399:   keys[0].field = 0;
1400:   keys[0].part  = 0;
1401:   keys[1].label = material;
1402:   keys[1].value = 2;
1403:   keys[1].field = 0;
1404:   keys[1].part  = 0;
1405:   keys[2].label = fault;
1406:   keys[2].value = 1;
1407:   keys[2].field = 1;
1408:   keys[2].part  = 0;
1409:   PetscCall(VecSet(locF, 0.));
1410:   PetscCall(DMPlexComputeResidualHybridByKey(dm, keys, cohesiveCells, 0.0, locX, NULL, 0.0, locF, user));
1411:   PetscCall(PetscObjectViewSynchronizedFromOptions((PetscObject)locF, (PetscObject)dm, "-local_residual_view"));
1412:   PetscCall(MatZeroEntries(J));
1413:   PetscCall(DMPlexComputeJacobianHybridByKey(dm, keys, cohesiveCells, 0.0, 0.0, locX, NULL, J, J, user));
1414:   PetscCall(MatAssemblyBegin(J, MAT_FINAL_ASSEMBLY));
1415:   PetscCall(MatAssemblyEnd(J, MAT_FINAL_ASSEMBLY));
1416:   PetscCall(MatViewFromOptions(J, NULL, "-local_jacobian_view"));

1418:   PetscCall(DMRestoreLocalVector(dm, &locX));
1419:   PetscCall(DMRestoreLocalVector(dm, &locF));
1420:   PetscCall(MatDestroy(&J));
1421:   PetscCall(ISDestroy(&cohesiveCells));

1423:   if (cMax < cEnd) {
1424:     PetscDS         ds;
1425:     PetscFE         fe;
1426:     PetscQuadrature quad;
1427:     IS             *perm;
1428:     const PetscInt *cone;
1429:     PetscInt        Na, a;

1431:     PetscCall(DMPlexGetCone(dm, cMax, &cone));
1432:     PetscCall(DMGetCellDS(dm, cMax, &ds, NULL));
1433:     PetscCall(PetscDSGetDiscretization(ds, 0, (PetscObject *)&fe));
1434:     PetscCall(PetscFEGetQuadrature(fe, &quad));
1435:     PetscCall(PetscQuadratureComputePermutations(quad, &Na, &perm));
1436:     for (a = 0; a < Na; ++a) PetscCall(ISDestroy(&perm[a]));
1437:     PetscCall(PetscFree(perm));
1438:   }
1439:   PetscFunctionReturn(PETSC_SUCCESS);
1440: }

1442: int main(int argc, char **argv)
1443: {
1444:   DM     dm;
1445:   AppCtx user;

1447:   PetscFunctionBeginUser;
1448:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
1449:   PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user));
1450:   PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm));
1451:   PetscCall(CreateDiscretization(dm, &user));
1452:   PetscCall(TestAssembly(dm, &user));
1453:   PetscCall(DMDestroy(&dm));
1454:   PetscCall(PetscFinalize());
1455:   return 0;
1456: }

1458: /*TEST

1460:   testset:
1461:     requires: triangle
1462:     args: -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1463:             -dm_plex_transform_active fault -dm_plex_save_transform -dm_plex_check_transform \
1464:           -dm_view ::ascii_info_detail -coarse_dm_view ::ascii_info_detail \
1465:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1466:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1467:     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"

1469:     test:
1470:       suffix: tri_0
1471:       args: -dm_plex_box_faces 1,1 -dm_plex_cohesive_label_fault 8
1472:     test:
1473:       suffix: tri_1
1474:       args: -dm_plex_box_faces 1,1 -dm_plex_cohesive_label_fault 8 \
1475:               -dm_plex_transform_extrude_use_tensor 0
1476:     test:
1477:       suffix: tri_2
1478:       args: -dm_plex_file_contents dat:tri_2_cv -dm_plex_cohesive_label_fault 11,15
1479:     test:
1480:       suffix: tri_2_perm
1481:       args: -dm_plex_file_contents dat:tri_2_cv -dm_plex_cohesive_label_fault 11,15 \
1482:             -dm_reorder_section -dm_reorder_section_type cohesive
1483:     # Note that the mesh is not parallel when the cohesive label is oriented
1484:     test:
1485:       suffix: tri_3
1486:       nsize: 2
1487:       args: -dm_plex_file_contents dat:tri_2_cv -dm_plex_cohesive_label_fault 11,15 \
1488:               -petscpartitioner_type shell -petscpartitioner_shell_sizes 2,2 \
1489:               -petscpartitioner_shell_points 0,3,1,2

1491:   testset:
1492:     requires: defined(PETSC_HAVE_EXECUTABLE_EXPORT)
1493:     nsize: 3
1494:     args: -dm_plex_file_contents dat:tri_3x3_cv -dm_plex_cohesive_label_fault 37,42,46 -petscpartitioner_type shell \
1495:           -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1496:             -dm_plex_transform_active fault -dm_plex_save_transform -dm_plex_check_transform \
1497:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1498:             -local_solution_view -local_residual_view
1499:     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"

1501:     # The owner of a shared fault vertex holds no fault edge next to it
1502:     test:
1503:       suffix: tri_4
1504:       args: -petscpartitioner_shell_sizes 4,5,9 -petscpartitioner_shell_points 6,7,10,15,0,1,3,8,16,2,4,5,9,11,12,13,14,17
1505:     # A process holds a shared fault vertex that it does not own, but no fault edge next to it
1506:     test:
1507:       suffix: tri_5
1508:       args: -petscpartitioner_shell_sizes 7,6,5 -petscpartitioner_shell_points 1,3,4,6,12,14,16,2,5,7,10,11,17,0,8,9,13,15
1509:     # With overlap, every process holds copies of fault edges that other processes own
1510:     test:
1511:       suffix: tri_6
1512:       args: -petscpartitioner_shell_sizes 7,6,5 -petscpartitioner_shell_points 1,3,4,6,12,14,16,2,5,7,10,11,17,0,8,9,13,15 \
1513:             -dm_distribute_overlap 1
1514:     # A process holds vertices of the fault, but no fault edge
1515:     test:
1516:       suffix: tri_7
1517:       args: -petscpartitioner_shell_sizes 6,6,6 -petscpartitioner_shell_points 1,3,6,8,13,17,5,7,9,10,15,16,0,2,4,11,12,14
1518:     # With overlap, the fault mesh made from the cohesive cells also has overlap
1519:     test:
1520:       suffix: tri_8
1521:       args: -petscpartitioner_shell_sizes 7,7,4 -petscpartitioner_shell_points 0,1,2,5,6,9,17,4,8,10,11,12,13,14,3,7,15,16 \
1522:             -dm_distribute_overlap 1

1524:   testset:
1525:     requires: triangle
1526:     args: -dm_plex_option_phases coh_,ref_ \
1527:             -coh_dm_refine 1 -coh_dm_plex_transform_type cohesive_extrude \
1528:               -coh_dm_plex_transform_active fault \
1529:             -ref_dm_refine 1 -ref_dm_plex_transform_type refine_regular \
1530:           -dm_view ::ascii_info_detail \
1531:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1532:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1533:     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"

1535:     test:
1536:       suffix: tri_0_ref
1537:       args: -dm_plex_box_faces 1,1 -dm_plex_cohesive_label_fault 8
1538:     test:
1539:       suffix: tri_2_ref
1540:       args: -dm_plex_file_contents dat:tri_2_cv -dm_plex_cohesive_label_fault 11,15
1541:     test:
1542:       suffix: tri_2_ref_perm
1543:       args: -dm_plex_file_contents dat:tri_2_cv -dm_plex_cohesive_label_fault 11,15 \
1544:             -dm_reorder_section -dm_reorder_section_type cohesive

1546:   testset:
1547:     args: -dm_plex_simplex 0 -dm_plex_box_faces 2,1 \
1548:           -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1549:             -dm_plex_transform_active fault -dm_plex_cohesive_label_fault 13 \
1550:             -dm_plex_save_transform -dm_plex_check_transform \
1551:           -dm_view ::ascii_info_detail -coarse_dm_view ::ascii_info_detail \
1552:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1553:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1554:     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"

1556:     test:
1557:       suffix: quad_0
1558:     test:
1559:       suffix: quad_1
1560:       args: -dm_plex_transform_extrude_use_tensor 0
1561:     test:
1562:       suffix: quad_2
1563:       nsize: 2
1564:       args: -petscpartitioner_type simple

1566:   test:
1567:     suffix: quad_3
1568:     nsize: 4
1569:     args: -test_num 1 \
1570:           -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1571:             -dm_plex_transform_active fault -dm_plex_save_transform -dm_plex_check_transform \
1572:           -dm_view ::ascii_info_detail -coarse_dm_view ::ascii_info_detail \
1573:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1574:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1575:     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"

1577:   test:
1578:     suffix: quad_4
1579:     args: -dm_plex_simplex 0 -dm_plex_box_faces 3,2 \
1580:           -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1581:             -dm_plex_transform_active fault -dm_plex_cohesive_label_fault 22,23 \
1582:             -dm_plex_save_transform -dm_plex_check_transform \
1583:           -dm_view ::ascii_info_detail -coarse_dm_view ::ascii_info_detail \
1584:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1585:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1586:     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"

1588:   test:
1589:     suffix: quad_5
1590:     args: -dm_plex_simplex 0 -dm_plex_box_faces 3,2 \
1591:             -dm_plex_cohesive_label_fault0 21 \
1592:             -dm_plex_cohesive_label_fault1 23 \
1593:           -f0_dm_refine 1 -f0_dm_plex_transform_type cohesive_extrude \
1594:             -f0_dm_plex_transform_active fault0  -f0_coarse_dm_view ::ascii_info_detail \
1595:           -f1_dm_refine 1 -f1_dm_plex_transform_type cohesive_extrude \
1596:             -f1_dm_plex_transform_active fault1  -f1_coarse_dm_view ::ascii_info_detail \
1597:           -dm_plex_save_transform -dm_plex_check_transform \
1598:           -dm_view ::ascii_info_detail
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: quad_6
1603:     args: -dm_plex_simplex 0 -dm_plex_box_faces 3,2 \
1604:             -dm_plex_cohesive_label_fault0 22,23 \
1605:             -dm_plex_cohesive_label_fault1 32 \
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: quad_6w
1616:     args: -dm_plex_simplex 0 -dm_plex_box_faces 3,2 \
1617:             -dm_plex_cohesive_label_fault0 22,23 \
1618:             -dm_plex_cohesive_label_fault1 32 \
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:             -f0_dm_plex_transform_cohesive_width 0.05 \
1622:           -f1_dm_refine 1 -f1_dm_plex_transform_type cohesive_extrude \
1623:             -f1_dm_plex_transform_active fault1  -f1_coarse_dm_view ::ascii_info_detail \
1624:             -f1_dm_plex_transform_cohesive_width 0.05 \
1625:           -dm_plex_save_transform -dm_plex_check_transform \
1626:           -dm_view ::ascii_info_detail
1627:     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"

1629:   testset:
1630:     args: -dm_plex_simplex 0 -dm_plex_box_faces 2,1 -dm_plex_cohesive_label_fault 13 \
1631:           -dm_plex_option_phases coh_,ref_ \
1632:             -coh_dm_refine 1 -coh_dm_plex_transform_type cohesive_extrude \
1633:               -coh_dm_plex_transform_active fault \
1634:             -ref_dm_refine 1 -ref_dm_plex_transform_type refine_regular \
1635:           -dm_view ::ascii_info_detail \
1636:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1637:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1638:     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"

1640:     test:
1641:       suffix: quad_0_ref

1643:   testset:
1644:     args: -dm_plex_dim 3 -dm_plex_shape doublet \
1645:           -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1646:             -dm_plex_transform_active fault -dm_plex_cohesive_label_fault 7 \
1647:             -dm_plex_save_transform -dm_plex_check_transform \
1648:           -dm_view ::ascii_info_detail -coarse_dm_view ::ascii_info_detail \
1649:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1650:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1651:     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"

1653:     test:
1654:       suffix: tet_0
1655:     test:
1656:       suffix: tet_1
1657:       nsize: 2
1658:       args: -petscpartitioner_type simple

1660:   testset:
1661:     args: -dm_plex_dim 3 -dm_plex_shape doublet -dm_plex_cohesive_label_fault 7 \
1662:           -dm_plex_option_phases coh_,ref_ \
1663:             -coh_dm_refine 1 -coh_dm_plex_transform_type cohesive_extrude \
1664:               -coh_dm_plex_transform_active fault \
1665:             -ref_dm_refine 1 -ref_dm_plex_transform_type refine_regular \
1666:           -dm_view ::ascii_info_detail \
1667:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1668:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1669:     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"

1671:     test:
1672:       suffix: tet_0_ref

1674:   testset:
1675:     args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_plex_box_faces 2,1,1 -dm_plex_box_upper 2,1,1 \
1676:           -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1677:             -dm_plex_transform_active fault -dm_plex_cohesive_label_fault 15 \
1678:             -dm_plex_save_transform -dm_plex_check_transform \
1679:           -dm_view ::ascii_info_detail -coarse_dm_view ::ascii_info_detail \
1680:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1681:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1682:     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"

1684:     test:
1685:       suffix: hex_0
1686:     test:
1687:       suffix: hex_1
1688:       nsize: 2
1689:       args: -petscpartitioner_type simple

1691:   test:
1692:     suffix: hex_2
1693:     nsize: 4
1694:     args: -test_num 2 \
1695:           -dm_refine 1 -dm_plex_transform_type cohesive_extrude \
1696:             -dm_plex_transform_active fault -dm_plex_save_transform -dm_plex_check_transform \
1697:           -dm_view ::ascii_info_detail -coarse_dm_view ::ascii_info_detail \
1698:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1699:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1700:     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"

1702:   test:
1703:     suffix: hex_3
1704:     args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_plex_box_faces 2,1,2 -dm_plex_box_upper 2.,1.,2. \
1705:             -dm_plex_cohesive_label_fault0 37,40 \
1706:             -dm_plex_cohesive_label_fault1 26 \
1707:           -f0_dm_refine 1 -f0_dm_plex_transform_type cohesive_extrude \
1708:             -f0_dm_plex_transform_active fault0  -f0_coarse_dm_view ::ascii_info_detail \
1709:           -f1_dm_refine 1 -f1_dm_plex_transform_type cohesive_extrude \
1710:             -f1_dm_plex_transform_active fault1  -f1_coarse_dm_view ::ascii_info_detail \
1711:           -dm_plex_save_transform -dm_plex_check_transform \
1712:           -dm_view ::ascii_info_detail
1713:     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"

1715:   test:
1716:     suffix: hex_4
1717:     args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_plex_box_faces 4,1,2 -dm_plex_box_upper 4.,1.,2. \
1718:             -dm_plex_cohesive_label_fault0 65,68 \
1719:             -dm_plex_cohesive_label_fault1 46 \
1720:           -f0_dm_refine 1 -f0_dm_plex_transform_type cohesive_extrude \
1721:             -f0_dm_plex_transform_active fault0  -f0_coarse_dm_view ::ascii_info_detail \
1722:           -f1_dm_refine 1 -f1_dm_plex_transform_type cohesive_extrude \
1723:             -f1_dm_plex_transform_active fault1  -f1_coarse_dm_view ::ascii_info_detail \
1724:           -dm_plex_save_transform -dm_plex_check_transform \
1725:           -dm_view ::ascii_info_detail
1726:     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"

1728:   testset:
1729:     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 \
1730:           -dm_plex_option_phases coh_,ref_ \
1731:             -coh_dm_refine 1 -coh_dm_plex_transform_type cohesive_extrude \
1732:               -coh_dm_plex_transform_active fault \
1733:             -ref_dm_refine 1 -ref_dm_plex_transform_type refine_regular \
1734:           -dm_view ::ascii_info_detail \
1735:           -displacement_petscspace_degree 1 -faulttraction_petscspace_degree 1 \
1736:             -local_section_view -local_solution_view -local_residual_view -local_jacobian_view
1737:     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"

1739:     test:
1740:       suffix: hex_0_ref

1742: TEST*/