1 #include <petsc/private/dmplextransformimpl.h> /*I "petscdmplextransform.h" I*/
2
DMPlexTransformView_Filter(DMPlexTransform tr,PetscViewer viewer)3 static PetscErrorCode DMPlexTransformView_Filter(DMPlexTransform tr, PetscViewer viewer)
4 {
5 PetscBool isascii;
6
7 PetscFunctionBegin;
8 PetscValidHeaderSpecific(tr, DMPLEXTRANSFORM_CLASSID, 1);
9 PetscValidHeaderSpecific(viewer, PETSC_VIEWER_CLASSID, 2);
10 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
11 if (isascii) {
12 const char *name;
13
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 }
21
DMPlexTransformSetUp_Filter(DMPlexTransform tr)22 static PetscErrorCode DMPlexTransformSetUp_Filter(DMPlexTransform tr)
23 {
24 DM dm;
25 DMLabel active;
26 PetscInt Nc;
27
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;
35
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;
45
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 }
60
DMPlexTransformDestroy_Filter(DMPlexTransform tr)61 static PetscErrorCode DMPlexTransformDestroy_Filter(DMPlexTransform tr)
62 {
63 DMPlexTransform_Filter *f = (DMPlexTransform_Filter *)tr->data;
64
65 PetscFunctionBegin;
66 PetscCall(PetscFree(f));
67 PetscFunctionReturn(PETSC_SUCCESS);
68 }
69
DMPlexTransformCellTransform_Filter(DMPlexTransform tr,DMPolytopeType source,PetscInt p,PetscInt * rt,PetscInt * Nt,DMPolytopeType * target[],PetscInt * size[],PetscInt * cone[],PetscInt * ornt[])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;
75
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 }
91
DMPlexTransformInitialize_Filter(DMPlexTransform tr)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 }
104
DMPlexTransformCreate_Filter(DMPlexTransform tr)105 PETSC_EXTERN PetscErrorCode DMPlexTransformCreate_Filter(DMPlexTransform tr)
106 {
107 DMPlexTransform_Filter *f;
108
109 PetscFunctionBegin;
110 PetscValidHeaderSpecific(tr, DMPLEXTRANSFORM_CLASSID, 1);
111 PetscCall(PetscNew(&f));
112 tr->data = f;
113
114 PetscCall(DMPlexTransformInitialize_Filter(tr));
115 PetscFunctionReturn(PETSC_SUCCESS);
116 }
117