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