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