1 #include <petsc/private/dmimpl.h> /*I "petscdm.h" I*/
2
3 PetscFunctionList DMGeomModelList = NULL;
4 PetscBool DMGeomModelRegisterAllCalled = PETSC_FALSE;
5
6 #if defined(PETSC_HAVE_EGADS)
7 PETSC_INTERN PetscErrorCode DMSnapToGeomModel_EGADS(DM, PetscInt, PetscInt, const PetscScalar[], PetscScalar[]);
8 #endif
9
DMSnapToGeomModelBall(DM dm,PetscInt p,PetscInt dE,const PetscScalar mcoords[],PetscScalar gcoords[])10 static PetscErrorCode DMSnapToGeomModelBall(DM dm, PetscInt p, PetscInt dE, const PetscScalar mcoords[], PetscScalar gcoords[])
11 {
12 PetscInt val;
13
14 PetscFunctionBeginUser;
15 PetscCall(DMGetLabelValue(dm, "marker", p, &val));
16 if (val >= 0) {
17 PetscReal norm = 0.;
18
19 for (PetscInt d = 0; d < dE; ++d) norm += PetscSqr(PetscRealPart(mcoords[d]));
20 norm = PetscSqrtReal(norm);
21 for (PetscInt d = 0; d < dE; ++d) gcoords[d] = mcoords[d] / norm;
22 } else {
23 for (PetscInt d = 0; d < dE; ++d) gcoords[d] = mcoords[d];
24 }
25 PetscFunctionReturn(PETSC_SUCCESS);
26 }
27
DMSnapToGeomModelCylinder(DM dm,PetscInt p,PetscInt dE,const PetscScalar mcoords[],PetscScalar gcoords[])28 static PetscErrorCode DMSnapToGeomModelCylinder(DM dm, PetscInt p, PetscInt dE, const PetscScalar mcoords[], PetscScalar gcoords[])
29 {
30 PetscReal gmin[3], gmax[3];
31 PetscInt val;
32
33 PetscFunctionBeginUser;
34 PetscCall(DMGetBoundingBox(dm, gmin, gmax));
35 PetscCall(DMGetLabelValue(dm, "generatrix", p, &val));
36 if (val >= 0) {
37 PetscReal norm = 0.;
38
39 for (PetscInt d = 0; d < dE - 1; ++d) norm += PetscSqr(PetscRealPart(mcoords[d]));
40 norm = PetscSqrtReal(norm);
41 for (PetscInt d = 0; d < dE - 1; ++d) gcoords[d] = mcoords[d] * gmax[0] / norm;
42 gcoords[dE - 1] = mcoords[dE - 1];
43 } else {
44 for (PetscInt d = 0; d < dE; ++d) gcoords[d] = mcoords[d];
45 }
46 PetscFunctionReturn(PETSC_SUCCESS);
47 }
48
49 /*@C
50 DMGeomModelRegisterAll - Registers all of the geometry model methods in the `DM` package.
51
52 Not Collective
53
54 Level: advanced
55
56 .seealso: `DM`, `DMGeomModelRegisterDestroy()`
57 @*/
DMGeomModelRegisterAll(void)58 PetscErrorCode DMGeomModelRegisterAll(void)
59 {
60 PetscFunctionBegin;
61 if (DMGeomModelRegisterAllCalled) PetscFunctionReturn(PETSC_SUCCESS);
62 DMGeomModelRegisterAllCalled = PETSC_TRUE;
63 PetscCall(DMGeomModelRegister("ball", DMSnapToGeomModelBall));
64 PetscCall(DMGeomModelRegister("cylinder", DMSnapToGeomModelCylinder));
65 #if defined(PETSC_HAVE_EGADS)
66 // FIXME: Brandon uses DMPlexSnapToGeomModel() here instead
67 PetscCall(DMGeomModelRegister("egads", DMSnapToGeomModel_EGADS));
68 #endif
69 PetscFunctionReturn(PETSC_SUCCESS);
70 }
71
72 /*@C
73 DMGeomModelRegister - Adds a geometry model to `DM`
74
75 Not Collective, No Fortran Support
76
77 Input Parameters:
78 + sname - name of a new user-defined geometry model
79 - fnc - geometry model function
80
81 Example Usage:
82 .vb
83 DMGeomModelRegister("my_geom_model", MySnapToGeomModel);
84 .ve
85
86 Then, your generator can be chosen with the procedural interface via
87 .vb
88 DMSetGeomModel(dm, "my_geom_model",...)
89 .ve
90 or at runtime via the option
91 .vb
92 -dm_geom_model my_geom_model
93 .ve
94
95 Level: advanced
96
97 Note:
98 `DMGeomModelRegister()` may be called multiple times to add several user-defined generators
99
100 .seealso: `DM`, `DMGeomModelRegisterAll()`, `DMPlexGeomModel()`, `DMGeomModelRegisterDestroy()`
101 @*/
DMGeomModelRegister(const char sname[],PetscErrorCode (* fnc)(DM,PetscInt,PetscInt,const PetscScalar[],PetscScalar[]))102 PetscErrorCode DMGeomModelRegister(const char sname[], PetscErrorCode (*fnc)(DM, PetscInt, PetscInt, const PetscScalar[], PetscScalar[]))
103 {
104 PetscFunctionBegin;
105 PetscCall(PetscFunctionListAdd(&DMGeomModelList, sname, fnc));
106 PetscFunctionReturn(PETSC_SUCCESS);
107 }
108
DMGeomModelRegisterDestroy(void)109 PetscErrorCode DMGeomModelRegisterDestroy(void)
110 {
111 PetscFunctionBegin;
112 PetscCall(PetscFunctionListDestroy(&DMGeomModelList));
113 DMGeomModelRegisterAllCalled = PETSC_FALSE;
114 PetscFunctionReturn(PETSC_SUCCESS);
115 }
116
117 /*@
118 DMSetSnapToGeomModel - Choose a geometry model for this `DM`.
119
120 Not Collective
121
122 Input Parameters:
123 + dm - The `DM` object
124 - name - A geometry model name, or `NULL` for the default
125
126 Level: intermediate
127
128 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMRefine()`, `DMPlexCreate()`, `DMSnapToGeomModel()`
129 @*/
DMSetSnapToGeomModel(DM dm,const char name[])130 PetscErrorCode DMSetSnapToGeomModel(DM dm, const char name[])
131 {
132 char geomname[PETSC_MAX_PATH_LEN];
133 PetscBool flg;
134
135 PetscFunctionBegin;
136 if (!name && dm->ops->snaptogeommodel) PetscFunctionReturn(PETSC_SUCCESS);
137 PetscCall(PetscOptionsGetString(((PetscObject)dm)->options, ((PetscObject)dm)->prefix, "-dm_geom_model", geomname, sizeof(geomname), &flg));
138 if (flg) name = geomname;
139 if (!name) {
140 PetscObject modelObj;
141
142 PetscCall(PetscObjectQuery((PetscObject)dm, "EGADS Model", &modelObj));
143 if (modelObj) name = "egads";
144 else {
145 PetscCall(PetscObjectQuery((PetscObject)dm, "EGADSlite Model", &modelObj));
146 if (modelObj) name = "egads";
147 }
148 }
149 if (!name) PetscFunctionReturn(PETSC_SUCCESS);
150
151 PetscCall(PetscFunctionListFind(DMGeomModelList, name, &dm->ops->snaptogeommodel));
152 PetscCheck(dm->ops->snaptogeommodel, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "Geometry model %s not registered; you may need to add --download-%s to your ./configure options", name, name);
153 PetscFunctionReturn(PETSC_SUCCESS);
154 }
155
156 /*@
157 DMSnapToGeomModel - Given a coordinate point 'mcoords' on the mesh point 'p', return the closest coordinate point 'gcoords' on the geometry model associated with that point.
158
159 Not Collective
160
161 Input Parameters:
162 + dm - The `DMPLEX` object
163 . p - The mesh point
164 . dE - The coordinate dimension
165 - mcoords - A coordinate point lying on the mesh point
166
167 Output Parameter:
168 . gcoords - The closest coordinate point on the geometry model associated with 'p' to the given point
169
170 Level: intermediate
171
172 Note:
173 Returns the original coordinates if no geometry model is found.
174
175 The coordinate dimension may be different from the coordinate dimension of the `dm`, for example if the transformation is extrusion.
176
177 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMRefine()`, `DMPlexCreate()`, `DMPlexSetRefinementUniform()`
178 @*/
DMSnapToGeomModel(DM dm,PetscInt p,PetscInt dE,const PetscScalar mcoords[],PetscScalar gcoords[])179 PetscErrorCode DMSnapToGeomModel(DM dm, PetscInt p, PetscInt dE, const PetscScalar mcoords[], PetscScalar gcoords[])
180 {
181 PetscFunctionBegin;
182 if (!dm->ops->snaptogeommodel)
183 for (PetscInt d = 0; d < dE; ++d) gcoords[d] = mcoords[d];
184 else PetscUseTypeMethod(dm, snaptogeommodel, p, dE, mcoords, gcoords);
185 PetscFunctionReturn(PETSC_SUCCESS);
186 }
187