Actual source code: plexref1d.c
1: #include <petsc/private/dmplextransformimpl.h>
3: static PetscErrorCode DMPlexTransformSetUp_1D(DMPlexTransform tr)
4: {
5: DM dm;
6: DMLabel active;
7: PetscInt pStart, pEnd, p;
9: DMPlexTransformGetDM(tr, &dm);
10: DMPlexTransformGetActive(tr, &active);
12: /* Calculate refineType for each cell */
13: DMLabelCreate(PETSC_COMM_SELF, "Refine Type", &tr->trType);
14: DMPlexGetChart(dm, &pStart, &pEnd);
15: for (p = pStart; p < pEnd; ++p) {
16: DMLabel trType = tr->trType;
17: DMPolytopeType ct;
18: PetscInt val;
20: DMPlexGetCellType(dm, p, &ct);
21: switch (ct) {
22: case DM_POLYTOPE_POINT: DMLabelSetValue(trType, p, 0); break;
23: case DM_POLYTOPE_SEGMENT:
24: case DM_POLYTOPE_POINT_PRISM_TENSOR:
25: DMLabelGetValue(active, p, &val);
26: if (val == 1) DMLabelSetValue(trType, p, val);
27: else DMLabelSetValue(trType, p, 2);
28: break;
29: default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Cannot handle points of type %s", DMPolytopeTypes[ct]);
30: }
31: }
32: return 0;
33: }
35: static PetscErrorCode DMPlexTransformGetSubcellOrientation_1D(DMPlexTransform tr, DMPolytopeType sct, PetscInt sp, PetscInt so, DMPolytopeType tct, PetscInt r, PetscInt o, PetscInt *rnew, PetscInt *onew)
36: {
37: PetscInt rt;
40: DMLabelGetValue(tr->trType, sp, &rt);
41: *rnew = r; *onew = o;
42: switch (rt) {
43: case 1:
44: DMPlexTransformGetSubcellOrientation_Regular(tr, sct, sp, so, tct, r, o, rnew, onew);
45: break;
46: default: DMPlexTransformGetSubcellOrientationIdentity(tr, sct, sp, so, tct, r, o, rnew, onew);
47: }
48: return 0;
49: }
51: static PetscErrorCode DMPlexTransformCellTransform_1D(DMPlexTransform tr, DMPolytopeType source, PetscInt p, PetscInt *rt, PetscInt *Nt, DMPolytopeType *target[], PetscInt *size[], PetscInt *cone[], PetscInt *ornt[])
52: {
53: DMLabel trType = tr->trType;
54: PetscInt val;
58: DMLabelGetValue(trType, p, &val);
59: if (rt) *rt = val;
60: switch (source) {
61: case DM_POLYTOPE_POINT:
62: DMPlexTransformCellRefine_Regular(tr, source, p, NULL, Nt, target, size, cone, ornt);
63: break;
64: case DM_POLYTOPE_POINT_PRISM_TENSOR:
65: case DM_POLYTOPE_SEGMENT:
66: if (val == 1) DMPlexTransformCellRefine_Regular(tr, source, p, NULL, Nt, target, size, cone, ornt);
67: else DMPlexTransformCellTransformIdentity(tr, source, p, NULL, Nt, target, size, cone, ornt);
68: break;
69: default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No refinement strategy for %s", DMPolytopeTypes[source]);
70: }
71: return 0;
72: }
74: static PetscErrorCode DMPlexTransformSetFromOptions_1D(PetscOptionItems *PetscOptionsObject, DMPlexTransform tr)
75: {
76: PetscInt cells[256], n = 256, i;
77: PetscBool flg;
80: PetscOptionsHead(PetscOptionsObject,"DMPlex Options");
81: PetscOptionsIntArray("-dm_plex_transform_1d_ref_cell", "Mark cells for refinement", "", cells, &n, &flg);
82: if (flg) {
83: DMLabel active;
85: DMLabelCreate(PETSC_COMM_SELF, "Adaptation Label", &active);
86: for (i = 0; i < n; ++i) DMLabelSetValue(active, cells[i], DM_ADAPT_REFINE);
87: DMPlexTransformSetActive(tr, active);
88: DMLabelDestroy(&active);
89: }
90: PetscOptionsTail();
91: return 0;
92: }
94: static PetscErrorCode DMPlexTransformView_1D(DMPlexTransform tr, PetscViewer viewer)
95: {
96: PetscBool isascii;
100: PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &isascii);
101: if (isascii) {
102: PetscViewerFormat format;
103: const char *name;
105: PetscObjectGetName((PetscObject) tr, &name);
106: PetscViewerASCIIPrintf(viewer, "1D refinement %s\n", name ? name : "");
107: PetscViewerGetFormat(viewer, &format);
108: if (format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
109: DMLabelView(tr->trType, viewer);
110: }
111: } else {
112: SETERRQ(PetscObjectComm((PetscObject) tr), PETSC_ERR_SUP, "Viewer type %s not yet supported for DMPlexTransform writing", ((PetscObject) viewer)->type_name);
113: }
114: return 0;
115: }
117: static PetscErrorCode DMPlexTransformDestroy_1D(DMPlexTransform tr)
118: {
119: PetscFree(tr->data);
120: return 0;
121: }
123: static PetscErrorCode DMPlexTransformInitialize_1D(DMPlexTransform tr)
124: {
125: tr->ops->view = DMPlexTransformView_1D;
126: tr->ops->setfromoptions = DMPlexTransformSetFromOptions_1D;
127: tr->ops->setup = DMPlexTransformSetUp_1D;
128: tr->ops->destroy = DMPlexTransformDestroy_1D;
129: tr->ops->celltransform = DMPlexTransformCellTransform_1D;
130: tr->ops->getsubcellorientation = DMPlexTransformGetSubcellOrientation_1D;
131: tr->ops->mapcoordinates = DMPlexTransformMapCoordinatesBarycenter_Internal;
132: return 0;
133: }
135: PETSC_EXTERN PetscErrorCode DMPlexTransformCreate_1D(DMPlexTransform tr)
136: {
137: DMPlexRefine_1D *f;
140: PetscNewLog(tr, &f);
141: tr->data = f;
143: DMPlexTransformInitialize_1D(tr);
144: return 0;
145: }