Actual source code: dmplextransformimpl.h
1: #pragma once
3: #include <petsc/private/dmpleximpl.h>
4: #include <petscdmplextransform.h>
6: PETSC_EXTERN PetscLogEvent DMPLEXTRANSFORM_SetUp;
7: PETSC_EXTERN PetscLogEvent DMPLEXTRANSFORM_Apply;
8: PETSC_EXTERN PetscLogEvent DMPLEXTRANSFORM_SetConeSizes;
9: PETSC_EXTERN PetscLogEvent DMPLEXTRANSFORM_SetCones;
10: PETSC_EXTERN PetscLogEvent DMPLEXTRANSFORM_CreateSF;
11: PETSC_EXTERN PetscLogEvent DMPLEXTRANSFORM_CreateLabels;
12: PETSC_EXTERN PetscLogEvent DMPLEXTRANSFORM_SetCoordinates;
14: typedef struct _p_DMPlexTransformOps *DMPlexTransformOps;
15: struct _p_DMPlexTransformOps {
16: PetscErrorCode (*view)(DMPlexTransform, PetscViewer);
17: PetscErrorCode (*setfromoptions)(DMPlexTransform, PetscOptionItems);
18: PetscErrorCode (*setup)(DMPlexTransform);
19: PetscErrorCode (*destroy)(DMPlexTransform);
20: PetscErrorCode (*setdimensions)(DMPlexTransform, DM, DM);
21: PetscErrorCode (*celltransform)(DMPlexTransform, DMPolytopeType, PetscInt, PetscInt *, PetscInt *, DMPolytopeType *[], PetscInt *[], PetscInt *[], PetscInt *[]);
22: PetscErrorCode (*ordersupports)(DMPlexTransform, DM, DM);
23: PetscErrorCode (*getsubcellorientation)(DMPlexTransform, DMPolytopeType, PetscInt, PetscInt, DMPolytopeType, PetscInt, PetscInt, PetscInt *, PetscInt *);
24: PetscErrorCode (*mapcoordinates)(DMPlexTransform, DMPolytopeType, DMPolytopeType, PetscInt, PetscInt, PetscInt, PetscInt, const PetscScalar[], PetscScalar[]);
25: };
27: struct _p_DMPlexTransform {
28: PETSCHEADER(struct _p_DMPlexTransformOps);
29: void *data;
31: DM dm; /* This is the DM for which the transform has been computed */
32: DMLabel active; /* If not NULL, indicates points that are participating in the transform */
33: DMLabel trType; /* If not NULL, this holds the transformation type for each point */
34: PetscBool setupcalled; /* true if setup has been called */
35: PetscInt *ctOrderOld; /* [i] = ct: An array with original cell types in depth order */
36: PetscInt *ctOrderInvOld; /* [ct] = i: An array with the ordinal numbers for each original cell type */
37: PetscInt *ctStart; /* [ct]: The number for the first cell of each polytope type in the original mesh */
38: PetscInt *ctOrderNew; /* [i] = ct: An array with produced cell types in depth order */
39: PetscInt *ctOrderInvNew; /* [ct] = i: An array with the ordinal numbers for each produced cell type */
40: PetscInt *ctStartNew; /* [ctNew]: The number for the first cell of each polytope type in the new mesh */
41: PetscInt *offset; /* [ct/rt][ctNew]: The offset from ctStartNew[ctNew] in the new point numbering of a point of type ctNew produced from an old point of type ct or refine type rt */
42: PetscInt depth; /* The depth of the transformed mesh */
43: PetscInt *depthStart; /* The starting point for each depth stratum */
44: PetscInt *depthEnd; /* The starting point for the next depth stratum */
45: PetscInt *trNv; /* The number of transformed vertices in the closure of a cell of each type */
46: PetscScalar **trVerts; /* The transformed vertex coordinates in the closure of a cell of each type */
47: PetscInt ****trSubVerts; /* The indices for vertices of subcell (rct, r) in a cell of each type */
48: PetscFE *coordFE; /* Finite element for each cell type, used for localized coordinate interpolation */
49: PetscFEGeom **refGeom; /* Geometry of the reference cell for each cell type */
50: PetscReal redFactor; // Used to scale maxCell in each direction
51: /* Label construction */
52: PetscBool labelMatchStrata; /* Flag to restrict labeled points to the same cell type as parents */
53: PetscInt labelReplicaInc; /* Multiplier to create new label values for replicas v = oldv + r * repInc */
54: };
56: typedef struct {
57: PetscInt dummy;
58: } DMPlexTransform_Filter;
60: typedef enum {
61: NORMAL_DEFAULT,
62: NORMAL_INPUT,
63: NORMAL_COMPUTE,
64: NORMAL_COMPUTE_BD
65: } PlexNormalAlg;
66: PETSC_EXTERN const char *const PlexNormalAlgs[];
68: typedef struct {
69: /* Inputs */
70: PetscInt dimEx; /* The dimension of the extruded mesh */
71: PetscInt cdim; /* The coordinate dimension of the input mesh */
72: PetscInt cdimEx; /* The coordinate dimension of the extruded mesh */
73: PetscInt layers; /* The number of extruded layers */
74: PetscReal thickness; /* The total thickness of the extruded layers */
75: PetscInt Nth; /* The number of specified thicknesses */
76: PetscReal *thicknesses; /* The input layer thicknesses */
77: PetscBool useTensor; /* Flag to create tensor cells */
78: PlexNormalAlg normalAlg; /* Algorithm to use for computing normal */
79: PetscReal normal[3]; /* Surface normal from input */
80: DM dmNormal; // DM for normal field
81: Vec vecNormal; // Normal at each vertex
82: PetscSimplePointFn *normalFunc; /* A function returning the normal at a given point */
83: PetscBool symmetric; /* Extrude layers symmetrically about the surface */
84: PetscBool periodic; /* Connect the extruded layer periodically to the beginning */
85: /* Calculated quantities */
86: PetscReal *layerPos; /* The position of each layer relative to the original surface, along the local normal direction */
87: PetscInt *Nt; /* The array of the number of target types */
88: DMPolytopeType **target; /* The array of target types */
89: PetscInt **size; /* The array of the number of each target type */
90: PetscInt **cone; /* The array of cones for each target cell */
91: PetscInt **ornt; /* The array of orientation for each target cell */
92: // Borrowed storage
93: const PetscInt *degree; // The root degree of all points in the original mesh
94: } DMPlexTransform_Extrude;
96: typedef struct {
97: PetscInt debug; // Debugging level
98: PetscBool useTensor; // Flag to create tensor cells
99: PetscReal width; // The width of a cohesive cell
100: PetscInt *Nt; // The array of the number of target types
101: DMPolytopeType **target; // The array of target types
102: PetscInt **size; // The array of the number of each target type
103: PetscInt **cone; // The array of cones for each target cell
104: PetscInt **ornt; // The array of orientation for each target cell
105: } DMPlexTransform_Cohesive;
107: typedef struct {
108: PetscInt dummy;
109: } DMPlexRefine_Regular;
111: typedef struct {
112: PetscInt dummy;
113: } DMPlexRefine_ToBox;
115: typedef struct {
116: PetscBool reflect; /* Flag to reflect the transformation */
117: } DMPlexRefine_ToSimplex;
119: typedef struct {
120: PetscInt dummy;
121: } DMPlexRefine_Alfeld;
123: typedef struct {
124: DMLabel splitPoints; /* List of edges to be bisected (1) and cells to be divided (2) */
125: PetscSection secEdgeLen; /* Section for edge length field */
126: PetscReal *edgeLen; /* Storage for edge length field */
127: } DMPlexRefine_SBR;
129: typedef struct {
130: PetscInt dummy;
131: } DMPlexRefine_1D;
133: typedef struct {
134: PetscInt n; /* The number of divisions to produce, so n = 1 gives 2 new cells */
135: PetscReal r; /* The factor increase for cell height */
136: PetscScalar *h; /* The computed cell heights, based on r */
137: PetscInt *Nt; /* The array of the number of target types */
138: DMPolytopeType **target; /* The array of target types */
139: PetscInt **size; /* The array of the number of each target type */
140: PetscInt **cone; /* The array of cones for each target cell */
141: PetscInt **ornt; /* The array of orientation for each target cell */
142: } DMPlexRefine_BL;
144: PetscErrorCode DMPlexTransformSetDimensions_Internal(DMPlexTransform, DM, DM);
145: PetscErrorCode DMPlexTransformMapCoordinatesBarycenter_Internal(DMPlexTransform, DMPolytopeType, DMPolytopeType, PetscInt, PetscInt, PetscInt, PetscInt, const PetscScalar[], PetscScalar[]);
146: PetscErrorCode DMPlexTransformGetSubcellOrientation_Regular(DMPlexTransform, DMPolytopeType, PetscInt, PetscInt, DMPolytopeType, PetscInt, PetscInt, PetscInt *, PetscInt *);
147: PetscErrorCode DMPlexTransformCellRefine_Regular(DMPlexTransform, DMPolytopeType, PetscInt, PetscInt *, PetscInt *, DMPolytopeType *[], PetscInt *[], PetscInt *[], PetscInt *[]);