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, 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};
 38:   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};
 39:   static PetscInt tseg_seg[]    = {0, -1, 0, 0, 0, 0, 0, -1};
 40:   static PetscInt tseg_quad[]   = {1, 2, 0, 2, 1, -3, 0, -3, 0, 0, 1, 0, 0, -1, 1, -1};
 41:   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,
 42:                                    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,
 43:                                    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};
 44:   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,
 45:                                    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,
 46:                                    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,
 47:                                    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,
 48:                                    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};
 49:   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,
 50:                                    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,
 51:                                    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,
 52:                                    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};
 53:   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,
 54:                                    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};
 55:   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,
 56:                                    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,
 57:                                    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,
 58:                                    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};
 59:   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,
 60:                                    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,
 61:                                    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};
 62:   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};
 63:   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,
 64:                                    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};
 65:   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,
 66:                                    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};
 67:   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};
 68:   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,
 69:                                    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,
 70:                                    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};
 71:   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,
 72:                                    -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,
 73:                                    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};

 76:   *rnew = r;
 77:   *onew = o;
 78:   if (!so) return 0;
 79:   if (convertTensor) {
 80:     switch (sct) {
 81:     case DM_POLYTOPE_POINT:
 82:     case DM_POLYTOPE_SEGMENT:
 83:     case DM_POLYTOPE_POINT_PRISM_TENSOR:
 84:     case DM_POLYTOPE_QUADRILATERAL:
 85:     case DM_POLYTOPE_HEXAHEDRON:
 86:       DMPlexTransformGetSubcellOrientation_Regular(tr, sct, sp, so, tct, r, o, rnew, onew);
 87:       break;
 88:     case DM_POLYTOPE_TRIANGLE:
 89:       switch (tct) {
 90:       case DM_POLYTOPE_POINT:
 91:         break;
 92:       case DM_POLYTOPE_SEGMENT:
 93:         *rnew = tri_seg[(so + 3) * 6 + r * 2];
 94:         *onew = DMPolytopeTypeComposeOrientation(tct, o, tri_seg[(so + 3) * 6 + r * 2 + 1]);
 95:         break;
 96:       case DM_POLYTOPE_QUADRILATERAL:
 97:         *rnew = tri_quad[(so + 3) * 6 + r * 2];
 98:         *onew = DMPolytopeTypeComposeOrientation(tct, o, tri_quad[(so + 3) * 6 + r * 2 + 1]);
 99:         break;
100:       default:
101:         SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
102:       }
103:       break;
104:     case DM_POLYTOPE_SEG_PRISM_TENSOR:
105:       switch (tct) {
106:       case DM_POLYTOPE_SEGMENT:
107:       case DM_POLYTOPE_POINT_PRISM_TENSOR:
108:         *rnew = tseg_seg[(so + 2) * 2 + r * 2];
109:         *onew = DMPolytopeTypeComposeOrientation(tct, o, tseg_seg[(so + 2) * 2 + r * 2 + 1]);
110:         break;
111:       case DM_POLYTOPE_QUADRILATERAL:
112:         *rnew = tseg_quad[(so + 2) * 4 + r * 2];
113:         *onew = DMPolytopeTypeComposeOrientation(tct, o, tseg_quad[(so + 2) * 4 + r * 2 + 1]);
114:         break;
115:       default:
116:         SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
117:       }
118:       break;
119:     case DM_POLYTOPE_TETRAHEDRON:
120:       switch (tct) {
121:       case DM_POLYTOPE_POINT:
122:         break;
123:       case DM_POLYTOPE_SEGMENT:
124:         *rnew = tet_seg[(so + 12) * 8 + r * 2];
125:         *onew = DMPolytopeTypeComposeOrientation(tct, o, tet_seg[(so + 12) * 8 + r * 2 + 1]);
126:         break;
127:       case DM_POLYTOPE_QUADRILATERAL:
128:         *rnew = tet_quad[(so + 12) * 12 + r * 2];
129:         *onew = DMPolytopeTypeComposeOrientation(tct, o, tet_quad[(so + 12) * 12 + r * 2 + 1]);
130:         break;
131:       case DM_POLYTOPE_HEXAHEDRON:
132:         *rnew = tet_hex[(so + 12) * 8 + r * 2];
133:         *onew = DMPolytopeTypeComposeOrientation(tct, o, tet_hex[(so + 12) * 8 + r * 2 + 1]);
134:         break;
135:       default:
136:         SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
137:       }
138:       break;
139:     case DM_POLYTOPE_TRI_PRISM:
140:       switch (tct) {
141:       case DM_POLYTOPE_POINT:
142:         break;
143:       case DM_POLYTOPE_SEGMENT:
144:         *rnew = trip_seg[(so + 6) * 10 + r * 2];
145:         *onew = DMPolytopeTypeComposeOrientation(tct, o, trip_seg[(so + 6) * 10 + r * 2 + 1]);
146:         break;
147:       case DM_POLYTOPE_QUADRILATERAL:
148:         *rnew = trip_quad[(so + 6) * 18 + r * 2];
149:         *onew = DMPolytopeTypeComposeOrientation(tct, o, trip_quad[(so + 6) * 18 + r * 2 + 1]);
150:         break;
151:       case DM_POLYTOPE_HEXAHEDRON:
152:         *rnew = trip_hex[(so + 6) * 12 + r * 2];
153:         *onew = DMPolytopeTypeComposeOrientation(tct, o, trip_hex[(so + 6) * 12 + r * 2 + 1]);
154:         break;
155:       default:
156:         SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
157:       }
158:       break;
159:     case DM_POLYTOPE_TRI_PRISM_TENSOR:
160:       switch (tct) {
161:       case DM_POLYTOPE_SEGMENT:
162:         *rnew = ctrip_seg[(so + 6) * 2 + r * 2];
163:         *onew = DMPolytopeTypeComposeOrientation(tct, o, ctrip_seg[(so + 6) * 2 + r * 2 + 1]);
164:         break;
165:       case DM_POLYTOPE_QUADRILATERAL:
166:         *rnew = ctrip_quad[(so + 6) * 6 + r * 2];
167:         *onew = DMPolytopeTypeComposeOrientation(tct, o, ctrip_quad[(so + 6) * 6 + r * 2 + 1]);
168:         break;
169:       case DM_POLYTOPE_HEXAHEDRON:
170:         *rnew = ctrip_hex[(so + 6) * 6 + r * 2];
171:         *onew = DMPolytopeTypeComposeOrientation(tct, o, ctrip_hex[(so + 6) * 6 + r * 2 + 1]);
172:         break;
173:       default:
174:         SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
175:       }
176:       break;
177:     case DM_POLYTOPE_QUAD_PRISM_TENSOR:
178:       switch (tct) {
179:       case DM_POLYTOPE_SEGMENT:
180:         *rnew = tquadp_seg[(so + 8) * 2 + r * 2];
181:         *onew = DMPolytopeTypeComposeOrientation(tct, o, tquadp_seg[(so + 8) * 2 + r * 2 + 1]);
182:         break;
183:       case DM_POLYTOPE_QUADRILATERAL:
184:         *rnew = tquadp_quad[(so + 8) * 8 + r * 2];
185:         *onew = DMPolytopeTypeComposeOrientation(tct, o, tquadp_quad[(so + 8) * 8 + r * 2 + 1]);
186:         break;
187:       case DM_POLYTOPE_HEXAHEDRON:
188:         *rnew = tquadp_hex[(so + 8) * 8 + r * 2];
189:         *onew = DMPolytopeTypeComposeOrientation(tct, o, tquadp_hex[(so + 8) * 8 + r * 2 + 1]);
190:         break;
191:       default:
192:         SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type %s is not produced by %s", DMPolytopeTypes[tct], DMPolytopeTypes[sct]);
193:       }
194:       break;
195:     default:
196:       SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported cell type %s", DMPolytopeTypes[sct]);
197:     }
198:   } else {
199:     switch (sct) {
200:     case DM_POLYTOPE_POINT:
201:     case DM_POLYTOPE_SEGMENT:
202:     case DM_POLYTOPE_POINT_PRISM_TENSOR:
203:     case DM_POLYTOPE_QUADRILATERAL:
204:     case DM_POLYTOPE_SEG_PRISM_TENSOR:
205:     case DM_POLYTOPE_HEXAHEDRON:
206:     case DM_POLYTOPE_QUAD_PRISM_TENSOR:
207:       DMPlexTransformGetSubcellOrientation_Regular(tr, sct, sp, so, tct, r, o, rnew, onew);
208:       break;
209:     default:
210:       SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported cell type %s", DMPolytopeTypes[sct]);
211:     }
212:   }
213:   return 0;
214: }

216: static PetscErrorCode DMPlexTransformCellRefine_ToBox(DMPlexTransform tr, DMPolytopeType source, PetscInt p, PetscInt *rt, PetscInt *Nt, DMPolytopeType *target[], PetscInt *size[], PetscInt *cone[], PetscInt *ornt[])
217: {
218:   PetscBool convertTensor = PETSC_TRUE;
219:   /* Change tensor edges to segments */
220:   static DMPolytopeType tedgeT[] = {DM_POLYTOPE_SEGMENT};
221:   static PetscInt       tedgeS[] = {1};
222:   static PetscInt       tedgeC[] = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 1, 1, 0};
223:   static PetscInt       tedgeO[] = {0, 0};
224:   /* Add 1 vertex, 3 edges inside every triangle, making 3 new quadrilaterals.
225:    2
226:    |\
227:    | \
228:    |  \
229:    |   \
230:    0    1
231:    |     \
232:    |      \
233:    2       1
234:    |\     / \
235:    | 2   1   \
236:    |  \ /     \
237:    1   |       0
238:    |   0        \
239:    |   |         \
240:    |   |          \
241:    0-0-0-----1-----1
242:   */
243:   static DMPolytopeType triT[] = {DM_POLYTOPE_POINT, DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL};
244:   static PetscInt       triS[] = {1, 3, 3};
245:   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};
246:   static PetscInt triO[] = {0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 0, 0, -1, 0, -1, 0, 0};
247:   /* Add 1 edge inside every tensor quad, making 2 new quadrilaterals
248:      2----2----1----3----3
249:      |         |         |
250:      |         |         |
251:      |         |         |
252:      4    A    6    B    5
253:      |         |         |
254:      |         |         |
255:      |         |         |
256:      0----0----0----1----1
257:   */
258:   static DMPolytopeType tsegT[] = {DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL};
259:   static PetscInt       tsegS[] = {1, 2};
260:   static PetscInt       tsegC[] = {DM_POLYTOPE_POINT, 1, 0, 0, DM_POLYTOPE_POINT, 1, 1, 0,
261:                                    /* TODO  Fix these */
262:                                    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};
263:   static PetscInt       tsegO[] = {0, 0, 0, 0, -1, -1, 0, 0, -1, -1};
264:   /* Add 6 triangles inside every cell, making 4 new hexs
265:      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
266:      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]
267:      called e3, e4, and e5. The faces of a tet, given in DMPlexGetRawFaces_Internal() are
268:        [v0, v1, v2], [v0, v3, v1], [v0, v2, v3], [v2, v1, v3]
269:      We make a new hex in each corner
270:        [v0, (e0, 0), (f0, 0), (e2, 0), (e3, 0), (f2, 0), (c0, 0), (f1, 0)]
271:        [v1, (e4, 0), (f3, 0), (e1, 0), (e0, 0), (f0, 0), (c0, 0), (f1, 0)]
272:        [v2, (e1, 0), (f3, 0), (e5, 0), (e2, 0), (f2, 0), (c0, 0), (f0, 0)]
273:        [v3, (e4, 0), (f1, 0), (e3, 0), (e5, 0), (f2, 0), (c0, 0), (f3, 0)]
274:      We create a new face for each edge
275:        [(e3, 0), (f2, 0), (c0, 0), (f1, 0)]
276:        [(f0, 0), (e0, 0), (f1, 0), (c0, 0)]
277:        [(e2, 0), (f0, 0), (c0, 0), (f2, 0)]
278:        [(f3, 0), (e4, 0), (f1, 0), (c0, 0)]
279:        [(e1, 0), (f3, 0), (c0, 0), (f0, 0)]
280:        [(e5, 0), (f3, 0), (c0, 0), (f2, 0)]
281:      I could write a program to generate these from the first hex by acting with the symmetry group to take one subcell into another.
282:    */
283:   static DMPolytopeType tetT[] = {DM_POLYTOPE_POINT, DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL, DM_POLYTOPE_HEXAHEDRON};
284:   static PetscInt       tetS[] = {1, 4, 6, 4};
285:   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};
286:   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};
287:   /* Add 3 quads inside every triangular prism, making 4 new prisms. */
288:   static DMPolytopeType tripT[] = {DM_POLYTOPE_POINT, DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL, DM_POLYTOPE_HEXAHEDRON};
289:   static PetscInt       tripS[] = {1, 5, 9, 6};
290:   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};
291:   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,
292:                              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};
293:   /* Add 3 tensor quads inside every tensor triangular prism, making 4 new tensor triangular prisms.
294:       2
295:       |\
296:       | \
297:       |  \
298:       0---1

300:       2

302:       0   1

304:       2
305:       |\
306:       | \
307:       |  \
308:       0---1
309:   */
310:   static DMPolytopeType ttriT[] = {DM_POLYTOPE_POINT_PRISM_TENSOR, DM_POLYTOPE_SEG_PRISM_TENSOR, DM_POLYTOPE_QUAD_PRISM_TENSOR};
311:   static PetscInt       ttriS[] = {1, 3, 3};
312:   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};
313:   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};
314:   /* TODO Add 3 quads inside every tensor triangular prism, making 4 new triangular prisms.
315:       2
316:       |\
317:       | \
318:       |  \
319:       0---1

321:       2

323:       0   1

325:       2
326:       |\
327:       | \
328:       |  \
329:       0---1
330:   */
331:   static DMPolytopeType ctripT[] = {DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL, DM_POLYTOPE_HEXAHEDRON};
332:   static PetscInt       ctripS[] = {1, 3, 3};
333:   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};
334:   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};
335:   static DMPolytopeType tquadpT[] = {DM_POLYTOPE_SEGMENT, DM_POLYTOPE_QUADRILATERAL, DM_POLYTOPE_HEXAHEDRON};
336:   static PetscInt       tquadpS[] = {1, 4, 4};
337:   static PetscInt       tquadpC[] = {
338:     DM_POLYTOPE_POINT,
339:     1,
340:     0,
341:     0,
342:     DM_POLYTOPE_POINT,
343:     1,
344:     1,
345:     0,
346:     DM_POLYTOPE_SEGMENT,
347:     1,
348:     0,
349:     0,
350:     DM_POLYTOPE_SEGMENT,
351:     0,
352:     0,
353:     DM_POLYTOPE_SEGMENT,
354:     1,
355:     1,
356:     0,
357:     DM_POLYTOPE_SEGMENT,
358:     1,
359:     2,
360:     0,
361:     DM_POLYTOPE_SEGMENT,
362:     1,
363:     0,
364:     1,
365:     DM_POLYTOPE_SEGMENT,
366:     0,
367:     0,
368:     DM_POLYTOPE_SEGMENT,
369:     1,
370:     1,
371:     1,
372:     DM_POLYTOPE_SEGMENT,
373:     1,
374:     3,
375:     0,
376:     DM_POLYTOPE_SEGMENT,
377:     1,
378:     0,
379:     2,
380:     DM_POLYTOPE_SEGMENT,
381:     0,
382:     0,
383:     DM_POLYTOPE_SEGMENT,
384:     1,
385:     1,
386:     2,
387:     DM_POLYTOPE_SEGMENT,
388:     1,
389:     4,
390:     0,
391:     DM_POLYTOPE_SEGMENT,
392:     1,
393:     0,
394:     3,
395:     DM_POLYTOPE_SEGMENT,
396:     0,
397:     0,
398:     DM_POLYTOPE_SEGMENT,
399:     1,
400:     1,
401:     3,
402:     DM_POLYTOPE_SEGMENT,
403:     1,
404:     5,
405:     0,
406:     DM_POLYTOPE_QUADRILATERAL,
407:     1,
408:     0,
409:     0,
410:     DM_POLYTOPE_QUADRILATERAL,
411:     1,
412:     1,
413:     0,
414:     DM_POLYTOPE_QUADRILATERAL,
415:     1,
416:     2,
417:     0,
418:     DM_POLYTOPE_QUADRILATERAL,
419:     0,
420:     3,
421:     DM_POLYTOPE_QUADRILATERAL,
422:     0,
423:     0,
424:     DM_POLYTOPE_QUADRILATERAL,
425:     1,
426:     5,
427:     1,
428:     DM_POLYTOPE_QUADRILATERAL,
429:     1,
430:     0,
431:     1,
432:     DM_POLYTOPE_QUADRILATERAL,
433:     1,
434:     1,
435:     1,
436:     DM_POLYTOPE_QUADRILATERAL,
437:     1,
438:     2,
439:     1,
440:     DM_POLYTOPE_QUADRILATERAL,
441:     0,
442:     1,
443:     DM_POLYTOPE_QUADRILATERAL,
444:     1,
445:     3,
446:     0,
447:     DM_POLYTOPE_QUADRILATERAL,
448:     0,
449:     0,
450:     DM_POLYTOPE_QUADRILATERAL,
451:     1,
452:     0,
453:     2,
454:     DM_POLYTOPE_QUADRILATERAL,
455:     1,
456:     1,
457:     2,
458:     DM_POLYTOPE_QUADRILATERAL,
459:     0,
460:     1,
461:     DM_POLYTOPE_QUADRILATERAL,
462:     1,
463:     4,
464:     0,
465:     DM_POLYTOPE_QUADRILATERAL,
466:     1,
467:     3,
468:     1,
469:     DM_POLYTOPE_QUADRILATERAL,
470:     0,
471:     2,
472:     DM_POLYTOPE_QUADRILATERAL,
473:     1,
474:     0,
475:     3,
476:     DM_POLYTOPE_QUADRILATERAL,
477:     1,
478:     1,
479:     3,
480:     DM_POLYTOPE_QUADRILATERAL,
481:     0,
482:     3,
483:     DM_POLYTOPE_QUADRILATERAL,
484:     1,
485:     4,
486:     1,
487:     DM_POLYTOPE_QUADRILATERAL,
488:     0,
489:     2,
490:     DM_POLYTOPE_QUADRILATERAL,
491:     1,
492:     5,
493:     0,
494:   };
495:   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};

498:   if (rt) *rt = 0;
499:   if (convertTensor) {
500:     switch (source) {
501:     case DM_POLYTOPE_POINT:
502:     case DM_POLYTOPE_SEGMENT:
503:     case DM_POLYTOPE_QUADRILATERAL:
504:     case DM_POLYTOPE_HEXAHEDRON:
505:       DMPlexTransformCellRefine_Regular(tr, source, p, rt, Nt, target, size, cone, ornt);
506:       break;
507:     case DM_POLYTOPE_POINT_PRISM_TENSOR:
508:       *Nt     = 1;
509:       *target = tedgeT;
510:       *size   = tedgeS;
511:       *cone   = tedgeC;
512:       *ornt   = tedgeO;
513:       break;
514:     case DM_POLYTOPE_SEG_PRISM_TENSOR:
515:       *Nt     = 2;
516:       *target = tsegT;
517:       *size   = tsegS;
518:       *cone   = tsegC;
519:       *ornt   = tsegO;
520:       break;
521:     case DM_POLYTOPE_TRI_PRISM_TENSOR:
522:       *Nt     = 3;
523:       *target = ctripT;
524:       *size   = ctripS;
525:       *cone   = ctripC;
526:       *ornt   = ctripO;
527:       break;
528:     case DM_POLYTOPE_TRIANGLE:
529:       *Nt     = 3;
530:       *target = triT;
531:       *size   = triS;
532:       *cone   = triC;
533:       *ornt   = triO;
534:       break;
535:     case DM_POLYTOPE_TETRAHEDRON:
536:       *Nt     = 4;
537:       *target = tetT;
538:       *size   = tetS;
539:       *cone   = tetC;
540:       *ornt   = tetO;
541:       break;
542:     case DM_POLYTOPE_TRI_PRISM:
543:       *Nt     = 4;
544:       *target = tripT;
545:       *size   = tripS;
546:       *cone   = tripC;
547:       *ornt   = tripO;
548:       break;
549:     case DM_POLYTOPE_QUAD_PRISM_TENSOR:
550:       *Nt     = 3;
551:       *target = tquadpT;
552:       *size   = tquadpS;
553:       *cone   = tquadpC;
554:       *ornt   = tquadpO;
555:       break;
556:     case DM_POLYTOPE_PYRAMID:
557:       *Nt     = 0;
558:       *target = NULL;
559:       *size   = NULL;
560:       *cone   = NULL;
561:       *ornt   = NULL;
562:       break;
563:     default:
564:       SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No refinement strategy for %s", DMPolytopeTypes[source]);
565:     }
566:   } else {
567:     switch (source) {
568:     case DM_POLYTOPE_POINT:
569:     case DM_POLYTOPE_POINT_PRISM_TENSOR:
570:     case DM_POLYTOPE_SEGMENT:
571:     case DM_POLYTOPE_QUADRILATERAL:
572:     case DM_POLYTOPE_SEG_PRISM_TENSOR:
573:     case DM_POLYTOPE_HEXAHEDRON:
574:     case DM_POLYTOPE_QUAD_PRISM_TENSOR:
575:       DMPlexTransformCellRefine_Regular(tr, source, p, rt, Nt, target, size, cone, ornt);
576:       break;
577:     case DM_POLYTOPE_TRIANGLE:
578:       *Nt     = 3;
579:       *target = triT;
580:       *size   = triS;
581:       *cone   = triC;
582:       *ornt   = triO;
583:       break;
584:     case DM_POLYTOPE_TETRAHEDRON:
585:       *Nt     = 4;
586:       *target = tetT;
587:       *size   = tetS;
588:       *cone   = tetC;
589:       *ornt   = tetO;
590:       break;
591:     case DM_POLYTOPE_TRI_PRISM:
592:       *Nt     = 4;
593:       *target = tripT;
594:       *size   = tripS;
595:       *cone   = tripC;
596:       *ornt   = tripO;
597:       break;
598:     case DM_POLYTOPE_TRI_PRISM_TENSOR:
599:       *Nt     = 3;
600:       *target = ttriT;
601:       *size   = ttriS;
602:       *cone   = ttriC;
603:       *ornt   = ttriO;
604:       break;
605:     case DM_POLYTOPE_PYRAMID:
606:       *Nt     = 0;
607:       *target = NULL;
608:       *size   = NULL;
609:       *cone   = NULL;
610:       *ornt   = NULL;
611:       break;
612:     default:
613:       SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No refinement strategy for %s", DMPolytopeTypes[source]);
614:     }
615:   }
616:   return 0;
617: }

619: static PetscErrorCode DMPlexTransformInitialize_ToBox(DMPlexTransform tr)
620: {
621:   tr->ops->view                  = DMPlexTransformView_ToBox;
622:   tr->ops->setup                 = DMPlexTransformSetUp_ToBox;
623:   tr->ops->destroy               = DMPlexTransformDestroy_ToBox;
624:   tr->ops->celltransform         = DMPlexTransformCellRefine_ToBox;
625:   tr->ops->getsubcellorientation = DMPlexTransformGetSubcellOrientation_ToBox;
626:   tr->ops->mapcoordinates        = DMPlexTransformMapCoordinatesBarycenter_Internal;
627:   return 0;
628: }

630: PETSC_EXTERN PetscErrorCode DMPlexTransformCreate_ToBox(DMPlexTransform tr)
631: {
632:   DMPlexRefine_ToBox *f;

635:   PetscNew(&f);
636:   tr->data = f;

638:   DMPlexTransformInitialize_ToBox(tr);
639:   return 0;
640: }