Actual source code: plexreftobox.c

  1: #include <petsc/private/dmplextransformimpl.h>

  3: static PetscErrorCode DMPlexTransformView_ToBox(DMPlexTransform tr, PetscViewer viewer)
  4: {
  5:   PetscBool      isascii;

  9:   PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &isascii);
 10:   if (isascii) {
 11:     const char *name;

 13:     PetscObjectGetName((PetscObject) tr, &name);
 14:     PetscViewerASCIIPrintf(viewer, "Transformation to box cells %s\n", name ? name : "");
 15:   } else {
 16:     SETERRQ(PetscObjectComm((PetscObject) tr), PETSC_ERR_SUP, "Viewer type %s not yet supported for DMPlexTransform writing", ((PetscObject) viewer)->type_name);
 17:   }
 18:   return 0;
 19: }

 21: static PetscErrorCode DMPlexTransformSetUp_ToBox(DMPlexTransform tr)
 22: {
 23:   return 0;
 24: }

 26: static PetscErrorCode DMPlexTransformDestroy_ToBox(DMPlexTransform tr)
 27: {
 28:   DMPlexRefine_ToBox *f = (DMPlexRefine_ToBox *) tr->data;

 30:   PetscFree(f);
 31:   return 0;
 32: }

 34: static PetscErrorCode DMPlexTransformGetSubcellOrientation_ToBox(DMPlexTransform tr, DMPolytopeType sct, PetscInt sp, PetscInt so, DMPolytopeType tct, PetscInt r, PetscInt o, PetscInt *rnew, PetscInt *onew)
 35: {
 36:   PetscBool      convertTensor = PETSC_TRUE;
 37:   static PetscInt tri_seg[]  = {0, 0, 2, 0, 1, 0,
 38:                                 2, 0, 1, 0, 0, 0,
 39:                                 1, 0, 0, 0, 2, 0,
 40:                                 0, 0, 1, 0, 2, 0,
 41:                                 1, 0, 2, 0, 0, 0,
 42:                                 2, 0, 0, 0, 1, 0};
 43:   static PetscInt tri_quad[] = {1, -3, 0, -3, 2, -4,
 44:                                 0, -2, 2, -2, 1, -2,
 45:                                 2, -1, 1, -4, 0, -1,
 46:                                 0,  0, 1,  0, 2,  0,
 47:                                 1,  1, 2,  2, 0,  1,
 48:                                 2,  3, 0,  3, 1,  2};
 49:   static PetscInt tseg_seg[]  = {0, -1,
 50:                                  0,  0,
 51:                                  0,  0,
 52:                                  0, -1};
 53:   static PetscInt tseg_quad[] = {1,  2, 0,  2,
 54:                                  1, -3, 0, -3,
 55:                                  0,  0, 1,  0,
 56:                                  0, -1, 1, -1};
 57:   static PetscInt tet_seg[]  = {3, 0, 2, 0, 0, 0, 1, 0,
 58:                                 3, 0, 1, 0, 2, 0, 0, 0,
 59:                                 3, 0, 0, 0, 1, 0, 2, 0,
 60:                                 2, 0, 3, 0, 1, 0, 0, 0,
 61:                                 2, 0, 0, 0, 3, 0, 1, 0,
 62:                                 2, 0, 1, 0, 0, 0, 3, 0,
 63:                                 1, 0, 0, 0, 2, 0, 3, 0,
 64:                                 1, 0, 3, 0, 0, 0, 2, 0,
 65:                                 1, 0, 2, 0, 3, 0, 0, 0,
 66:                                 0, 0, 1, 0, 3, 0, 2, 0,
 67:                                 0, 0, 2, 0, 1, 0, 3, 0,
 68:                                 0, 0, 3, 0, 2, 0, 1, 0,
 69:                                 0, 0, 1, 0, 2, 0, 3, 0,
 70:                                 0, 0, 3, 0, 1, 0, 2, 0,
 71:                                 0, 0, 2, 0, 3, 0, 1, 0,
 72:                                 1, 0, 0, 0, 3, 0, 2, 0,
 73:                                 1, 0, 2, 0, 0, 0, 3, 0,
 74:                                 1, 0, 3, 0, 2, 0, 0, 0,
 75:                                 2, 0, 3, 0, 0, 0, 1, 0,
 76:                                 2, 0, 1, 0, 3, 0, 0, 0,
 77:                                 2, 0, 0, 0, 1, 0, 3, 0,
 78:                                 3, 0, 2, 0, 1, 0, 0, 0,
 79:                                 3, 0, 0, 0, 2, 0, 1, 0,
 80:                                 3, 0, 1, 0, 0, 0, 2, 0};
 81:   static PetscInt tet_quad[] = {2, 0, 5, -3, 4, 0, 0, 3, 3, 1, 1, 1,
 82:                                 0, 0, 3, 0, 5, 0, 1, 0, 4, -2, 2, 0,
 83:                                 1, 1, 4, -3, 3, -3, 2, 3, 5, -2, 0, 0,
 84:                                 3, 1, 5, 3, 0, 0, 4, 3, 2, 0, 1, -3,
 85:                                 4, 0, 2, 3, 5, -2, 1, -4, 0, -2, 3, 1,
 86:                                 1, -3, 0, -3, 2, -2, 3, 0, 5, 0, 4, 0,
 87:                                 2, -2, 1, -4, 0, -2, 4, -3, 3, -3, 5, 0,
 88:                                 4, -2, 3, -4, 1, 1, 5, 3, 0, 0, 2, -2,
 89:                                 5, 0, 0, 3, 3, 1, 2, -3, 1, -3, 4, -2,
 90:                                 3, -3, 1, 0, 4, -2, 0, -3, 2, -2, 5, -2,
 91:                                 0, -2, 2, -3, 1, -3, 5, -3, 4, 0, 3, -3,
 92:                                 5, -2, 4, 3, 2, 0, 3, -4, 1, 1, 0, -2,
 93:                                 0, 0, 1, 0, 2, 0, 3, 0, 4, 0, 5, 0,
 94:                                 3, 1, 4, 3, 1, -3, 5, 3, 2, -2, 0, 0,
 95:                                 5, 0, 2, -3, 4, -2, 0, 3, 1, 1, 3, 1,
 96:                                 4, 0, 1, -4, 3, 1, 2, 3, 0, 0, 5, -2,
 97:                                 2, 0, 0, 3, 1, 1, 5, -3, 3, -3, 4, 0,
 98:                                 5, -2, 3, -4, 0, -2, 4, 3, 1, -3, 2, 0,
 99:                                 4, -2, 5, 3, 2, -2, 3, -4, 0, -2, 1, 1,
100:                                 3, -3, 0, -3, 5, -2, 1, 0, 2, 0, 4, -2,
101:                                 1, 1, 2, 3, 0, 0, 4, -3, 5, 0, 3, -3,
102:                                 0, -2, 5, -3, 3, -3, 2, -3, 4, -2, 1, -3,
103:                                 2, -2, 4, -3, 5, 0, 1, -4, 3, 1, 0, -2,
104:                                 1, -3, 3, 0, 4, 0, 0, -3, 5, -2, 2, -2};
105:   static PetscInt tet_hex[]  = {2, -2, 3, -2, 1, -10, 0, -13,
106:                                 3, -10, 1, -13, 2, -10, 0, -10,
107:                                 1, -2, 2, -13, 3, -13, 0, -2,
108:                                 3, -13, 2, -10, 0, -2, 1, -2,
109:                                 2, -13, 0, -10, 3, -2, 1, -13,
110:                                 0, -13, 3, -10, 2, -2, 1, -10,
111:                                 0, -10, 1, -2, 3, -10, 2, -10,
112:                                 1, -10, 3, -13, 0, -13, 2, -2,
113:                                 3, -2, 0, -2, 1, -13, 2, -13,
114:                                 1, -13, 0, -13, 2, -13, 3, -2,
115:                                 0, -2, 2, -2, 1, -2, 3, -13,
116:                                 2, -10, 1, -10, 0, -10, 3, -10,
117:                                 0, 0, 1, 0, 2, 0, 3, 0,
118:                                 1, 17, 2, 17, 0, 17, 3, 16,
119:                                 2, 16, 0, 16, 1, 16, 3, 17,
120:                                 1, 16, 0, 17, 3, 17, 2, 16,
121:                                 0, 16, 3, 16, 1, 0, 2, 17,
122:                                 3, 0, 1, 17, 0, 0, 2, 0,
123:                                 2, 17, 3, 0, 0, 16, 1, 0,
124:                                 3, 17, 0, 0, 2, 16, 1, 16,
125:                                 0, 17, 2, 0, 3, 16, 1, 17,
126:                                 3, 16, 2, 16, 1, 17, 0, 17,
127:                                 2, 0, 1, 16, 3, 0, 0, 0,
128:                                 1, 0, 3, 17, 2, 17, 0, 16};
129:   static PetscInt trip_seg[]  = {1, 0, 0, 0, 3, 0, 4, 0, 2, 0,
130:                                  1, 0, 0, 0, 4, 0, 2, 0, 3, 0,
131:                                  1, 0, 0, 0, 2, 0, 3, 0, 4, 0,
132:                                  0, 0, 1, 0, 3, 0, 2, 0, 4, 0,
133:                                  0, 0, 1, 0, 4, 0, 3, 0, 2, 0,
134:                                  0, 0, 1, 0, 2, 0, 4, 0, 3, 0,
135:                                  0, 0, 1, 0, 2, 0, 3, 0, 4, 0,
136:                                  0, 0, 1, 0, 4, 0, 2, 0, 3, 0,
137:                                  0, 0, 1, 0, 3, 0, 4, 0, 2, 0,
138:                                  1, 0, 0, 0, 2, 0, 4, 0, 3, 0,
139:                                  1, 0, 0, 0, 4, 0, 3, 0, 2, 0,
140:                                  1, 0, 0, 0, 3, 0, 2, 0, 4, 0};
141:   static PetscInt trip_quad[] = {1, 1, 2, 2, 0, 1, 7, -1, 8, -1, 6, -1, 4, -1, 5, -1, 3, -1,
142:                                  2, 3, 0, 3, 1, 2, 8, -1, 6, -1, 7, -1, 5, -1, 3, -1, 4, -1,
143:                                  0, 0, 1, 0, 2, 0, 6, -1, 7, -1, 8, -1, 3, -1, 4, -1, 5, -1,
144:                                  2, -1, 1, -4, 0, -1, 4, 0, 3, 0, 5, 0, 7, 0, 6, 0, 8, 0,
145:                                  0, -2, 2, -2, 1, -2, 5, 0, 4, 0, 3, 0, 8, 0, 7, 0, 6, 0,
146:                                  1, -3, 0, -3, 2, -4, 3, 0, 5, 0, 4, 0, 6, 0, 8, 0, 7, 0,
147:                                  0, 0, 1, 0, 2, 0, 3, 0, 4, 0, 5, 0, 6, 0, 7, 0, 8, 0,
148:                                  2, 3, 0, 3, 1, 2, 5, 0, 3, 0, 4, 0, 8, 0, 6, 0, 7, 0,
149:                                  1, 1, 2, 2, 0, 1, 4, 0, 5, 0, 3, 0, 7, 0, 8, 0, 6, 0,
150:                                  1, -3, 0, -3, 2, -4, 6, -1, 8, -1, 7, -1, 3, -1, 5, -1, 4, -1,
151:                                  0, -2, 2, -2, 1, -2, 8, -1, 7, -1, 6, -1, 5, -1, 4, -1, 3, -1,
152:                                  2, -1, 1, -4, 0, -1, 7, -1, 6, -1, 8, -1, 4, -1, 3, -1, 5, -1};
153:   static PetscInt trip_hex[]  = {4, -12, 5, -6, 3, -12, 1, -12, 2, -6, 0, -12,
154:                                  5, -11, 3, -11, 4, -6, 2, -11, 0, -11, 1, -6,
155:                                  3, -9, 4, -9, 5, -9, 0, -9, 1, -9, 2, -9,
156:                                  2, -3, 1, -4, 0, -3, 5, -3, 4, -4, 3, -3,
157:                                  0, -2, 2, -2, 1, -2, 3, -2, 5, -2, 4, -2,
158:                                  1, -1, 0, -1, 2, -4, 4, -1, 3, -1, 5, -4,
159:                                  0, 0, 1, 0, 2, 0, 3, 0, 4, 0, 5, 0,
160:                                  2, 1, 0, 1, 1, 2, 5, 1, 3, 1, 4, 2,
161:                                  1, 3, 2, 2, 0, 3, 4, 3, 5, 2, 3, 3,
162:                                  4, 8, 3, 8, 5, 11, 1, 8, 0, 8, 2, 11,
163:                                  3, 10, 5, 10, 4, 10, 0, 10, 2, 10, 1, 10,
164:                                  5, 5, 4, 11, 3, 5, 2, 5, 1, 11, 0, 5};
165:   static PetscInt ctrip_seg[]  = {0, -1,
166:                                   0, -1,
167:                                   0, -1,
168:                                   0, 0,
169:                                   0, 0,
170:                                   0, 0,
171:                                   0, 0,
172:                                   0, 0,
173:                                   0, 0,
174:                                   0, -1,
175:                                   0, -1,
176:                                   0, -1};
177:   static PetscInt ctrip_quad[] = {0, -1, 2, -1, 1, -1,
178:                                   2, -1, 1, -1, 0, -1,
179:                                   1, -1, 0, -1, 2, -1,
180:                                   0, 0, 2, 0, 1, 0,
181:                                   2, 0, 1, 0, 0, 0,
182:                                   1, 0, 0, 0, 2, 0,
183:                                   0, 0, 1, 0, 2, 0,
184:                                   1, 0, 2, 0, 0, 0,
185:                                   2, 0, 0, 0, 1, 0,
186:                                   0, -1, 1, -1, 2, -1,
187:                                   1, -1, 2, -1, 0, -1,
188:                                   2, -1, 0, -1, 1, -1};
189:   static PetscInt ctrip_hex[]  = {1, 8, 0, 8, 2, 11,
190:                                   0, 10, 2, 10, 1, 10,
191:                                   2, 5, 1, 11, 0, 5,
192:                                   1, -1, 0, -1, 2, -4,
193:                                   0, -2, 2, -2, 1, -2,
194:                                   2, -3, 1, -4, 0, -3,
195:                                   0, 0, 1, 0, 2, 0,
196:                                   1, 3, 2, 2, 0, 3,
197:                                   2, 1, 0, 1, 1, 2,
198:                                   0, -9, 1, -9, 2, -9,
199:                                   1, -12, 2, -6, 0, -12,
200:                                   2, -11, 0, -11, 1, -6};
201:   static PetscInt tquadp_seg[]  = {0, -1,
202:                                    0, -1,
203:                                    0, -1,
204:                                    0, -1,
205:                                    0,  0,
206:                                    0,  0,
207:                                    0,  0,
208:                                    0,  0,
209:                                    0,  0,
210:                                    0,  0,
211:                                    0,  0,
212:                                    0,  0,
213:                                    0, -1,
214:                                    0, -1,
215:                                    0, -1,
216:                                    0, -1};
217:   static PetscInt tquadp_quad[] = {1, -1, 0, -1, 3, -1, 2, -1,
218:                                    0, -1, 3, -1, 2, -1, 1, -1,
219:                                    3, -1, 2, -1, 1, -1, 0, -1,
220:                                    2, -1, 1, -1, 0, -1, 3, -1,
221:                                    1,  0, 0,  0, 3,  0, 2,  0,
222:                                    0,  0, 3,  0, 2,  0, 1,  0,
223:                                    3,  0, 2,  0, 1,  0, 0,  0,
224:                                    2,  0, 1,  0, 0,  0, 3,  0,
225:                                    0,  0, 1,  0, 2,  0, 3,  0,
226:                                    1,  0, 2,  0, 3,  0, 0,  0,
227:                                    2,  0, 3,  0, 0,  0, 1,  0,
228:                                    3,  0, 0,  0, 1,  0, 2,  0,
229:                                    0, -1, 1, -1, 2, -1, 3, -1,
230:                                    1, -1, 2, -1, 3, -1, 0, -1,
231:                                    2, -1, 3, -1, 0, -1, 1, -1,
232:                                    3, -1, 0, -1, 1, -1, 2, -1};
233:   static PetscInt tquadp_hex[]  = {2, 11,  1,  11, 0,  11, 3,  11,
234:                                    1,  8,  0,   8, 3,   8, 2,   8,
235:                                    0, 10,  3,  10, 2,  10, 1,  10,
236:                                    3,  5,  2,   5, 1,   5, 0,   5,
237:                                    2, -4,  1,  -4, 0,  -4, 3,  -4,
238:                                    1, -1,  0,  -1, 3,  -1, 2,  -1,
239:                                    0, -2,  3,  -2, 2,  -2, 1,  -2,
240:                                    3, -3,  2,  -3, 1,  -3, 0,  -3,
241:                                    0,   0, 1,   0, 2,   0, 3,   0,
242:                                    1,   3, 2,   3, 3,   3, 0,   3,
243:                                    2,   2, 3,   2, 0,   2, 1,   2,
244:                                    3,   1, 0,   1, 1,   1, 2,   1,
245:                                    0,  -9, 1,  -9, 2,  -9, 3,  -9,
246:                                    1, -12, 2, -12, 3, -12, 0, -12,
247:                                    2,  -6, 3,  -6, 0,  -6, 1,  -6,
248:                                    3, -11, 0, -11, 1, -11, 2, -11};

251:   *rnew = r; *onew = o;
252:   if (!so) return 0;
253:   if (convertTensor) {
254:     switch (sct) {
255:       case DM_POLYTOPE_POINT:
256:       case DM_POLYTOPE_SEGMENT:
257:       case DM_POLYTOPE_POINT_PRISM_TENSOR:
258:       case DM_POLYTOPE_QUADRILATERAL:
259:       case DM_POLYTOPE_HEXAHEDRON:
260:         DMPlexTransformGetSubcellOrientation_Regular(tr, sct, sp, so, tct, r, o, rnew, onew);
261:         break;
262:       case DM_POLYTOPE_TRIANGLE:
263:       switch (tct) {
264:         case DM_POLYTOPE_POINT: break;
265:         case DM_POLYTOPE_SEGMENT:
266:           *rnew = tri_seg[(so+3)*6 + r*2];
267:           *onew = DMPolytopeTypeComposeOrientation(tct, o, tri_seg[(so+3)*6 + r*2 + 1]);
268:           break;
269:         case DM_POLYTOPE_QUADRILATERAL:
270:           *rnew = tri_quad[(so+3)*6 + r*2];
271:           *onew = DMPolytopeTypeComposeOrientation(tct, o, tri_quad[(so+3)*6 + r*2 + 1]);
272:           break;
273:         default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
274:       }
275:       break;
276:       case DM_POLYTOPE_SEG_PRISM_TENSOR:
277:       switch (tct) {
278:         case DM_POLYTOPE_SEGMENT:
279:         case DM_POLYTOPE_POINT_PRISM_TENSOR:
280:           *rnew = tseg_seg[(so+2)*2 + r*2];
281:           *onew = DMPolytopeTypeComposeOrientation(tct, o, tseg_seg[(so+2)*2 + r*2 + 1]);
282:           break;
283:         case DM_POLYTOPE_QUADRILATERAL:
284:           *rnew = tseg_quad[(so+2)*4 + r*2];
285:           *onew = DMPolytopeTypeComposeOrientation(tct, o, tseg_quad[(so+2)*4 + r*2 + 1]);
286:           break;
287:         default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
288:       }
289:       break;
290:       case DM_POLYTOPE_TETRAHEDRON:
291:       switch (tct) {
292:         case DM_POLYTOPE_POINT: break;
293:         case DM_POLYTOPE_SEGMENT:
294:           *rnew = tet_seg[(so+12)*8 + r*2];
295:           *onew = DMPolytopeTypeComposeOrientation(tct, o, tet_seg[(so+12)*8 + r*2 + 1]);
296:           break;
297:         case DM_POLYTOPE_QUADRILATERAL:
298:           *rnew = tet_quad[(so+12)*12 + r*2];
299:           *onew = DMPolytopeTypeComposeOrientation(tct, o, tet_quad[(so+12)*12 + r*2 + 1]);
300:           break;
301:         case DM_POLYTOPE_HEXAHEDRON:
302:           *rnew = tet_hex[(so+12)*8 + r*2];
303:           *onew = DMPolytopeTypeComposeOrientation(tct, o, tet_hex[(so+12)*8 + r*2 + 1]);
304:           break;
305:         default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
306:       }
307:       break;
308:       case DM_POLYTOPE_TRI_PRISM:
309:       switch (tct) {
310:         case DM_POLYTOPE_POINT: break;
311:         case DM_POLYTOPE_SEGMENT:
312:           *rnew = trip_seg[(so+6)*10 + r*2];
313:           *onew = DMPolytopeTypeComposeOrientation(tct, o, trip_seg[(so+6)*10 + r*2 + 1]);
314:           break;
315:         case DM_POLYTOPE_QUADRILATERAL:
316:           *rnew = trip_quad[(so+6)*18 + r*2];
317:           *onew = DMPolytopeTypeComposeOrientation(tct, o, trip_quad[(so+6)*18 + r*2 + 1]);
318:           break;
319:         case DM_POLYTOPE_HEXAHEDRON:
320:           *rnew = trip_hex[(so+6)*12 + r*2];
321:           *onew = DMPolytopeTypeComposeOrientation(tct, o, trip_hex[(so+6)*12 + r*2 + 1]);
322:           break;
323:         default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
324:       }
325:       break;
326:       case DM_POLYTOPE_TRI_PRISM_TENSOR:
327:       switch (tct) {
328:         case DM_POLYTOPE_SEGMENT:
329:           *rnew = ctrip_seg[(so+6)*2 + r*2];
330:           *onew = DMPolytopeTypeComposeOrientation(tct, o, ctrip_seg[(so+6)*2 + r*2 + 1]);
331:           break;
332:         case DM_POLYTOPE_QUADRILATERAL:
333:           *rnew = ctrip_quad[(so+6)*6 + r*2];
334:           *onew = DMPolytopeTypeComposeOrientation(tct, o, ctrip_quad[(so+6)*6 + r*2 + 1]);
335:           break;
336:         case DM_POLYTOPE_HEXAHEDRON:
337:           *rnew = ctrip_hex[(so+6)*6 + r*2];
338:           *onew = DMPolytopeTypeComposeOrientation(tct, o, ctrip_hex[(so+6)*6 + r*2 + 1]);
339:           break;
340:         default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
341:       }
342:       break;
343:       case DM_POLYTOPE_QUAD_PRISM_TENSOR:
344:       switch (tct) {
345:         case DM_POLYTOPE_SEGMENT:
346:           *rnew = tquadp_seg[(so+8)*2 + r*2];
347:           *onew = DMPolytopeTypeComposeOrientation(tct, o, tquadp_seg[(so+8)*2 + r*2 + 1]);
348:           break;
349:         case DM_POLYTOPE_QUADRILATERAL:
350:           *rnew = tquadp_quad[(so+8)*8 + r*2];
351:           *onew = DMPolytopeTypeComposeOrientation(tct, o, tquadp_quad[(so+8)*8 + r*2 + 1]);
352:           break;
353:         case DM_POLYTOPE_HEXAHEDRON:
354:           *rnew = tquadp_hex[(so+8)*8 + r*2];
355:           *onew = DMPolytopeTypeComposeOrientation(tct, o, tquadp_hex[(so+8)*8 + r*2 + 1]);
356:           break;
357:         default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
358:       }
359:       break;
360:       default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported cell type %s", DMPolytopeTypes[sct]);
361:     }
362:   } else {
363:     switch (sct) {
364:       case DM_POLYTOPE_POINT:
365:       case DM_POLYTOPE_SEGMENT:
366:       case DM_POLYTOPE_POINT_PRISM_TENSOR:
367:       case DM_POLYTOPE_QUADRILATERAL:
368:       case DM_POLYTOPE_SEG_PRISM_TENSOR:
369:       case DM_POLYTOPE_HEXAHEDRON:
370:       case DM_POLYTOPE_QUAD_PRISM_TENSOR:
371:         DMPlexTransformGetSubcellOrientation_Regular(tr, sct, sp, so, tct, r, o, rnew, onew);
372:         break;
373:       default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported cell type %s", DMPolytopeTypes[sct]);
374:     }
375:   }
376:   return 0;
377: }

379: static PetscErrorCode DMPlexTransformCellRefine_ToBox(DMPlexTransform tr, DMPolytopeType source, PetscInt p, PetscInt *rt, PetscInt *Nt, DMPolytopeType *target[], PetscInt *size[], PetscInt *cone[], PetscInt *ornt[])
380: {
381:   PetscBool      convertTensor = PETSC_TRUE;
382:   /* Change tensor edges to segments */
383:   static DMPolytopeType tedgeT[]  = {DM_POLYTOPE_SEGMENT};
384:   static PetscInt       tedgeS[]  = {1};
385:   static PetscInt       tedgeC[]  = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 1, 1, 0};
386:   static PetscInt       tedgeO[]  = {                         0,                          0};
387:   /* Add 1 vertex, 3 edges inside every triangle, making 3 new quadrilaterals.
388:    2
389:    |\
390:    | \
391:    |  \
392:    |   \
393:    0    1
394:    |     \
395:    |      \
396:    2       1
397:    |\     / \
398:    | 2   1   \
399:    |  \ /     \
400:    1   |       0
401:    |   0        \
402:    |   |         \
403:    |   |          \
404:    0-0-0-----1-----1
405:   */
406:   static DMPolytopeType triT[]    = {DM_POLYTOPE_POINT, DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL};
407:   static PetscInt       triS[]    = {1, 3, 3};
408:   static PetscInt       triC[]    = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 0,    0,
409:                                      DM_POLYTOPE_POINT, 1, 1, 0, DM_POLYTOPE_POINT, 0,    0,
410:                                      DM_POLYTOPE_POINT, 1, 2, 0, DM_POLYTOPE_POINT, 0,    0,
411:                                      DM_POLYTOPE_SEGMENT, 1, 0, 0, DM_POLYTOPE_SEGMENT, 0,    0, DM_POLYTOPE_SEGMENT, 0,    2, DM_POLYTOPE_SEGMENT, 1, 2, 1,
412:                                      DM_POLYTOPE_SEGMENT, 1, 0, 1, DM_POLYTOPE_SEGMENT, 1, 1, 0, DM_POLYTOPE_SEGMENT, 0,    1, DM_POLYTOPE_SEGMENT, 0,    0,
413:                                      DM_POLYTOPE_SEGMENT, 0,    2, DM_POLYTOPE_SEGMENT, 0,    1, DM_POLYTOPE_SEGMENT, 1, 1, 1, DM_POLYTOPE_SEGMENT, 1, 2, 0};
414:   static PetscInt       triO[]    = {0, 0,
415:                                      0, 0,
416:                                      0, 0,
417:                                      0,  0, -1,  0,
418:                                      0,  0,  0, -1,
419:                                      0, -1,  0,  0};
420:   /* Add 1 edge inside every tensor quad, making 2 new quadrilaterals
421:      2----2----1----3----3
422:      |         |         |
423:      |         |         |
424:      |         |         |
425:      4    A    6    B    5
426:      |         |         |
427:      |         |         |
428:      |         |         |
429:      0----0----0----1----1
430:   */
431:   static DMPolytopeType tsegT[]  = {DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL};
432:   static PetscInt       tsegS[]  = {1, 2};
433:   static PetscInt       tsegC[]  = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 1, 1, 0,
434:                                      /* TODO  Fix these */
435:                                      DM_POLYTOPE_SEGMENT, 1, 0, 0, DM_POLYTOPE_SEGMENT, 0,    0, DM_POLYTOPE_SEGMENT, 1, 1, 0, DM_POLYTOPE_SEGMENT, 1, 2, 0,
436:                                      DM_POLYTOPE_SEGMENT, 1, 0, 1, DM_POLYTOPE_SEGMENT, 1, 3, 0, DM_POLYTOPE_SEGMENT, 1, 1, 1, DM_POLYTOPE_SEGMENT, 0,    0};
437:   static PetscInt       tsegO[]  = {0, 0,
438:                                     0, 0, -1, -1,
439:                                     0, 0, -1, -1};
440:   /* Add 6 triangles inside every cell, making 4 new hexs
441:      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
442:      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]
443:      called e3, e4, and e5. The faces of a tet, given in DMPlexGetRawFaces_Internal() are
444:        [v0, v1, v2], [v0, v3, v1], [v0, v2, v3], [v2, v1, v3]
445:      We make a new hex in each corner
446:        [v0, (e0, 0), (f0, 0), (e2, 0), (e3, 0), (f2, 0), (c0, 0), (f1, 0)]
447:        [v1, (e4, 0), (f3, 0), (e1, 0), (e0, 0), (f0, 0), (c0, 0), (f1, 0)]
448:        [v2, (e1, 0), (f3, 0), (e5, 0), (e2, 0), (f2, 0), (c0, 0), (f0, 0)]
449:        [v3, (e4, 0), (f1, 0), (e3, 0), (e5, 0), (f2, 0), (c0, 0), (f3, 0)]
450:      We create a new face for each edge
451:        [(e3, 0), (f2, 0), (c0, 0), (f1, 0)]
452:        [(f0, 0), (e0, 0), (f1, 0), (c0, 0)]
453:        [(e2, 0), (f0, 0), (c0, 0), (f2, 0)]
454:        [(f3, 0), (e4, 0), (f1, 0), (c0, 0)]
455:        [(e1, 0), (f3, 0), (c0, 0), (f0, 0)]
456:        [(e5, 0), (f3, 0), (c0, 0), (f2, 0)]
457:      I could write a program to generate these from the first hex by acting with the symmetry group to take one subcell into another.
458:    */
459:   static DMPolytopeType tetT[]    = {DM_POLYTOPE_POINT, DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL, DM_POLYTOPE_HEXAHEDRON};
460:   static PetscInt       tetS[]    = {1, 4, 6, 4};
461:   static PetscInt       tetC[]    = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 0, 0,
462:                                      DM_POLYTOPE_POINT, 1, 1, 0, DM_POLYTOPE_POINT, 0, 0,
463:                                      DM_POLYTOPE_POINT, 1, 2, 0, DM_POLYTOPE_POINT, 0, 0,
464:                                      DM_POLYTOPE_POINT, 1, 3, 0, DM_POLYTOPE_POINT, 0, 0,
465:                                      DM_POLYTOPE_SEGMENT, 1, 2, 2, DM_POLYTOPE_SEGMENT, 0,    2, DM_POLYTOPE_SEGMENT, 0,    1, DM_POLYTOPE_SEGMENT, 1, 1, 0,
466:                                      DM_POLYTOPE_SEGMENT, 1, 0, 0, DM_POLYTOPE_SEGMENT, 1, 1, 2, DM_POLYTOPE_SEGMENT, 0,    1, DM_POLYTOPE_SEGMENT, 0,    0,
467:                                      DM_POLYTOPE_SEGMENT, 1, 0, 2, DM_POLYTOPE_SEGMENT, 0,    0, DM_POLYTOPE_SEGMENT, 0,    2, DM_POLYTOPE_SEGMENT, 1, 2, 0,
468:                                      DM_POLYTOPE_SEGMENT, 1, 3, 1, DM_POLYTOPE_SEGMENT, 1, 1, 1, DM_POLYTOPE_SEGMENT, 0,    1, DM_POLYTOPE_SEGMENT, 0,    3,
469:                                      DM_POLYTOPE_SEGMENT, 1, 3, 0, DM_POLYTOPE_SEGMENT, 0,    3, DM_POLYTOPE_SEGMENT, 0,    0, DM_POLYTOPE_SEGMENT, 1, 0, 1,
470:                                      DM_POLYTOPE_SEGMENT, 1, 3, 2, DM_POLYTOPE_SEGMENT, 0,    3, DM_POLYTOPE_SEGMENT, 0,    2, DM_POLYTOPE_SEGMENT, 1, 2, 1,
471:                                      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,
472:                                      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,
473:                                      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,
474:                                      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};
475:   static PetscInt       tetO[]    = {0, 0,
476:                                      0, 0,
477:                                      0, 0,
478:                                      0, 0,
479:                                      0,  0, -1, -1,
480:                                     -1,  0,  0, -1,
481:                                      0,  0, -1, -1,
482:                                     -1,  0,  0, -1,
483:                                      0,  0, -1, -1,
484:                                      0,  0, -1, -1,
485:                                      0,  0,  0,  0,  0,  0,
486:                                      1, -3,  1,  0,  0,  3,
487:                                      0, -2,  1, -3,  0,  3,
488:                                      1, -2,  3, -4, -2,  3};
489:   /* Add 3 quads inside every triangular prism, making 4 new prisms. */
490:   static DMPolytopeType tripT[]   = {DM_POLYTOPE_POINT, DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL, DM_POLYTOPE_HEXAHEDRON};
491:   static PetscInt       tripS[]   = {1, 5, 9, 6};
492:   static PetscInt       tripC[]   = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 0, 0,
493:                                      DM_POLYTOPE_POINT, 1, 1, 0, DM_POLYTOPE_POINT, 0, 0,
494:                                      DM_POLYTOPE_POINT, 1, 2, 0, DM_POLYTOPE_POINT, 0, 0,
495:                                      DM_POLYTOPE_POINT, 1, 3, 0, DM_POLYTOPE_POINT, 0, 0,
496:                                      DM_POLYTOPE_POINT, 1, 4, 0, DM_POLYTOPE_POINT, 0, 0,
497:                                      DM_POLYTOPE_SEGMENT, 1, 2, 3, DM_POLYTOPE_SEGMENT, 0,    2, DM_POLYTOPE_SEGMENT, 0,    4, DM_POLYTOPE_SEGMENT, 1, 4, 1,
498:                                      DM_POLYTOPE_SEGMENT, 1, 2, 1, DM_POLYTOPE_SEGMENT, 1, 3, 3, DM_POLYTOPE_SEGMENT, 0,    3, DM_POLYTOPE_SEGMENT, 0,    2,
499:                                      DM_POLYTOPE_SEGMENT, 0,    4, DM_POLYTOPE_SEGMENT, 0,    3, DM_POLYTOPE_SEGMENT, 1, 3, 1, DM_POLYTOPE_SEGMENT, 1, 4, 3,
500:                                      DM_POLYTOPE_SEGMENT, 1, 0, 2, DM_POLYTOPE_SEGMENT, 0,    0, DM_POLYTOPE_SEGMENT, 0,    2, DM_POLYTOPE_SEGMENT, 1, 2, 0,
501:                                      DM_POLYTOPE_SEGMENT, 1, 0, 1, DM_POLYTOPE_SEGMENT, 0,    0, DM_POLYTOPE_SEGMENT, 0,    3, DM_POLYTOPE_SEGMENT, 1, 3, 0,
502:                                      DM_POLYTOPE_SEGMENT, 1, 0, 0, DM_POLYTOPE_SEGMENT, 0,    0, DM_POLYTOPE_SEGMENT, 0,    4, DM_POLYTOPE_SEGMENT, 1, 4, 0,
503:                                      DM_POLYTOPE_SEGMENT, 0,    2, DM_POLYTOPE_SEGMENT, 0,    1, DM_POLYTOPE_SEGMENT, 1, 1, 0, DM_POLYTOPE_SEGMENT, 1, 2, 2,
504:                                      DM_POLYTOPE_SEGMENT, 0,    3, DM_POLYTOPE_SEGMENT, 0,    1, DM_POLYTOPE_SEGMENT, 1, 1, 1, DM_POLYTOPE_SEGMENT, 1, 3, 2,
505:                                      DM_POLYTOPE_SEGMENT, 0,    4, DM_POLYTOPE_SEGMENT, 0,    1, DM_POLYTOPE_SEGMENT, 1, 1, 2, DM_POLYTOPE_SEGMENT, 1, 4, 2,
506:                                      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,
507:                                      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,
508:                                      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,
509:                                      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,
510:                                      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,
511:                                      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};
512:   static PetscInt       tripO[]   = {0, 0,
513:                                      0, 0,
514:                                      0, 0,
515:                                      0, 0,
516:                                      0, 0,
517:                                      0,  0, -1, -1,
518:                                     -1,  0,  0, -1,
519:                                      0, -1, -1,  0,
520:                                      0,  0, -1, -1,
521:                                      0,  0, -1, -1,
522:                                      0,  0, -1, -1,
523:                                      0, -1, -1,  0,
524:                                      0, -1, -1,  0,
525:                                      0, -1, -1,  0,
526:                                      0,  0,  0, -3,  0,  1,
527:                                      0,  0,  0,  0,  0, -2,
528:                                      0,  0,  0,  0, -3,  1,
529:                                     -2,  0,  0, -3,  0,  1,
530:                                     -2,  0,  0,  0,  0, -2,
531:                                     -2,  0,  0,  0, -3,  1};
532:   /* Add 3 tensor quads inside every tensor triangular prism, making 4 new tensor triangular prisms.
533:       2
534:       |\
535:       | \
536:       |  \
537:       0---1

539:       2

541:       0   1

543:       2
544:       |\
545:       | \
546:       |  \
547:       0---1
548:   */
549:   static DMPolytopeType ttriT[]  = {DM_POLYTOPE_POINT_PRISM_TENSOR, DM_POLYTOPE_SEG_PRISM_TENSOR, DM_POLYTOPE_QUAD_PRISM_TENSOR};
550:   static PetscInt       ttriS[]  = {1, 3, 3};
551:   static PetscInt       ttriC[]  = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 1, 1, 0,
552:                                      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,
553:                                      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,
554:                                      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,
555:                                      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,
556:                                      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,
557:                                      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};
558:   static PetscInt       ttriO[]  = {0, 0,
559:                                      0, 0, 0, 0,
560:                                      0, 0, 0, 0,
561:                                      0, 0, 0, 0,
562:                                      0, 0, 0,  0, -1, 0,
563:                                      0, 0, 0,  0,  0, -1,
564:                                      0, 0, 0, -1,  0, 0};
565:   /* TODO Add 3 quads inside every tensor triangular prism, making 4 new triangular prisms.
566:       2
567:       |\
568:       | \
569:       |  \
570:       0---1

572:       2

574:       0   1

576:       2
577:       |\
578:       | \
579:       |  \
580:       0---1
581:   */
582:   static DMPolytopeType ctripT[]  = {DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL, DM_POLYTOPE_HEXAHEDRON};
583:   static PetscInt       ctripS[]  = {1, 3, 3};
584:   static PetscInt       ctripC[]  = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 1, 1, 0,
585:                                      DM_POLYTOPE_SEGMENT, 1, 0, 0, DM_POLYTOPE_SEGMENT, 0,    0, DM_POLYTOPE_SEGMENT, 1, 1, 0, DM_POLYTOPE_SEGMENT, 1, 2, 0,
586:                                      DM_POLYTOPE_SEGMENT, 1, 0, 1, DM_POLYTOPE_SEGMENT, 0,    0, DM_POLYTOPE_SEGMENT, 1, 1, 1, DM_POLYTOPE_SEGMENT, 1, 3, 0,
587:                                      DM_POLYTOPE_SEGMENT, 1, 0, 2, DM_POLYTOPE_SEGMENT, 0,    0, DM_POLYTOPE_SEGMENT, 1, 1, 2, DM_POLYTOPE_SEGMENT, 1, 4, 0,
588:                                      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,
589:                                      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,
590:                                      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};
591:   static PetscInt       ctripO[]  = {0, 0,
592:                                      0, 0, -1, -1,
593:                                      0, 0, -1, -1,
594:                                      0, 0, -1, -1,
595:                                     -2, 0, 0, -3,  0,  1,
596:                                     -2, 0, 0,  0,  0, -2,
597:                                     -2, 0, 0,  0, -3,  1};
598:   static DMPolytopeType tquadpT[] = {DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL, DM_POLYTOPE_HEXAHEDRON};
599:   static PetscInt       tquadpS[] = {1, 4, 4};
600:   static PetscInt       tquadpC[] = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 1, 1, 0,
601:                                      DM_POLYTOPE_SEGMENT, 1, 0, 0, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 1, 1, 0, DM_POLYTOPE_SEGMENT, 1, 2, 0,
602:                                      DM_POLYTOPE_SEGMENT, 1, 0, 1, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 1, 1, 1, DM_POLYTOPE_SEGMENT, 1, 3, 0,
603:                                      DM_POLYTOPE_SEGMENT, 1, 0, 2, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 1, 1, 2, DM_POLYTOPE_SEGMENT, 1, 4, 0,
604:                                      DM_POLYTOPE_SEGMENT, 1, 0, 3, DM_POLYTOPE_SEGMENT, 0, 0, DM_POLYTOPE_SEGMENT, 1, 1, 3, DM_POLYTOPE_SEGMENT, 1, 5, 0,
605:                                      DM_POLYTOPE_QUADRILATERAL, 1, 0, 0, DM_POLYTOPE_QUADRILATERAL, 1, 1, 0, DM_POLYTOPE_QUADRILATERAL, 1, 2, 0, DM_POLYTOPE_QUADRILATERAL, 0,    3, DM_POLYTOPE_QUADRILATERAL, 0,    0, DM_POLYTOPE_QUADRILATERAL, 1, 5, 1,
606:                                      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,
607:                                      DM_POLYTOPE_QUADRILATERAL, 1, 0, 2, DM_POLYTOPE_QUADRILATERAL, 1, 1, 2, DM_POLYTOPE_QUADRILATERAL, 0,    1, DM_POLYTOPE_QUADRILATERAL, 1, 4, 0, DM_POLYTOPE_QUADRILATERAL, 1, 3, 1, DM_POLYTOPE_QUADRILATERAL, 0,    2,
608:                                      DM_POLYTOPE_QUADRILATERAL, 1, 0, 3, DM_POLYTOPE_QUADRILATERAL, 1, 1, 3, DM_POLYTOPE_QUADRILATERAL, 0,    3, DM_POLYTOPE_QUADRILATERAL, 1, 4, 1, DM_POLYTOPE_QUADRILATERAL, 0,    2, DM_POLYTOPE_QUADRILATERAL, 1, 5, 0,
609:   };
610:   static PetscInt       tquadpO[] = {0,  0,
611:                                      0,  0, -1, -1,
612:                                      0,  0, -1, -1,
613:                                      0,  0, -1, -1,
614:                                      0,  0, -1, -1,
615:                                     -2,  0,  0, -3,  0,  1,
616:                                     -2,  0,  0,  0,  0, -2,
617:                                     -2,  0, -3,  0,  0,  1,
618:                                     -2,  0,  0,  0, -3,  1};

621:   if (rt) *rt = 0;
622:   if (convertTensor) {
623:     switch (source) {
624:       case DM_POLYTOPE_POINT:
625:       case DM_POLYTOPE_SEGMENT:
626:       case DM_POLYTOPE_QUADRILATERAL:
627:       case DM_POLYTOPE_HEXAHEDRON:
628:         DMPlexTransformCellRefine_Regular(tr, source, p, rt, Nt, target, size, cone, ornt);
629:         break;
630:       case DM_POLYTOPE_POINT_PRISM_TENSOR: *Nt = 1; *target = tedgeT;  *size = tedgeS;  *cone = tedgeC;  *ornt = tedgeO;  break;
631:       case DM_POLYTOPE_SEG_PRISM_TENSOR:   *Nt = 2; *target = tsegT;   *size = tsegS;   *cone = tsegC;   *ornt = tsegO;   break;
632:       case DM_POLYTOPE_TRI_PRISM_TENSOR:   *Nt = 3; *target = ctripT;  *size = ctripS;  *cone = ctripC;  *ornt = ctripO;  break;
633:       case DM_POLYTOPE_TRIANGLE:           *Nt = 3; *target = triT;    *size = triS;    *cone = triC;    *ornt = triO;    break;
634:       case DM_POLYTOPE_TETRAHEDRON:        *Nt = 4; *target = tetT;    *size = tetS;    *cone = tetC;    *ornt = tetO;    break;
635:       case DM_POLYTOPE_TRI_PRISM:          *Nt = 4; *target = tripT;   *size = tripS;   *cone = tripC;   *ornt = tripO;   break;
636:       case DM_POLYTOPE_QUAD_PRISM_TENSOR:  *Nt = 3; *target = tquadpT; *size = tquadpS; *cone = tquadpC; *ornt = tquadpO; break;
637:       case DM_POLYTOPE_PYRAMID:            *Nt = 0; *target = NULL;    *size = NULL;    *cone = NULL;    *ornt = NULL;    break;
638:       default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No refinement strategy for %s", DMPolytopeTypes[source]);
639:     }
640:   } else {
641:     switch (source) {
642:       case DM_POLYTOPE_POINT:
643:       case DM_POLYTOPE_POINT_PRISM_TENSOR:
644:       case DM_POLYTOPE_SEGMENT:
645:       case DM_POLYTOPE_QUADRILATERAL:
646:       case DM_POLYTOPE_SEG_PRISM_TENSOR:
647:       case DM_POLYTOPE_HEXAHEDRON:
648:       case DM_POLYTOPE_QUAD_PRISM_TENSOR:
649:         DMPlexTransformCellRefine_Regular(tr, source, p, rt, Nt, target, size, cone, ornt);
650:         break;
651:       case DM_POLYTOPE_TRIANGLE:           *Nt = 3; *target = triT;    *size = triS;    *cone = triC;    *ornt = triO;    break;
652:       case DM_POLYTOPE_TETRAHEDRON:        *Nt = 4; *target = tetT;    *size = tetS;    *cone = tetC;    *ornt = tetO;    break;
653:       case DM_POLYTOPE_TRI_PRISM:          *Nt = 4; *target = tripT;   *size = tripS;   *cone = tripC;   *ornt = tripO;   break;
654:       case DM_POLYTOPE_TRI_PRISM_TENSOR:   *Nt = 3; *target = ttriT;   *size = ttriS;   *cone = ttriC;   *ornt = ttriO;   break;
655:       case DM_POLYTOPE_PYRAMID:            *Nt = 0; *target = NULL;    *size = NULL;    *cone = NULL;    *ornt = NULL;    break;
656:       default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No refinement strategy for %s", DMPolytopeTypes[source]);
657:     }
658:   }
659:   return 0;
660: }

662: static PetscErrorCode DMPlexTransformInitialize_ToBox(DMPlexTransform tr)
663: {
664:   tr->ops->view    = DMPlexTransformView_ToBox;
665:   tr->ops->setup   = DMPlexTransformSetUp_ToBox;
666:   tr->ops->destroy = DMPlexTransformDestroy_ToBox;
667:   tr->ops->celltransform = DMPlexTransformCellRefine_ToBox;
668:   tr->ops->getsubcellorientation = DMPlexTransformGetSubcellOrientation_ToBox;
669:   tr->ops->mapcoordinates = DMPlexTransformMapCoordinatesBarycenter_Internal;
670:   return 0;
671: }

673: PETSC_EXTERN PetscErrorCode DMPlexTransformCreate_ToBox(DMPlexTransform tr)
674: {
675:   DMPlexRefine_ToBox *f;

678:   PetscNewLog(tr, &f);
679:   tr->data = f;

681:   DMPlexTransformInitialize_ToBox(tr);
682:   return 0;
683: }