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