1 #include <petsc/private/dmimpl.h> /*I "petscdm.h" I*/ 2 3 PETSC_EXTERN PetscErrorCode DMIsForest(DM, PetscBool *); 4 5 DMGeneratorFunctionList DMGenerateList = NULL; 6 PetscBool DMGenerateRegisterAllCalled = PETSC_FALSE; 7 8 #if defined(PETSC_HAVE_TRIANGLE) 9 PETSC_EXTERN PetscErrorCode DMPlexGenerate_Triangle(DM, PetscBool, DM *); 10 PETSC_EXTERN PetscErrorCode DMPlexRefine_Triangle(DM, PetscReal *, DM *); 11 #endif 12 #if defined(PETSC_HAVE_TETGEN) 13 PETSC_EXTERN PetscErrorCode DMPlexGenerate_Tetgen(DM, PetscBool, DM *); 14 PETSC_EXTERN PetscErrorCode DMPlexRefine_Tetgen(DM, double *, DM *); 15 #endif 16 #if defined(PETSC_HAVE_CTETGEN) 17 PETSC_EXTERN PetscErrorCode DMPlexGenerate_CTetgen(DM, PetscBool, DM *); 18 PETSC_EXTERN PetscErrorCode DMPlexRefine_CTetgen(DM, double *, DM *); 19 #endif 20 #if defined(PETSC_HAVE_PRAGMATIC) 21 PETSC_EXTERN PetscErrorCode DMAdaptMetric_Pragmatic_Plex(DM, Vec, DMLabel, DMLabel, DM *); 22 #endif 23 #if defined(PETSC_HAVE_MMG) 24 PETSC_EXTERN PetscErrorCode DMAdaptMetric_Mmg_Plex(DM, Vec, DMLabel, DMLabel, DM *); 25 #endif 26 #if defined(PETSC_HAVE_PARMMG) 27 PETSC_EXTERN PetscErrorCode DMAdaptMetric_ParMmg_Plex(DM, Vec, DMLabel, DMLabel, DM *); 28 #endif 29 PETSC_EXTERN PetscErrorCode DMPlexTransformAdaptLabel(DM, Vec, DMLabel, DMLabel, DM *); 30 PETSC_EXTERN PetscErrorCode DMAdaptLabel_Forest(DM, Vec, DMLabel, DMLabel, DM *); 31 32 /*@C 33 DMGenerateRegisterAll - Registers all of the mesh generation methods in the `DM` package. 34 35 Not Collective 36 37 Level: advanced 38 39 .seealso: `DM`, `DMGenerateRegisterDestroy()` 40 @*/ 41 PetscErrorCode DMGenerateRegisterAll(void) 42 { 43 PetscFunctionBegin; 44 if (DMGenerateRegisterAllCalled) PetscFunctionReturn(PETSC_SUCCESS); 45 DMGenerateRegisterAllCalled = PETSC_TRUE; 46 #if defined(PETSC_HAVE_TRIANGLE) 47 PetscCall(DMGenerateRegister("triangle", DMPlexGenerate_Triangle, DMPlexRefine_Triangle, NULL, 1)); 48 #endif 49 #if defined(PETSC_HAVE_CTETGEN) 50 PetscCall(DMGenerateRegister("ctetgen", DMPlexGenerate_CTetgen, DMPlexRefine_CTetgen, NULL, 2)); 51 #endif 52 #if defined(PETSC_HAVE_TETGEN) 53 PetscCall(DMGenerateRegister("tetgen", DMPlexGenerate_Tetgen, DMPlexRefine_Tetgen, NULL, 2)); 54 #endif 55 #if defined(PETSC_HAVE_PRAGMATIC) 56 PetscCall(DMGenerateRegister("pragmatic", NULL, NULL, DMAdaptMetric_Pragmatic_Plex, -1)); 57 #endif 58 #if defined(PETSC_HAVE_MMG) 59 PetscCall(DMGenerateRegister("mmg", NULL, NULL, DMAdaptMetric_Mmg_Plex, -1)); 60 #endif 61 #if defined(PETSC_HAVE_PARMMG) 62 PetscCall(DMGenerateRegister("parmmg", NULL, NULL, DMAdaptMetric_ParMmg_Plex, -1)); 63 #endif 64 PetscCall(DMGenerateRegister("cellrefiner", NULL, NULL, DMPlexTransformAdaptLabel, -1)); 65 PetscCall(DMGenerateRegister("forest", NULL, NULL, DMAdaptLabel_Forest, -1)); 66 PetscFunctionReturn(PETSC_SUCCESS); 67 } 68 69 /*@C 70 DMGenerateRegister - Adds a grid generator to `DM` 71 72 Not Collective, No Fortran Support 73 74 Input Parameters: 75 + sname - name of a new user-defined grid generator 76 . fnc - generator function 77 . rfnc - refinement function 78 . alfnc - adapt by label function 79 - dim - dimension of boundary of domain 80 81 Example Usage: 82 .vb 83 DMGenerateRegister("my_generator", MyGeneratorCreate, MyGeneratorRefiner, MyGeneratorAdaptor, dim); 84 .ve 85 86 Then, your generator can be chosen with the procedural interface via 87 .vb 88 DMGenerate(dm, "my_generator",...) 89 .ve 90 or at runtime via the option 91 .vb 92 -dm_generator my_generator 93 .ve 94 95 Level: advanced 96 97 Note: 98 `DMGenerateRegister()` may be called multiple times to add several user-defined generators 99 100 .seealso: `DM`, `DMGenerateRegisterAll()`, `DMPlexGenerate()`, `DMGenerateRegisterDestroy()` 101 @*/ 102 PetscErrorCode DMGenerateRegister(const char sname[], PetscErrorCode (*fnc)(DM, PetscBool, DM *), PetscErrorCode (*rfnc)(DM, PetscReal *, DM *), PetscErrorCode (*alfnc)(DM, Vec, DMLabel, DMLabel, DM *), PetscInt dim) 103 { 104 DMGeneratorFunctionList entry; 105 106 PetscFunctionBegin; 107 PetscCall(PetscNew(&entry)); 108 PetscCall(PetscStrallocpy(sname, &entry->name)); 109 entry->generate = fnc; 110 entry->refine = rfnc; 111 entry->adapt = alfnc; 112 entry->dim = dim; 113 entry->next = NULL; 114 if (!DMGenerateList) DMGenerateList = entry; 115 else { 116 DMGeneratorFunctionList fl = DMGenerateList; 117 while (fl->next) fl = fl->next; 118 fl->next = entry; 119 } 120 PetscFunctionReturn(PETSC_SUCCESS); 121 } 122 123 extern PetscBool DMGenerateRegisterAllCalled; 124 125 PetscErrorCode DMGenerateRegisterDestroy(void) 126 { 127 DMGeneratorFunctionList next, fl; 128 129 PetscFunctionBegin; 130 next = fl = DMGenerateList; 131 while (next) { 132 next = fl ? fl->next : NULL; 133 if (fl) PetscCall(PetscFree(fl->name)); 134 PetscCall(PetscFree(fl)); 135 fl = next; 136 } 137 DMGenerateList = NULL; 138 DMGenerateRegisterAllCalled = PETSC_FALSE; 139 PetscFunctionReturn(PETSC_SUCCESS); 140 } 141 142 /*@ 143 DMAdaptLabel - Adapt a `DM` based on a `DMLabel` with values interpreted as coarsening and refining flags. Specific implementations of `DM` maybe have 144 specialized flags, but all implementations should accept flag values `DM_ADAPT_DETERMINE`, `DM_ADAPT_KEEP`, `DM_ADAPT_REFINE`, and, 145 `DM_ADAPT_COARSEN`. 146 147 Collective 148 149 Input Parameters: 150 + dm - the pre-adaptation `DM` object 151 - label - label with the flags 152 153 Output Parameters: 154 . dmAdapt - the adapted `DM` object: may be `NULL` if an adapted `DM` could not be produced. 155 156 Level: intermediate 157 158 .seealso: `DM`, `DMAdaptMetric()`, `DMCoarsen()`, `DMRefine()` 159 @*/ 160 PetscErrorCode DMAdaptLabel(DM dm, DMLabel label, DM *dmAdapt) 161 { 162 DMGeneratorFunctionList fl; 163 char adaptname[PETSC_MAX_PATH_LEN]; 164 const char *name; 165 PetscInt dim; 166 PetscBool flg, isForest, found = PETSC_FALSE; 167 168 PetscFunctionBegin; 169 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 170 if (label) PetscValidHeaderSpecific(label, DMLABEL_CLASSID, 2); 171 PetscAssertPointer(dmAdapt, 3); 172 *dmAdapt = NULL; 173 PetscCall(DMGetDimension(dm, &dim)); 174 PetscCall(DMIsForest(dm, &isForest)); 175 name = isForest ? "forest" : "cellrefiner"; 176 PetscCall(PetscOptionsGetString(((PetscObject)dm)->options, ((PetscObject)dm)->prefix, "-dm_adaptor", adaptname, sizeof(adaptname), &flg)); 177 if (flg) name = adaptname; 178 179 fl = DMGenerateList; 180 while (fl) { 181 PetscCall(PetscStrcmp(fl->name, name, &flg)); 182 if (flg) { 183 PetscCall((*fl->adapt)(dm, NULL, label, NULL, dmAdapt)); 184 found = PETSC_TRUE; 185 } 186 fl = fl->next; 187 } 188 PetscCheck(found, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Grid adaptor %s not registered; you may need to add --download-%s to your ./configure options", name, name); 189 if (*dmAdapt) { 190 void *ctx; 191 192 (*dmAdapt)->prealloc_only = dm->prealloc_only; /* maybe this should go .... */ 193 PetscCall(PetscFree((*dmAdapt)->vectype)); 194 PetscCall(PetscStrallocpy(dm->vectype, (char **)&(*dmAdapt)->vectype)); 195 PetscCall(PetscFree((*dmAdapt)->mattype)); 196 PetscCall(PetscStrallocpy(dm->mattype, (char **)&(*dmAdapt)->mattype)); 197 PetscCall(DMGetApplicationContext(dm, &ctx)); 198 PetscCall(DMSetApplicationContext(*dmAdapt, ctx)); 199 } 200 PetscFunctionReturn(PETSC_SUCCESS); 201 } 202 203 /*@ 204 DMAdaptMetric - Generates a mesh adapted to the specified metric field. 205 206 Input Parameters: 207 + dm - The DM object 208 . metric - The metric to which the mesh is adapted, defined vertex-wise. 209 . bdLabel - Label for boundary tags, which will be preserved in the output mesh. `bdLabel` should be `NULL` if there is no such label, and should be different from "_boundary_". 210 - rgLabel - Label for cell tags, which will be preserved in the output mesh. `rgLabel` should be `NULL` if there is no such label, and should be different from "_regions_". 211 212 Output Parameter: 213 . dmAdapt - Pointer to the `DM` object containing the adapted mesh 214 215 Note: 216 The label in the adapted mesh will be registered under the name of the input `DMLabel` object 217 218 Level: advanced 219 220 .seealso: `DMAdaptLabel()`, `DMCoarsen()`, `DMRefine()` 221 @*/ 222 PetscErrorCode DMAdaptMetric(DM dm, Vec metric, DMLabel bdLabel, DMLabel rgLabel, DM *dmAdapt) 223 { 224 DMGeneratorFunctionList fl; 225 char adaptname[PETSC_MAX_PATH_LEN]; 226 const char *name; 227 const char *const adaptors[3] = {"pragmatic", "mmg", "parmmg"}; 228 PetscInt dim; 229 PetscBool flg, found = PETSC_FALSE; 230 231 PetscFunctionBegin; 232 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 233 PetscValidHeaderSpecific(metric, VEC_CLASSID, 2); 234 if (bdLabel) PetscValidHeaderSpecific(bdLabel, DMLABEL_CLASSID, 3); 235 if (rgLabel) PetscValidHeaderSpecific(rgLabel, DMLABEL_CLASSID, 4); 236 PetscAssertPointer(dmAdapt, 5); 237 *dmAdapt = NULL; 238 PetscCall(DMGetDimension(dm, &dim)); 239 PetscCall(PetscOptionsGetString(((PetscObject)dm)->options, ((PetscObject)dm)->prefix, "-dm_adaptor", adaptname, sizeof(adaptname), &flg)); 240 241 /* Default to Mmg in serial and ParMmg in parallel */ 242 if (flg) name = adaptname; 243 else { 244 MPI_Comm comm; 245 PetscMPIInt size; 246 247 PetscCall(PetscObjectGetComm((PetscObject)dm, &comm)); 248 PetscCallMPI(MPI_Comm_size(comm, &size)); 249 if (size == 1) name = adaptors[1]; 250 else name = adaptors[2]; 251 } 252 253 fl = DMGenerateList; 254 while (fl) { 255 PetscCall(PetscStrcmp(fl->name, name, &flg)); 256 if (flg) { 257 PetscCall((*fl->adapt)(dm, metric, bdLabel, rgLabel, dmAdapt)); 258 found = PETSC_TRUE; 259 } 260 fl = fl->next; 261 } 262 PetscCheck(found, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Grid adaptor %s not registered; you may need to add --download-%s to your ./configure options", name, name); 263 if (*dmAdapt) { 264 (*dmAdapt)->prealloc_only = dm->prealloc_only; /* maybe this should go .... */ 265 PetscCall(PetscFree((*dmAdapt)->vectype)); 266 PetscCall(PetscStrallocpy(dm->vectype, (char **)&(*dmAdapt)->vectype)); 267 PetscCall(PetscFree((*dmAdapt)->mattype)); 268 PetscCall(PetscStrallocpy(dm->mattype, (char **)&(*dmAdapt)->mattype)); 269 } 270 PetscFunctionReturn(PETSC_SUCCESS); 271 } 272