Actual source code: plexreftobox.c

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

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

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

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

 23: static PetscErrorCode DMPlexTransformSetUp_ToBox(DMPlexTransform tr)
 24: {
 26:   return(0);
 27: }

 29: static PetscErrorCode DMPlexTransformDestroy_ToBox(DMPlexTransform tr)
 30: {
 31:   DMPlexRefine_ToBox *f = (DMPlexRefine_ToBox *) tr->data;
 32:   PetscErrorCode          ierr;

 35:   PetscFree(f);
 36:   return(0);
 37: }

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

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

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

546:       2

548:       0   1

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

579:       2

581:       0   1

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

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

669: static PetscErrorCode DMPlexTransformInitialize_ToBox(DMPlexTransform tr)
670: {
672:   tr->ops->view    = DMPlexTransformView_ToBox;
673:   tr->ops->setup   = DMPlexTransformSetUp_ToBox;
674:   tr->ops->destroy = DMPlexTransformDestroy_ToBox;
675:   tr->ops->celltransform = DMPlexTransformCellRefine_ToBox;
676:   tr->ops->getsubcellorientation = DMPlexTransformGetSubcellOrientation_ToBox;
677:   tr->ops->mapcoordinates = DMPlexTransformMapCoordinatesBarycenter_Internal;
678:   return(0);
679: }

681: PETSC_EXTERN PetscErrorCode DMPlexTransformCreate_ToBox(DMPlexTransform tr)
682: {
683:   DMPlexRefine_ToBox *f;
684:   PetscErrorCode          ierr;

688:   PetscNewLog(tr, &f);
689:   tr->data = f;

691:   DMPlexTransformInitialize_ToBox(tr);
692:   return(0);
693: }