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: PetscFunctionBegin;
10: PetscCall(DMPlexTransformGetDM(tr, &dm));
11: PetscCall(DMPlexTransformGetActive(tr, &active));
12: PetscCheck(active, PetscObjectComm((PetscObject)tr), PETSC_ERR_ARG_WRONGSTATE, "DMPlexTransform must have an adaptation label in order to use 1D algorithm");
13: /* Calculate refineType for each cell */
14: PetscCall(DMLabelCreate(PETSC_COMM_SELF, "Refine Type", &tr->trType));
15: PetscCall(DMPlexGetChart(dm, &pStart, &pEnd));
16: for (p = pStart; p < pEnd; ++p) {
17: DMLabel trType = tr->trType;
18: DMPolytopeType ct;
19: PetscInt val;
21: PetscCall(DMPlexGetCellType(dm, p, &ct));
22: switch (ct) {
23: case DM_POLYTOPE_POINT:
24: PetscCall(DMLabelSetValue(trType, p, 0));
25: break;
26: case DM_POLYTOPE_SEGMENT:
27: case DM_POLYTOPE_POINT_PRISM_TENSOR:
28: PetscCall(DMLabelGetValue(active, p, &val));
29: if (val == 1) PetscCall(DMLabelSetValue(trType, p, val));
30: else PetscCall(DMLabelSetValue(trType, p, 2));
31: break;
32: default:
33: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Cannot handle points of type %s", DMPolytopeTypes[ct]);
34: }
35: }
36: PetscFunctionReturn(PETSC_SUCCESS);
37: }
39: static PetscErrorCode DMPlexTransformGetSubcellOrientation_1D(DMPlexTransform tr, DMPolytopeType sct, PetscInt sp, PetscInt so, DMPolytopeType tct, PetscInt r, PetscInt o, PetscInt *rnew, PetscInt *onew)
40: {
41: PetscInt rt;
43: PetscFunctionBeginHot;
44: PetscCall(DMLabelGetValue(tr->trType, sp, &rt));
45: *rnew = r;
46: *onew = o;
47: switch (rt) {
48: case 1:
49: PetscCall(DMPlexTransformGetSubcellOrientation_Regular(tr, sct, sp, so, tct, r, o, rnew, onew));
50: break;
51: default:
52: PetscCall(DMPlexTransformGetSubcellOrientationIdentity(tr, sct, sp, so, tct, r, o, rnew, onew));
53: }
54: PetscFunctionReturn(PETSC_SUCCESS);
55: }
57: static PetscErrorCode DMPlexTransformCellTransform_1D(DMPlexTransform tr, DMPolytopeType source, PetscInt p, PetscInt *rt, PetscInt *Nt, DMPolytopeType *target[], PetscInt *size[], PetscInt *cone[], PetscInt *ornt[])
58: {
59: DMLabel trType = tr->trType;
60: PetscInt val;
62: PetscFunctionBeginHot;
63: PetscCheck(p >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Point argument is invalid");
64: PetscCall(DMLabelGetValue(trType, p, &val));
65: if (rt) *rt = val;
66: switch (source) {
67: case DM_POLYTOPE_POINT:
68: PetscCall(DMPlexTransformCellRefine_Regular(tr, source, p, NULL, Nt, target, size, cone, ornt));
69: break;
70: case DM_POLYTOPE_POINT_PRISM_TENSOR:
71: case DM_POLYTOPE_SEGMENT:
72: if (val == 1) PetscCall(DMPlexTransformCellRefine_Regular(tr, source, p, NULL, Nt, target, size, cone, ornt));
73: else PetscCall(DMPlexTransformCellTransformIdentity(tr, source, p, NULL, Nt, target, size, cone, ornt));
74: break;
75: default:
76: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No refinement strategy for %s", DMPolytopeTypes[source]);
77: }
78: PetscFunctionReturn(PETSC_SUCCESS);
79: }
81: static PetscErrorCode DMPlexTransformSetFromOptions_1D(DMPlexTransform tr, PetscOptionItems *PetscOptionsObject)
82: {
83: PetscInt cells[256], n = 256, i;
84: PetscBool flg;
86: PetscFunctionBegin;
87: PetscOptionsHeadBegin(PetscOptionsObject, "DMPlex Options");
88: PetscCall(PetscOptionsIntArray("-dm_plex_transform_1d_ref_cell", "Mark cells for refinement", "", cells, &n, &flg));
89: if (flg) {
90: DMLabel active;
92: PetscCall(DMLabelCreate(PETSC_COMM_SELF, "Adaptation Label", &active));
93: for (i = 0; i < n; ++i) PetscCall(DMLabelSetValue(active, cells[i], DM_ADAPT_REFINE));
94: PetscCall(DMPlexTransformSetActive(tr, active));
95: PetscCall(DMLabelDestroy(&active));
96: }
97: PetscOptionsHeadEnd();
98: PetscFunctionReturn(PETSC_SUCCESS);
99: }
101: static PetscErrorCode DMPlexTransformView_1D(DMPlexTransform tr, PetscViewer viewer)
102: {
103: PetscBool isascii;
105: PetscFunctionBegin;
108: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
109: if (isascii) {
110: PetscViewerFormat format;
111: const char *name;
113: PetscCall(PetscObjectGetName((PetscObject)tr, &name));
114: PetscCall(PetscViewerASCIIPrintf(viewer, "1D refinement %s\n", name ? name : ""));
115: PetscCall(PetscViewerGetFormat(viewer, &format));
116: if (format == PETSC_VIEWER_ASCII_INFO_DETAIL) PetscCall(DMLabelView(tr->trType, viewer));
117: } else {
118: SETERRQ(PetscObjectComm((PetscObject)tr), PETSC_ERR_SUP, "Viewer type %s not yet supported for DMPlexTransform writing", ((PetscObject)viewer)->type_name);
119: }
120: PetscFunctionReturn(PETSC_SUCCESS);
121: }
123: static PetscErrorCode DMPlexTransformDestroy_1D(DMPlexTransform tr)
124: {
125: PetscFunctionBegin;
126: PetscCall(PetscFree(tr->data));
127: PetscFunctionReturn(PETSC_SUCCESS);
128: }
130: static PetscErrorCode DMPlexTransformInitialize_1D(DMPlexTransform tr)
131: {
132: PetscFunctionBegin;
133: tr->ops->view = DMPlexTransformView_1D;
134: tr->ops->setfromoptions = DMPlexTransformSetFromOptions_1D;
135: tr->ops->setup = DMPlexTransformSetUp_1D;
136: tr->ops->destroy = DMPlexTransformDestroy_1D;
137: tr->ops->setdimensions = DMPlexTransformSetDimensions_Internal;
138: tr->ops->celltransform = DMPlexTransformCellTransform_1D;
139: tr->ops->getsubcellorientation = DMPlexTransformGetSubcellOrientation_1D;
140: tr->ops->mapcoordinates = DMPlexTransformMapCoordinatesBarycenter_Internal;
141: PetscFunctionReturn(PETSC_SUCCESS);
142: }
144: PETSC_EXTERN PetscErrorCode DMPlexTransformCreate_1D(DMPlexTransform tr)
145: {
146: DMPlexRefine_1D *f;
148: PetscFunctionBegin;
150: PetscCall(PetscNew(&f));
151: tr->data = f;
153: PetscCall(DMPlexTransformInitialize_1D(tr));
154: PetscFunctionReturn(PETSC_SUCCESS);
155: }