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