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