Actual source code: plextrfilter.c
1: #include <petsc/private/dmplextransformimpl.h>
3: static PetscErrorCode DMPlexTransformView_Filter(DMPlexTransform tr, PetscViewer viewer)
4: {
5: PetscBool isascii;
7: PetscFunctionBegin;
10: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
11: if (isascii) {
12: const char *name;
14: PetscCall(PetscObjectGetName((PetscObject)tr, &name));
15: PetscCall(PetscViewerASCIIPrintf(viewer, "Filter transformation %s\n", name ? name : ""));
16: } else {
17: SETERRQ(PetscObjectComm((PetscObject)tr), PETSC_ERR_SUP, "Viewer type %s not yet supported for DMPlexTransform writing", ((PetscObject)viewer)->type_name);
18: }
19: PetscFunctionReturn(PETSC_SUCCESS);
20: }
22: static PetscErrorCode DMPlexTransformSetUp_Filter(DMPlexTransform tr)
23: {
24: DM dm;
25: DMLabel active;
26: PetscInt Nc;
28: PetscFunctionBegin;
29: PetscCall(DMPlexTransformGetDM(tr, &dm));
30: PetscCall(DMPlexTransformGetActive(tr, &active));
31: if (active) {
32: IS filterIS;
33: const PetscInt *filterCells;
34: PetscInt c;
36: PetscCall(DMLabelCreate(PETSC_COMM_SELF, "Filter Type", &tr->trType));
37: PetscCall(DMLabelGetStratumIS(active, DM_ADAPT_REFINE, &filterIS));
38: PetscCall(DMLabelGetStratumSize(active, DM_ADAPT_REFINE, &Nc));
39: if (filterIS) PetscCall(ISGetIndices(filterIS, &filterCells));
40: for (c = 0; c < Nc; ++c) {
41: const PetscInt cell = filterCells[c];
42: PetscInt *closure = NULL;
43: DMPolytopeType ct;
44: PetscInt Ncl, cl;
46: PetscCall(DMPlexGetTransitiveClosure(dm, cell, PETSC_TRUE, &Ncl, &closure));
47: for (cl = 0; cl < Ncl * 2; cl += 2) {
48: PetscCall(DMPlexGetCellType(dm, closure[cl], &ct));
49: PetscCall(DMLabelSetValue(tr->trType, closure[cl], ct));
50: }
51: PetscCall(DMPlexRestoreTransitiveClosure(dm, cell, PETSC_TRUE, &Ncl, &closure));
52: }
53: if (filterIS) {
54: PetscCall(ISRestoreIndices(filterIS, &filterCells));
55: PetscCall(ISDestroy(&filterIS));
56: }
57: }
58: PetscFunctionReturn(PETSC_SUCCESS);
59: }
61: static PetscErrorCode DMPlexTransformDestroy_Filter(DMPlexTransform tr)
62: {
63: DMPlexTransform_Filter *f = (DMPlexTransform_Filter *)tr->data;
65: PetscFunctionBegin;
66: PetscCall(PetscFree(f));
67: PetscFunctionReturn(PETSC_SUCCESS);
68: }
70: static PetscErrorCode DMPlexTransformCellTransform_Filter(DMPlexTransform tr, DMPolytopeType source, PetscInt p, PetscInt *rt, PetscInt *Nt, DMPolytopeType *target[], PetscInt *size[], PetscInt *cone[], PetscInt *ornt[])
71: {
72: PetscFunctionBeginHot;
73: if (tr->trType && p >= 0) {
74: PetscInt val;
76: PetscCall(DMLabelGetValue(tr->trType, p, &val));
77: if (val >= 0) {
78: if (rt) *rt = val;
79: PetscCall(DMPlexTransformCellTransformIdentity(tr, source, p, NULL, Nt, target, size, cone, ornt));
80: PetscFunctionReturn(PETSC_SUCCESS);
81: }
82: }
83: if (rt) *rt = -1;
84: *Nt = 0;
85: *target = NULL;
86: *size = NULL;
87: *cone = NULL;
88: *ornt = NULL;
89: PetscFunctionReturn(PETSC_SUCCESS);
90: }
92: static PetscErrorCode DMPlexTransformInitialize_Filter(DMPlexTransform tr)
93: {
94: PetscFunctionBegin;
95: tr->ops->view = DMPlexTransformView_Filter;
96: tr->ops->setup = DMPlexTransformSetUp_Filter;
97: tr->ops->destroy = DMPlexTransformDestroy_Filter;
98: tr->ops->setdimensions = DMPlexTransformSetDimensions_Internal;
99: tr->ops->celltransform = DMPlexTransformCellTransform_Filter;
100: tr->ops->getsubcellorientation = DMPlexTransformGetSubcellOrientationIdentity;
101: tr->ops->mapcoordinates = DMPlexTransformMapCoordinatesBarycenter_Internal;
102: PetscFunctionReturn(PETSC_SUCCESS);
103: }
105: PETSC_EXTERN PetscErrorCode DMPlexTransformCreate_Filter(DMPlexTransform tr)
106: {
107: DMPlexTransform_Filter *f;
109: PetscFunctionBegin;
111: PetscCall(PetscNew(&f));
112: tr->data = f;
114: PetscCall(DMPlexTransformInitialize_Filter(tr));
115: PetscFunctionReturn(PETSC_SUCCESS);
116: }