Actual source code: plexreftobox.c
1: #include <petsc/private/dmplextransformimpl.h>
3: static PetscErrorCode DMPlexTransformView_ToBox(DMPlexTransform tr, PetscViewer viewer)
4: {
5: PetscBool isascii;
7: PetscFunctionBegin;
10: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
11: if (isascii) {
12: const char *name;
14: PetscCall(PetscObjectGetName((PetscObject)tr, &name));
15: PetscCall(PetscViewerASCIIPrintf(viewer, "Transformation to box cells %s\n", name ? name : ""));
16: } else {
17: SETERRQ(PetscObjectComm((PetscObject)tr), PETSC_ERR_SUP, "Viewer type %s not yet supported for DMPlexTransform writing", ((PetscObject)viewer)->type_name);
18: }
19: PetscFunctionReturn(PETSC_SUCCESS);
20: }
22: static PetscErrorCode DMPlexTransformSetUp_ToBox(DMPlexTransform tr)
23: {
24: PetscFunctionBegin;
25: PetscFunctionReturn(PETSC_SUCCESS);
26: }
28: static PetscErrorCode DMPlexTransformDestroy_ToBox(DMPlexTransform tr)
29: {
30: DMPlexRefine_ToBox *f = (DMPlexRefine_ToBox *)tr->data;
32: PetscFunctionBegin;
33: PetscCall(PetscFree(f));
34: PetscFunctionReturn(PETSC_SUCCESS);
35: }
37: static PetscErrorCode DMPlexTransformGetSubcellOrientation_ToBox(DMPlexTransform tr, DMPolytopeType sct, PetscInt sp, PetscInt so, DMPolytopeType tct, PetscInt r, PetscInt o, PetscInt *rnew, PetscInt *onew)
38: {
39: PetscBool convertTensor = PETSC_TRUE;
40: static PetscInt tri_seg[] = {0, 0, 2, 0, 1, 0, 2, 0, 1, 0, 0, 0, 1, 0, 0, 0, 2, 0, 0, 0, 1, 0, 2, 0, 1, 0, 2, 0, 0, 0, 2, 0, 0, 0, 1, 0};
41: static PetscInt tri_quad[] = {1, -3, 0, -3, 2, -4, 0, -2, 2, -2, 1, -2, 2, -1, 1, -4, 0, -1, 0, 0, 1, 0, 2, 0, 1, 1, 2, 2, 0, 1, 2, 3, 0, 3, 1, 2};
42: static PetscInt tseg_seg[] = {0, -1, 0, 0, 0, 0, 0, -1};
43: static PetscInt tseg_quad[] = {1, 2, 0, 2, 1, -3, 0, -3, 0, 0, 1, 0, 0, -1, 1, -1};
44: static PetscInt tet_seg[] = {3, 0, 2, 0, 0, 0, 1, 0, 3, 0, 1, 0, 2, 0, 0, 0, 3, 0, 0, 0, 1, 0, 2, 0, 2, 0, 3, 0, 1, 0, 0, 0, 2, 0, 0, 0, 3, 0, 1, 0, 2, 0, 1, 0, 0, 0, 3, 0, 1, 0, 0, 0, 2, 0, 3, 0, 1, 0, 3, 0, 0, 0, 2, 0,
45: 1, 0, 2, 0, 3, 0, 0, 0, 0, 0, 1, 0, 3, 0, 2, 0, 0, 0, 2, 0, 1, 0, 3, 0, 0, 0, 3, 0, 2, 0, 1, 0, 0, 0, 1, 0, 2, 0, 3, 0, 0, 0, 3, 0, 1, 0, 2, 0, 0, 0, 2, 0, 3, 0, 1, 0, 1, 0, 0, 0, 3, 0, 2, 0,
46: 1, 0, 2, 0, 0, 0, 3, 0, 1, 0, 3, 0, 2, 0, 0, 0, 2, 0, 3, 0, 0, 0, 1, 0, 2, 0, 1, 0, 3, 0, 0, 0, 2, 0, 0, 0, 1, 0, 3, 0, 3, 0, 2, 0, 1, 0, 0, 0, 3, 0, 0, 0, 2, 0, 1, 0, 3, 0, 1, 0, 0, 0, 2, 0};
47: static PetscInt tet_quad[] = {2, 0, 5, -3, 4, 0, 0, 3, 3, 1, 1, 1, 0, 0, 3, 0, 5, 0, 1, 0, 4, -2, 2, 0, 1, 1, 4, -3, 3, -3, 2, 3, 5, -2, 0, 0, 3, 1, 5, 3, 0, 0, 4, 3, 2, 0, 1, -3, 4, 0, 2, 3, 5, -2, 1, -4, 0, -2,
48: 3, 1, 1, -3, 0, -3, 2, -2, 3, 0, 5, 0, 4, 0, 2, -2, 1, -4, 0, -2, 4, -3, 3, -3, 5, 0, 4, -2, 3, -4, 1, 1, 5, 3, 0, 0, 2, -2, 5, 0, 0, 3, 3, 1, 2, -3, 1, -3, 4, -2, 3, -3, 1, 0, 4, -2, 0, -3,
49: 2, -2, 5, -2, 0, -2, 2, -3, 1, -3, 5, -3, 4, 0, 3, -3, 5, -2, 4, 3, 2, 0, 3, -4, 1, 1, 0, -2, 0, 0, 1, 0, 2, 0, 3, 0, 4, 0, 5, 0, 3, 1, 4, 3, 1, -3, 5, 3, 2, -2, 0, 0, 5, 0, 2, -3, 4, -2,
50: 0, 3, 1, 1, 3, 1, 4, 0, 1, -4, 3, 1, 2, 3, 0, 0, 5, -2, 2, 0, 0, 3, 1, 1, 5, -3, 3, -3, 4, 0, 5, -2, 3, -4, 0, -2, 4, 3, 1, -3, 2, 0, 4, -2, 5, 3, 2, -2, 3, -4, 0, -2, 1, 1, 3, -3, 0, -3,
51: 5, -2, 1, 0, 2, 0, 4, -2, 1, 1, 2, 3, 0, 0, 4, -3, 5, 0, 3, -3, 0, -2, 5, -3, 3, -3, 2, -3, 4, -2, 1, -3, 2, -2, 4, -3, 5, 0, 1, -4, 3, 1, 0, -2, 1, -3, 3, 0, 4, 0, 0, -3, 5, -2, 2, -2};
52: static PetscInt tet_hex[] = {2, -2, 3, -2, 1, -10, 0, -13, 3, -10, 1, -13, 2, -10, 0, -10, 1, -2, 2, -13, 3, -13, 0, -2, 3, -13, 2, -10, 0, -2, 1, -2, 2, -13, 0, -10, 3, -2, 1, -13, 0, -13, 3, -10, 2, -2, 1, -10,
53: 0, -10, 1, -2, 3, -10, 2, -10, 1, -10, 3, -13, 0, -13, 2, -2, 3, -2, 0, -2, 1, -13, 2, -13, 1, -13, 0, -13, 2, -13, 3, -2, 0, -2, 2, -2, 1, -2, 3, -13, 2, -10, 1, -10, 0, -10, 3, -10,
54: 0, 0, 1, 0, 2, 0, 3, 0, 1, 17, 2, 17, 0, 17, 3, 16, 2, 16, 0, 16, 1, 16, 3, 17, 1, 16, 0, 17, 3, 17, 2, 16, 0, 16, 3, 16, 1, 0, 2, 17, 3, 0, 1, 17, 0, 0, 2, 0,
55: 2, 17, 3, 0, 0, 16, 1, 0, 3, 17, 0, 0, 2, 16, 1, 16, 0, 17, 2, 0, 3, 16, 1, 17, 3, 16, 2, 16, 1, 17, 0, 17, 2, 0, 1, 16, 3, 0, 0, 0, 1, 0, 3, 17, 2, 17, 0, 16};
56: static PetscInt trip_seg[] = {1, 0, 0, 0, 3, 0, 4, 0, 2, 0, 1, 0, 0, 0, 4, 0, 2, 0, 3, 0, 1, 0, 0, 0, 2, 0, 3, 0, 4, 0, 0, 0, 1, 0, 3, 0, 2, 0, 4, 0, 0, 0, 1, 0, 4, 0, 3, 0, 2, 0, 0, 0, 1, 0, 2, 0, 4, 0, 3, 0,
57: 0, 0, 1, 0, 2, 0, 3, 0, 4, 0, 0, 0, 1, 0, 4, 0, 2, 0, 3, 0, 0, 0, 1, 0, 3, 0, 4, 0, 2, 0, 1, 0, 0, 0, 2, 0, 4, 0, 3, 0, 1, 0, 0, 0, 4, 0, 3, 0, 2, 0, 1, 0, 0, 0, 3, 0, 2, 0, 4, 0};
58: static PetscInt trip_quad[] = {1, 1, 2, 2, 0, 1, 7, -1, 8, -1, 6, -1, 4, -1, 5, -1, 3, -1, 2, 3, 0, 3, 1, 2, 8, -1, 6, -1, 7, -1, 5, -1, 3, -1, 4, -1, 0, 0, 1, 0, 2, 0, 6, -1, 7, -1, 8, -1, 3, -1, 4, -1, 5, -1,
59: 2, -1, 1, -4, 0, -1, 4, 0, 3, 0, 5, 0, 7, 0, 6, 0, 8, 0, 0, -2, 2, -2, 1, -2, 5, 0, 4, 0, 3, 0, 8, 0, 7, 0, 6, 0, 1, -3, 0, -3, 2, -4, 3, 0, 5, 0, 4, 0, 6, 0, 8, 0, 7, 0,
60: 0, 0, 1, 0, 2, 0, 3, 0, 4, 0, 5, 0, 6, 0, 7, 0, 8, 0, 2, 3, 0, 3, 1, 2, 5, 0, 3, 0, 4, 0, 8, 0, 6, 0, 7, 0, 1, 1, 2, 2, 0, 1, 4, 0, 5, 0, 3, 0, 7, 0, 8, 0, 6, 0,
61: 1, -3, 0, -3, 2, -4, 6, -1, 8, -1, 7, -1, 3, -1, 5, -1, 4, -1, 0, -2, 2, -2, 1, -2, 8, -1, 7, -1, 6, -1, 5, -1, 4, -1, 3, -1, 2, -1, 1, -4, 0, -1, 7, -1, 6, -1, 8, -1, 4, -1, 3, -1, 5, -1};
62: static PetscInt trip_hex[] = {4, -12, 5, -6, 3, -12, 1, -12, 2, -6, 0, -12, 5, -11, 3, -11, 4, -6, 2, -11, 0, -11, 1, -6, 3, -9, 4, -9, 5, -9, 0, -9, 1, -9, 2, -9, 2, -3, 1, -4, 0, -3, 5, -3, 4, -4, 3, -3,
63: 0, -2, 2, -2, 1, -2, 3, -2, 5, -2, 4, -2, 1, -1, 0, -1, 2, -4, 4, -1, 3, -1, 5, -4, 0, 0, 1, 0, 2, 0, 3, 0, 4, 0, 5, 0, 2, 1, 0, 1, 1, 2, 5, 1, 3, 1, 4, 2,
64: 1, 3, 2, 2, 0, 3, 4, 3, 5, 2, 3, 3, 4, 8, 3, 8, 5, 11, 1, 8, 0, 8, 2, 11, 3, 10, 5, 10, 4, 10, 0, 10, 2, 10, 1, 10, 5, 5, 4, 11, 3, 5, 2, 5, 1, 11, 0, 5};
65: static PetscInt ctrip_seg[] = {0, -1, 0, -1, 0, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, -1, 0, -1};
66: static PetscInt ctrip_quad[] = {0, -1, 2, -1, 1, -1, 2, -1, 1, -1, 0, -1, 1, -1, 0, -1, 2, -1, 0, 0, 2, 0, 1, 0, 2, 0, 1, 0, 0, 0, 1, 0, 0, 0, 2, 0,
67: 0, 0, 1, 0, 2, 0, 1, 0, 2, 0, 0, 0, 2, 0, 0, 0, 1, 0, 0, -1, 1, -1, 2, -1, 1, -1, 2, -1, 0, -1, 2, -1, 0, -1, 1, -1};
68: static PetscInt ctrip_hex[] = {1, 8, 0, 8, 2, 11, 0, 10, 2, 10, 1, 10, 2, 5, 1, 11, 0, 5, 1, -1, 0, -1, 2, -4, 0, -2, 2, -2, 1, -2, 2, -3, 1, -4, 0, -3,
69: 0, 0, 1, 0, 2, 0, 1, 3, 2, 2, 0, 3, 2, 1, 0, 1, 1, 2, 0, -9, 1, -9, 2, -9, 1, -12, 2, -6, 0, -12, 2, -11, 0, -11, 1, -6};
70: static PetscInt tquadp_seg[] = {0, -1, 0, -1, 0, -1, 0, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, -1, 0, -1, 0, -1};
71: static PetscInt tquadp_quad[] = {1, -1, 0, -1, 3, -1, 2, -1, 0, -1, 3, -1, 2, -1, 1, -1, 3, -1, 2, -1, 1, -1, 0, -1, 2, -1, 1, -1, 0, -1, 3, -1, 1, 0, 0, 0, 3, 0, 2, 0, 0, 0, 3,
72: 0, 2, 0, 1, 0, 3, 0, 2, 0, 1, 0, 0, 0, 2, 0, 1, 0, 0, 0, 3, 0, 0, 0, 1, 0, 2, 0, 3, 0, 1, 0, 2, 0, 3, 0, 0, 0, 2, 0, 3, 0, 0, 0,
73: 1, 0, 3, 0, 0, 0, 1, 0, 2, 0, 0, -1, 1, -1, 2, -1, 3, -1, 1, -1, 2, -1, 3, -1, 0, -1, 2, -1, 3, -1, 0, -1, 1, -1, 3, -1, 0, -1, 1, -1, 2, -1};
74: static PetscInt tquadp_hex[] = {2, 11, 1, 11, 0, 11, 3, 11, 1, 8, 0, 8, 3, 8, 2, 8, 0, 10, 3, 10, 2, 10, 1, 10, 3, 5, 2, 5, 1, 5, 0, 5, 2, -4, 1, -4, 0, -4, 3, -4, 1, -1, 0,
75: -1, 3, -1, 2, -1, 0, -2, 3, -2, 2, -2, 1, -2, 3, -3, 2, -3, 1, -3, 0, -3, 0, 0, 1, 0, 2, 0, 3, 0, 1, 3, 2, 3, 3, 3, 0, 3, 2, 2, 3, 2, 0, 2,
76: 1, 2, 3, 1, 0, 1, 1, 1, 2, 1, 0, -9, 1, -9, 2, -9, 3, -9, 1, -12, 2, -12, 3, -12, 0, -12, 2, -6, 3, -6, 0, -6, 1, -6, 3, -11, 0, -11, 1, -11, 2, -11};
78: PetscFunctionBeginHot;
79: *rnew = r;
80: *onew = o;
81: if (!so) PetscFunctionReturn(PETSC_SUCCESS);
82: if (convertTensor) {
83: switch (sct) {
84: case DM_POLYTOPE_POINT:
85: case DM_POLYTOPE_SEGMENT:
86: case DM_POLYTOPE_POINT_PRISM_TENSOR:
87: case DM_POLYTOPE_QUADRILATERAL:
88: case DM_POLYTOPE_HEXAHEDRON:
89: PetscCall(DMPlexTransformGetSubcellOrientation_Regular(tr, sct, sp, so, tct, r, o, rnew, onew));
90: break;
91: case DM_POLYTOPE_TRIANGLE:
92: switch (tct) {
93: case DM_POLYTOPE_POINT:
94: break;
95: case DM_POLYTOPE_SEGMENT:
96: *rnew = tri_seg[(so + 3) * 6 + r * 2];
97: *onew = DMPolytopeTypeComposeOrientation(tct, o, tri_seg[(so + 3) * 6 + r * 2 + 1]);
98: break;
99: case DM_POLYTOPE_QUADRILATERAL:
100: *rnew = tri_quad[(so + 3) * 6 + r * 2];
101: *onew = DMPolytopeTypeComposeOrientation(tct, o, tri_quad[(so + 3) * 6 + r * 2 + 1]);
102: break;
103: default:
104: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
105: }
106: break;
107: case DM_POLYTOPE_SEG_PRISM_TENSOR:
108: switch (tct) {
109: case DM_POLYTOPE_SEGMENT:
110: case DM_POLYTOPE_POINT_PRISM_TENSOR:
111: *rnew = tseg_seg[(so + 2) * 2 + r * 2];
112: *onew = DMPolytopeTypeComposeOrientation(tct, o, tseg_seg[(so + 2) * 2 + r * 2 + 1]);
113: break;
114: case DM_POLYTOPE_QUADRILATERAL:
115: *rnew = tseg_quad[(so + 2) * 4 + r * 2];
116: *onew = DMPolytopeTypeComposeOrientation(tct, o, tseg_quad[(so + 2) * 4 + r * 2 + 1]);
117: break;
118: default:
119: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
120: }
121: break;
122: case DM_POLYTOPE_TETRAHEDRON:
123: switch (tct) {
124: case DM_POLYTOPE_POINT:
125: break;
126: case DM_POLYTOPE_SEGMENT:
127: *rnew = tet_seg[(so + 12) * 8 + r * 2];
128: *onew = DMPolytopeTypeComposeOrientation(tct, o, tet_seg[(so + 12) * 8 + r * 2 + 1]);
129: break;
130: case DM_POLYTOPE_QUADRILATERAL:
131: *rnew = tet_quad[(so + 12) * 12 + r * 2];
132: *onew = DMPolytopeTypeComposeOrientation(tct, o, tet_quad[(so + 12) * 12 + r * 2 + 1]);
133: break;
134: case DM_POLYTOPE_HEXAHEDRON:
135: *rnew = tet_hex[(so + 12) * 8 + r * 2];
136: *onew = DMPolytopeTypeComposeOrientation(tct, o, tet_hex[(so + 12) * 8 + r * 2 + 1]);
137: break;
138: default:
139: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
140: }
141: break;
142: case DM_POLYTOPE_TRI_PRISM:
143: switch (tct) {
144: case DM_POLYTOPE_POINT:
145: break;
146: case DM_POLYTOPE_SEGMENT:
147: *rnew = trip_seg[(so + 6) * 10 + r * 2];
148: *onew = DMPolytopeTypeComposeOrientation(tct, o, trip_seg[(so + 6) * 10 + r * 2 + 1]);
149: break;
150: case DM_POLYTOPE_QUADRILATERAL:
151: *rnew = trip_quad[(so + 6) * 18 + r * 2];
152: *onew = DMPolytopeTypeComposeOrientation(tct, o, trip_quad[(so + 6) * 18 + r * 2 + 1]);
153: break;
154: case DM_POLYTOPE_HEXAHEDRON:
155: *rnew = trip_hex[(so + 6) * 12 + r * 2];
156: *onew = DMPolytopeTypeComposeOrientation(tct, o, trip_hex[(so + 6) * 12 + r * 2 + 1]);
157: break;
158: default:
159: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
160: }
161: break;
162: case DM_POLYTOPE_TRI_PRISM_TENSOR:
163: switch (tct) {
164: case DM_POLYTOPE_SEGMENT:
165: *rnew = ctrip_seg[(so + 6) * 2 + r * 2];
166: *onew = DMPolytopeTypeComposeOrientation(tct, o, ctrip_seg[(so + 6) * 2 + r * 2 + 1]);
167: break;
168: case DM_POLYTOPE_QUADRILATERAL:
169: *rnew = ctrip_quad[(so + 6) * 6 + r * 2];
170: *onew = DMPolytopeTypeComposeOrientation(tct, o, ctrip_quad[(so + 6) * 6 + r * 2 + 1]);
171: break;
172: case DM_POLYTOPE_HEXAHEDRON:
173: *rnew = ctrip_hex[(so + 6) * 6 + r * 2];
174: *onew = DMPolytopeTypeComposeOrientation(tct, o, ctrip_hex[(so + 6) * 6 + r * 2 + 1]);
175: break;
176: default:
177: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
178: }
179: break;
180: case DM_POLYTOPE_QUAD_PRISM_TENSOR:
181: switch (tct) {
182: case DM_POLYTOPE_SEGMENT:
183: *rnew = tquadp_seg[(so + 8) * 2 + r * 2];
184: *onew = DMPolytopeTypeComposeOrientation(tct, o, tquadp_seg[(so + 8) * 2 + r * 2 + 1]);
185: break;
186: case DM_POLYTOPE_QUADRILATERAL:
187: *rnew = tquadp_quad[(so + 8) * 8 + r * 2];
188: *onew = DMPolytopeTypeComposeOrientation(tct, o, tquadp_quad[(so + 8) * 8 + r * 2 + 1]);
189: break;
190: case DM_POLYTOPE_HEXAHEDRON:
191: *rnew = tquadp_hex[(so + 8) * 8 + r * 2];
192: *onew = DMPolytopeTypeComposeOrientation(tct, o, tquadp_hex[(so + 8) * 8 + r * 2 + 1]);
193: break;
194: default:
195: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
196: }
197: break;
198: default:
199: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported cell type %s", DMPolytopeTypes[sct]);
200: }
201: } else {
202: switch (sct) {
203: case DM_POLYTOPE_POINT:
204: case DM_POLYTOPE_SEGMENT:
205: case DM_POLYTOPE_POINT_PRISM_TENSOR:
206: case DM_POLYTOPE_QUADRILATERAL:
207: case DM_POLYTOPE_SEG_PRISM_TENSOR:
208: case DM_POLYTOPE_HEXAHEDRON:
209: case DM_POLYTOPE_QUAD_PRISM_TENSOR:
210: PetscCall(DMPlexTransformGetSubcellOrientation_Regular(tr, sct, sp, so, tct, r, o, rnew, onew));
211: break;
212: default:
213: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported cell type %s", DMPolytopeTypes[sct]);
214: }
215: }
216: PetscFunctionReturn(PETSC_SUCCESS);
217: }
219: static PetscErrorCode DMPlexTransformCellRefine_ToBox(DMPlexTransform tr, DMPolytopeType source, PetscInt p, PetscInt *rt, PetscInt *Nt, DMPolytopeType *target[], PetscInt *size[], PetscInt *cone[], PetscInt *ornt[])
220: {
221: PetscBool convertTensor = PETSC_TRUE;
222: /* Change tensor edges to segments */
223: static DMPolytopeType tedgeT[] = {DM_POLYTOPE_SEGMENT};
224: static PetscInt tedgeS[] = {1};
225: static PetscInt tedgeC[] = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 1, 1, 0};
226: static PetscInt tedgeO[] = {0, 0};
227: /* Add 1 vertex, 3 edges inside every triangle, making 3 new quadrilaterals.
228: 2
229: |\
230: | \
231: | \
232: | \
233: 0 1
234: | \
235: | \
236: 2 1
237: |\ / \
238: | 2 1 \
239: | \ / \
240: 1 | 0
241: | 0 \
242: | | \
243: | | \
244: 0-0-0-----1-----1
245: */
246: static DMPolytopeType triT[] = {DM_POLYTOPE_POINT, DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL};
247: static PetscInt triS[] = {1, 3, 3};
248: static PetscInt triC[] = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 0, 0, DM_POLYTOPE_POINT, 1, 1, 0, DM_POLYTOPE_POINT, 0, 0, DM_POLYTOPE_POINT, 1, 2, 0, DM_POLYTOPE_POINT, 0, 0, DM_POLYTOPE_SEGMENT, 1, 0, 0, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 0, 2, DM_POLYTOPE_SEGMENT, 1, 2, 1, DM_POLYTOPE_SEGMENT, 1, 0, 1, DM_POLYTOPE_SEGMENT, 1, 1, 0, DM_POLYTOPE_SEGMENT, 0, 1, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 0, 2, DM_POLYTOPE_SEGMENT, 0, 1, DM_POLYTOPE_SEGMENT, 1, 1, 1, DM_POLYTOPE_SEGMENT, 1, 2, 0};
249: static PetscInt triO[] = {0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 0, 0, -1, 0, -1, 0, 0};
250: /* Add 1 edge inside every tensor quad, making 2 new quadrilaterals
251: 2----2----1----3----3
252: | | |
253: | | |
254: | | |
255: 4 A 6 B 5
256: | | |
257: | | |
258: | | |
259: 0----0----0----1----1
260: */
261: static DMPolytopeType tsegT[] = {DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL};
262: static PetscInt tsegS[] = {1, 2};
263: static PetscInt tsegC[] = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 1, 1, 0,
264: /* TODO Fix these */
265: DM_POLYTOPE_SEGMENT, 1, 0, 0, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 1, 1, 0, DM_POLYTOPE_SEGMENT, 1, 2, 0, DM_POLYTOPE_SEGMENT, 1, 0, 1, DM_POLYTOPE_SEGMENT, 1, 3, 0, DM_POLYTOPE_SEGMENT, 1, 1, 1, DM_POLYTOPE_SEGMENT, 0, 0};
266: static PetscInt tsegO[] = {0, 0, 0, 0, -1, -1, 0, 0, -1, -1};
267: /* Add 6 triangles inside every cell, making 4 new hexs
268: The vertices of our reference tet are [(-1, -1, -1), (-1, 1, -1), (1, -1, -1), (-1, -1, 1)], which we call [v0, v1, v2, v3]. The first
269: three edges are [v0, v1], [v1, v2], [v2, v0] called e0, e1, and e2, and then three edges to the top point [v0, v3], [v1, v3], [v2, v3]
270: called e3, e4, and e5. The faces of a tet, given in DMPlexGetRawFaces_Internal() are
271: [v0, v1, v2], [v0, v3, v1], [v0, v2, v3], [v2, v1, v3]
272: We make a new hex in each corner
273: [v0, (e0, 0), (f0, 0), (e2, 0), (e3, 0), (f2, 0), (c0, 0), (f1, 0)]
274: [v1, (e4, 0), (f3, 0), (e1, 0), (e0, 0), (f0, 0), (c0, 0), (f1, 0)]
275: [v2, (e1, 0), (f3, 0), (e5, 0), (e2, 0), (f2, 0), (c0, 0), (f0, 0)]
276: [v3, (e4, 0), (f1, 0), (e3, 0), (e5, 0), (f2, 0), (c0, 0), (f3, 0)]
277: We create a new face for each edge
278: [(e3, 0), (f2, 0), (c0, 0), (f1, 0)]
279: [(f0, 0), (e0, 0), (f1, 0), (c0, 0)]
280: [(e2, 0), (f0, 0), (c0, 0), (f2, 0)]
281: [(f3, 0), (e4, 0), (f1, 0), (c0, 0)]
282: [(e1, 0), (f3, 0), (c0, 0), (f0, 0)]
283: [(e5, 0), (f3, 0), (c0, 0), (f2, 0)]
284: I could write a program to generate these from the first hex by acting with the symmetry group to take one subcell into another.
285: */
286: static DMPolytopeType tetT[] = {DM_POLYTOPE_POINT, DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL, DM_POLYTOPE_HEXAHEDRON};
287: static PetscInt tetS[] = {1, 4, 6, 4};
288: static PetscInt tetC[] = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 0, 0, DM_POLYTOPE_POINT, 1, 1, 0, DM_POLYTOPE_POINT, 0, 0, DM_POLYTOPE_POINT, 1, 2, 0, DM_POLYTOPE_POINT, 0, 0, DM_POLYTOPE_POINT, 1, 3, 0, DM_POLYTOPE_POINT, 0, 0, DM_POLYTOPE_SEGMENT, 1, 2, 2, DM_POLYTOPE_SEGMENT, 0, 2, DM_POLYTOPE_SEGMENT, 0, 1, DM_POLYTOPE_SEGMENT, 1, 1, 0, DM_POLYTOPE_SEGMENT, 1, 0, 0, DM_POLYTOPE_SEGMENT, 1, 1, 2, DM_POLYTOPE_SEGMENT, 0, 1, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 1, 0, 2, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 0, 2, DM_POLYTOPE_SEGMENT, 1, 2, 0, DM_POLYTOPE_SEGMENT, 1, 3, 1, DM_POLYTOPE_SEGMENT, 1, 1, 1, DM_POLYTOPE_SEGMENT, 0, 1, DM_POLYTOPE_SEGMENT, 0, 3, DM_POLYTOPE_SEGMENT, 1, 3, 0, DM_POLYTOPE_SEGMENT, 0, 3, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 1, 0, 1, DM_POLYTOPE_SEGMENT, 1, 3, 2, DM_POLYTOPE_SEGMENT, 0, 3, DM_POLYTOPE_SEGMENT, 0, 2, DM_POLYTOPE_SEGMENT, 1, 2, 1, DM_POLYTOPE_QUADRILATERAL, 1, 0, 0, DM_POLYTOPE_QUADRILATERAL, 0, 0, DM_POLYTOPE_QUADRILATERAL, 1, 2, 0, DM_POLYTOPE_QUADRILATERAL, 0, 1, DM_POLYTOPE_QUADRILATERAL, 0, 2, DM_POLYTOPE_QUADRILATERAL, 1, 1, 0, DM_POLYTOPE_QUADRILATERAL, 1, 3, 1, DM_POLYTOPE_QUADRILATERAL, 0, 1, DM_POLYTOPE_QUADRILATERAL, 1, 0, 1, DM_POLYTOPE_QUADRILATERAL, 0, 3, DM_POLYTOPE_QUADRILATERAL, 0, 4, DM_POLYTOPE_QUADRILATERAL, 1, 1, 2, DM_POLYTOPE_QUADRILATERAL, 1, 3, 0, DM_POLYTOPE_QUADRILATERAL, 0, 2, DM_POLYTOPE_QUADRILATERAL, 1, 2, 1, DM_POLYTOPE_QUADRILATERAL, 0, 4, DM_POLYTOPE_QUADRILATERAL, 0, 5, DM_POLYTOPE_QUADRILATERAL, 1, 0, 2, DM_POLYTOPE_QUADRILATERAL, 1, 1, 1, DM_POLYTOPE_QUADRILATERAL, 0, 5, DM_POLYTOPE_QUADRILATERAL, 1, 2, 2, DM_POLYTOPE_QUADRILATERAL, 0, 3, DM_POLYTOPE_QUADRILATERAL, 0, 0, DM_POLYTOPE_QUADRILATERAL, 1, 3, 2};
289: static PetscInt tetO[] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, -1, -1, 0, 0, -1, 0, 0, -1, -1, -1, 0, 0, -1, 0, 0, -1, -1, 0, 0, -1, -1, 0, 0, 0, 0, 0, 0, 1, -3, 1, 0, 0, 3, 0, -2, 1, -3, 0, 3, 1, -2, 3, -4, -2, 3};
290: /* Add 3 quads inside every triangular prism, making 4 new prisms. */
291: static DMPolytopeType tripT[] = {DM_POLYTOPE_POINT, DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL, DM_POLYTOPE_HEXAHEDRON};
292: static PetscInt tripS[] = {1, 5, 9, 6};
293: static PetscInt tripC[] = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 0, 0, DM_POLYTOPE_POINT, 1, 1, 0, DM_POLYTOPE_POINT, 0, 0, DM_POLYTOPE_POINT, 1, 2, 0, DM_POLYTOPE_POINT, 0, 0, DM_POLYTOPE_POINT, 1, 3, 0, DM_POLYTOPE_POINT, 0, 0, DM_POLYTOPE_POINT, 1, 4, 0, DM_POLYTOPE_POINT, 0, 0, DM_POLYTOPE_SEGMENT, 1, 2, 3, DM_POLYTOPE_SEGMENT, 0, 2, DM_POLYTOPE_SEGMENT, 0, 4, DM_POLYTOPE_SEGMENT, 1, 4, 1, DM_POLYTOPE_SEGMENT, 1, 2, 1, DM_POLYTOPE_SEGMENT, 1, 3, 3, DM_POLYTOPE_SEGMENT, 0, 3, DM_POLYTOPE_SEGMENT, 0, 2, DM_POLYTOPE_SEGMENT, 0, 4, DM_POLYTOPE_SEGMENT, 0, 3, DM_POLYTOPE_SEGMENT, 1, 3, 1, DM_POLYTOPE_SEGMENT, 1, 4, 3, DM_POLYTOPE_SEGMENT, 1, 0, 2, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 0, 2, DM_POLYTOPE_SEGMENT, 1, 2, 0, DM_POLYTOPE_SEGMENT, 1, 0, 1, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 0, 3, DM_POLYTOPE_SEGMENT, 1, 3, 0, DM_POLYTOPE_SEGMENT, 1, 0, 0, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 0, 4, DM_POLYTOPE_SEGMENT, 1, 4, 0, DM_POLYTOPE_SEGMENT, 0, 2, DM_POLYTOPE_SEGMENT, 0, 1, DM_POLYTOPE_SEGMENT, 1, 1, 0, DM_POLYTOPE_SEGMENT, 1, 2, 2, DM_POLYTOPE_SEGMENT, 0, 3, DM_POLYTOPE_SEGMENT, 0, 1, DM_POLYTOPE_SEGMENT, 1, 1, 1, DM_POLYTOPE_SEGMENT, 1, 3, 2, DM_POLYTOPE_SEGMENT, 0, 4, DM_POLYTOPE_SEGMENT, 0, 1, DM_POLYTOPE_SEGMENT, 1, 1, 2, DM_POLYTOPE_SEGMENT, 1, 4, 2, DM_POLYTOPE_QUADRILATERAL, 1, 0, 0, DM_POLYTOPE_QUADRILATERAL, 0, 0, DM_POLYTOPE_QUADRILATERAL, 1, 2, 0, DM_POLYTOPE_QUADRILATERAL, 0, 5, DM_POLYTOPE_QUADRILATERAL, 0, 3, DM_POLYTOPE_QUADRILATERAL, 1, 4, 1, DM_POLYTOPE_QUADRILATERAL, 1, 0, 2, DM_POLYTOPE_QUADRILATERAL, 0, 1, DM_POLYTOPE_QUADRILATERAL, 1, 2, 1, DM_POLYTOPE_QUADRILATERAL, 0, 4, DM_POLYTOPE_QUADRILATERAL, 1, 3, 0, DM_POLYTOPE_QUADRILATERAL, 0, 3, DM_POLYTOPE_QUADRILATERAL, 1, 0, 1, DM_POLYTOPE_QUADRILATERAL, 0, 2, DM_POLYTOPE_QUADRILATERAL, 0, 5, DM_POLYTOPE_QUADRILATERAL, 1, 3, 1, DM_POLYTOPE_QUADRILATERAL, 0, 4, DM_POLYTOPE_QUADRILATERAL, 1, 4, 0, DM_POLYTOPE_QUADRILATERAL, 0, 0, DM_POLYTOPE_QUADRILATERAL, 1, 1, 0, DM_POLYTOPE_QUADRILATERAL, 1, 2, 3, DM_POLYTOPE_QUADRILATERAL, 0, 8, DM_POLYTOPE_QUADRILATERAL, 0, 6, DM_POLYTOPE_QUADRILATERAL, 1, 4, 2, DM_POLYTOPE_QUADRILATERAL, 0, 1, DM_POLYTOPE_QUADRILATERAL, 1, 1, 1, DM_POLYTOPE_QUADRILATERAL, 1, 2, 2, DM_POLYTOPE_QUADRILATERAL, 0, 7, DM_POLYTOPE_QUADRILATERAL, 1, 3, 3, DM_POLYTOPE_QUADRILATERAL, 0, 6, DM_POLYTOPE_QUADRILATERAL, 0, 2, DM_POLYTOPE_QUADRILATERAL, 1, 1, 2, DM_POLYTOPE_QUADRILATERAL, 0, 8, DM_POLYTOPE_QUADRILATERAL, 1, 3, 2, DM_POLYTOPE_QUADRILATERAL, 0, 7, DM_POLYTOPE_QUADRILATERAL, 1, 4, 3};
294: static PetscInt tripO[] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, -1, -1, 0, 0, -1, 0, -1, -1, 0, 0, 0, -1, -1, 0, 0, -1, -1, 0, 0, -1, -1, 0, -1, -1, 0, 0, -1, -1,
295: 0, 0, -1, -1, 0, 0, 0, 0, -3, 0, 1, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, -3, 1, -2, 0, 0, -3, 0, 1, -2, 0, 0, 0, 0, -2, -2, 0, 0, 0, -3, 1};
296: /* Add 3 tensor quads inside every tensor triangular prism, making 4 new tensor triangular prisms.
297: 2
298: |\
299: | \
300: | \
301: 0---1
303: 2
305: 0 1
307: 2
308: |\
309: | \
310: | \
311: 0---1
312: */
313: static DMPolytopeType ttriT[] = {DM_POLYTOPE_POINT_PRISM_TENSOR, DM_POLYTOPE_SEG_PRISM_TENSOR, DM_POLYTOPE_QUAD_PRISM_TENSOR};
314: static PetscInt ttriS[] = {1, 3, 3};
315: static PetscInt ttriC[] = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 1, 1, 0, DM_POLYTOPE_SEGMENT, 1, 0, 0, DM_POLYTOPE_SEGMENT, 1, 1, 0, DM_POLYTOPE_POINT_PRISM_TENSOR, 1, 2, 0, DM_POLYTOPE_POINT_PRISM_TENSOR, 0, 0, DM_POLYTOPE_SEGMENT, 1, 0, 1, DM_POLYTOPE_SEGMENT, 1, 1, 1, DM_POLYTOPE_POINT_PRISM_TENSOR, 1, 3, 0, DM_POLYTOPE_POINT_PRISM_TENSOR, 0, 0, DM_POLYTOPE_SEGMENT, 1, 0, 2, DM_POLYTOPE_SEGMENT, 1, 1, 2, DM_POLYTOPE_POINT_PRISM_TENSOR, 1, 4, 0, DM_POLYTOPE_POINT_PRISM_TENSOR, 0, 0, DM_POLYTOPE_QUADRILATERAL, 1, 0, 0, DM_POLYTOPE_QUADRILATERAL, 1, 1, 0, DM_POLYTOPE_SEG_PRISM_TENSOR, 1, 2, 0, DM_POLYTOPE_SEG_PRISM_TENSOR, 0, 0, DM_POLYTOPE_SEG_PRISM_TENSOR, 0, 2, DM_POLYTOPE_SEG_PRISM_TENSOR, 1, 4, 1, DM_POLYTOPE_QUADRILATERAL, 1, 0, 1, DM_POLYTOPE_QUADRILATERAL, 1, 1, 1, DM_POLYTOPE_SEG_PRISM_TENSOR, 1, 2, 1, DM_POLYTOPE_SEG_PRISM_TENSOR, 1, 3, 0, DM_POLYTOPE_SEG_PRISM_TENSOR, 0, 1, DM_POLYTOPE_SEG_PRISM_TENSOR, 0, 0, DM_POLYTOPE_QUADRILATERAL, 1, 0, 2, DM_POLYTOPE_QUADRILATERAL, 1, 1, 2, DM_POLYTOPE_SEG_PRISM_TENSOR, 0, 2, DM_POLYTOPE_SEG_PRISM_TENSOR, 0, 1, DM_POLYTOPE_SEG_PRISM_TENSOR, 1, 3, 1, DM_POLYTOPE_SEG_PRISM_TENSOR, 1, 4, 0};
316: static PetscInt ttriO[] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 0, 0, 0, 0, -1, 0, 0, 0, -1, 0, 0};
317: /* TODO Add 3 quads inside every tensor triangular prism, making 4 new triangular prisms.
318: 2
319: |\
320: | \
321: | \
322: 0---1
324: 2
326: 0 1
328: 2
329: |\
330: | \
331: | \
332: 0---1
333: */
334: static DMPolytopeType ctripT[] = {DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL, DM_POLYTOPE_HEXAHEDRON};
335: static PetscInt ctripS[] = {1, 3, 3};
336: static PetscInt ctripC[] = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 1, 1, 0, DM_POLYTOPE_SEGMENT, 1, 0, 0, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 1, 1, 0, DM_POLYTOPE_SEGMENT, 1, 2, 0, DM_POLYTOPE_SEGMENT, 1, 0, 1, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 1, 1, 1, DM_POLYTOPE_SEGMENT, 1, 3, 0, DM_POLYTOPE_SEGMENT, 1, 0, 2, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 1, 1, 2, DM_POLYTOPE_SEGMENT, 1, 4, 0, DM_POLYTOPE_QUADRILATERAL, 1, 0, 0, DM_POLYTOPE_QUADRILATERAL, 1, 1, 0, DM_POLYTOPE_QUADRILATERAL, 1, 2, 0, DM_POLYTOPE_QUADRILATERAL, 0, 2, DM_POLYTOPE_QUADRILATERAL, 0, 0, DM_POLYTOPE_QUADRILATERAL, 1, 4, 1, DM_POLYTOPE_QUADRILATERAL, 1, 0, 1, DM_POLYTOPE_QUADRILATERAL, 1, 1, 1, DM_POLYTOPE_QUADRILATERAL, 1, 2, 1, DM_POLYTOPE_QUADRILATERAL, 0, 1, DM_POLYTOPE_QUADRILATERAL, 1, 3, 0, DM_POLYTOPE_QUADRILATERAL, 0, 0, DM_POLYTOPE_QUADRILATERAL, 1, 0, 2, DM_POLYTOPE_QUADRILATERAL, 1, 1, 2, DM_POLYTOPE_QUADRILATERAL, 0, 2, DM_POLYTOPE_QUADRILATERAL, 1, 3, 1, DM_POLYTOPE_QUADRILATERAL, 0, 1, DM_POLYTOPE_QUADRILATERAL, 1, 4, 0};
337: static PetscInt ctripO[] = {0, 0, 0, 0, -1, -1, 0, 0, -1, -1, 0, 0, -1, -1, -2, 0, 0, -3, 0, 1, -2, 0, 0, 0, 0, -2, -2, 0, 0, 0, -3, 1};
338: static DMPolytopeType tquadpT[] = {DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL, DM_POLYTOPE_HEXAHEDRON};
339: static PetscInt tquadpS[] = {1, 4, 4};
340: static PetscInt tquadpC[] = {
341: DM_POLYTOPE_POINT,
342: 1,
343: 0,
344: 0,
345: DM_POLYTOPE_POINT,
346: 1,
347: 1,
348: 0,
349: DM_POLYTOPE_SEGMENT,
350: 1,
351: 0,
352: 0,
353: DM_POLYTOPE_SEGMENT,
354: 0,
355: 0,
356: DM_POLYTOPE_SEGMENT,
357: 1,
358: 1,
359: 0,
360: DM_POLYTOPE_SEGMENT,
361: 1,
362: 2,
363: 0,
364: DM_POLYTOPE_SEGMENT,
365: 1,
366: 0,
367: 1,
368: DM_POLYTOPE_SEGMENT,
369: 0,
370: 0,
371: DM_POLYTOPE_SEGMENT,
372: 1,
373: 1,
374: 1,
375: DM_POLYTOPE_SEGMENT,
376: 1,
377: 3,
378: 0,
379: DM_POLYTOPE_SEGMENT,
380: 1,
381: 0,
382: 2,
383: DM_POLYTOPE_SEGMENT,
384: 0,
385: 0,
386: DM_POLYTOPE_SEGMENT,
387: 1,
388: 1,
389: 2,
390: DM_POLYTOPE_SEGMENT,
391: 1,
392: 4,
393: 0,
394: DM_POLYTOPE_SEGMENT,
395: 1,
396: 0,
397: 3,
398: DM_POLYTOPE_SEGMENT,
399: 0,
400: 0,
401: DM_POLYTOPE_SEGMENT,
402: 1,
403: 1,
404: 3,
405: DM_POLYTOPE_SEGMENT,
406: 1,
407: 5,
408: 0,
409: DM_POLYTOPE_QUADRILATERAL,
410: 1,
411: 0,
412: 0,
413: DM_POLYTOPE_QUADRILATERAL,
414: 1,
415: 1,
416: 0,
417: DM_POLYTOPE_QUADRILATERAL,
418: 1,
419: 2,
420: 0,
421: DM_POLYTOPE_QUADRILATERAL,
422: 0,
423: 3,
424: DM_POLYTOPE_QUADRILATERAL,
425: 0,
426: 0,
427: DM_POLYTOPE_QUADRILATERAL,
428: 1,
429: 5,
430: 1,
431: DM_POLYTOPE_QUADRILATERAL,
432: 1,
433: 0,
434: 1,
435: DM_POLYTOPE_QUADRILATERAL,
436: 1,
437: 1,
438: 1,
439: DM_POLYTOPE_QUADRILATERAL,
440: 1,
441: 2,
442: 1,
443: DM_POLYTOPE_QUADRILATERAL,
444: 0,
445: 1,
446: DM_POLYTOPE_QUADRILATERAL,
447: 1,
448: 3,
449: 0,
450: DM_POLYTOPE_QUADRILATERAL,
451: 0,
452: 0,
453: DM_POLYTOPE_QUADRILATERAL,
454: 1,
455: 0,
456: 2,
457: DM_POLYTOPE_QUADRILATERAL,
458: 1,
459: 1,
460: 2,
461: DM_POLYTOPE_QUADRILATERAL,
462: 0,
463: 1,
464: DM_POLYTOPE_QUADRILATERAL,
465: 1,
466: 4,
467: 0,
468: DM_POLYTOPE_QUADRILATERAL,
469: 1,
470: 3,
471: 1,
472: DM_POLYTOPE_QUADRILATERAL,
473: 0,
474: 2,
475: DM_POLYTOPE_QUADRILATERAL,
476: 1,
477: 0,
478: 3,
479: DM_POLYTOPE_QUADRILATERAL,
480: 1,
481: 1,
482: 3,
483: DM_POLYTOPE_QUADRILATERAL,
484: 0,
485: 3,
486: DM_POLYTOPE_QUADRILATERAL,
487: 1,
488: 4,
489: 1,
490: DM_POLYTOPE_QUADRILATERAL,
491: 0,
492: 2,
493: DM_POLYTOPE_QUADRILATERAL,
494: 1,
495: 5,
496: 0,
497: };
498: static PetscInt tquadpO[] = {0, 0, 0, 0, -1, -1, 0, 0, -1, -1, 0, 0, -1, -1, 0, 0, -1, -1, -2, 0, 0, -3, 0, 1, -2, 0, 0, 0, 0, -2, -2, 0, -3, 0, 0, 1, -2, 0, 0, 0, -3, 1};
500: PetscFunctionBeginHot;
501: if (rt) *rt = 0;
502: if (convertTensor) {
503: switch (source) {
504: case DM_POLYTOPE_POINT:
505: case DM_POLYTOPE_SEGMENT:
506: case DM_POLYTOPE_QUADRILATERAL:
507: case DM_POLYTOPE_HEXAHEDRON:
508: PetscCall(DMPlexTransformCellRefine_Regular(tr, source, p, rt, Nt, target, size, cone, ornt));
509: break;
510: case DM_POLYTOPE_POINT_PRISM_TENSOR:
511: *Nt = 1;
512: *target = tedgeT;
513: *size = tedgeS;
514: *cone = tedgeC;
515: *ornt = tedgeO;
516: break;
517: case DM_POLYTOPE_SEG_PRISM_TENSOR:
518: *Nt = 2;
519: *target = tsegT;
520: *size = tsegS;
521: *cone = tsegC;
522: *ornt = tsegO;
523: break;
524: case DM_POLYTOPE_TRI_PRISM_TENSOR:
525: *Nt = 3;
526: *target = ctripT;
527: *size = ctripS;
528: *cone = ctripC;
529: *ornt = ctripO;
530: break;
531: case DM_POLYTOPE_TRIANGLE:
532: *Nt = 3;
533: *target = triT;
534: *size = triS;
535: *cone = triC;
536: *ornt = triO;
537: break;
538: case DM_POLYTOPE_TETRAHEDRON:
539: *Nt = 4;
540: *target = tetT;
541: *size = tetS;
542: *cone = tetC;
543: *ornt = tetO;
544: break;
545: case DM_POLYTOPE_TRI_PRISM:
546: *Nt = 4;
547: *target = tripT;
548: *size = tripS;
549: *cone = tripC;
550: *ornt = tripO;
551: break;
552: case DM_POLYTOPE_QUAD_PRISM_TENSOR:
553: *Nt = 3;
554: *target = tquadpT;
555: *size = tquadpS;
556: *cone = tquadpC;
557: *ornt = tquadpO;
558: break;
559: case DM_POLYTOPE_PYRAMID:
560: *Nt = 0;
561: *target = NULL;
562: *size = NULL;
563: *cone = NULL;
564: *ornt = NULL;
565: break;
566: default:
567: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No refinement strategy for %s", DMPolytopeTypes[source]);
568: }
569: } else {
570: switch (source) {
571: case DM_POLYTOPE_POINT:
572: case DM_POLYTOPE_POINT_PRISM_TENSOR:
573: case DM_POLYTOPE_SEGMENT:
574: case DM_POLYTOPE_QUADRILATERAL:
575: case DM_POLYTOPE_SEG_PRISM_TENSOR:
576: case DM_POLYTOPE_HEXAHEDRON:
577: case DM_POLYTOPE_QUAD_PRISM_TENSOR:
578: PetscCall(DMPlexTransformCellRefine_Regular(tr, source, p, rt, Nt, target, size, cone, ornt));
579: break;
580: case DM_POLYTOPE_TRIANGLE:
581: *Nt = 3;
582: *target = triT;
583: *size = triS;
584: *cone = triC;
585: *ornt = triO;
586: break;
587: case DM_POLYTOPE_TETRAHEDRON:
588: *Nt = 4;
589: *target = tetT;
590: *size = tetS;
591: *cone = tetC;
592: *ornt = tetO;
593: break;
594: case DM_POLYTOPE_TRI_PRISM:
595: *Nt = 4;
596: *target = tripT;
597: *size = tripS;
598: *cone = tripC;
599: *ornt = tripO;
600: break;
601: case DM_POLYTOPE_TRI_PRISM_TENSOR:
602: *Nt = 3;
603: *target = ttriT;
604: *size = ttriS;
605: *cone = ttriC;
606: *ornt = ttriO;
607: break;
608: case DM_POLYTOPE_PYRAMID:
609: *Nt = 0;
610: *target = NULL;
611: *size = NULL;
612: *cone = NULL;
613: *ornt = NULL;
614: break;
615: default:
616: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No refinement strategy for %s", DMPolytopeTypes[source]);
617: }
618: }
619: PetscFunctionReturn(PETSC_SUCCESS);
620: }
622: static PetscErrorCode DMPlexTransformInitialize_ToBox(DMPlexTransform tr)
623: {
624: PetscFunctionBegin;
625: tr->ops->view = DMPlexTransformView_ToBox;
626: tr->ops->setup = DMPlexTransformSetUp_ToBox;
627: tr->ops->destroy = DMPlexTransformDestroy_ToBox;
628: tr->ops->setdimensions = DMPlexTransformSetDimensions_Internal;
629: tr->ops->celltransform = DMPlexTransformCellRefine_ToBox;
630: tr->ops->getsubcellorientation = DMPlexTransformGetSubcellOrientation_ToBox;
631: tr->ops->mapcoordinates = DMPlexTransformMapCoordinatesBarycenter_Internal;
632: PetscFunctionReturn(PETSC_SUCCESS);
633: }
635: PETSC_EXTERN PetscErrorCode DMPlexTransformCreate_ToBox(DMPlexTransform tr)
636: {
637: DMPlexRefine_ToBox *f;
639: PetscFunctionBegin;
641: PetscCall(PetscNew(&f));
642: tr->data = f;
644: PetscCall(DMPlexTransformInitialize_ToBox(tr));
645: PetscFunctionReturn(PETSC_SUCCESS);
646: }