1 #include <petsc/private/dmpleximpl.h> /*I "petscdmplex.h" I*/ 2 3 PetscErrorCode DMPlexInvertCell_Internal(PetscInt dim, PetscInt numCorners, PetscInt cone[]) 4 { 5 int tmpc; 6 7 PetscFunctionBegin; 8 if (dim != 3) PetscFunctionReturn(0); 9 switch (numCorners) { 10 case 4: 11 tmpc = cone[0]; 12 cone[0] = cone[1]; 13 cone[1] = tmpc; 14 break; 15 case 8: 16 tmpc = cone[1]; 17 cone[1] = cone[3]; 18 cone[3] = tmpc; 19 break; 20 default: break; 21 } 22 PetscFunctionReturn(0); 23 } 24 25 /*@C 26 DMPlexInvertCell - This flips tetrahedron and hexahedron orientation since Plex stores them internally with outward normals. Other cells are left untouched. 27 28 Input Parameters: 29 + numCorners - The number of vertices in a cell 30 - cone - The incoming cone 31 32 Output Parameter: 33 . cone - The inverted cone (in-place) 34 35 Level: developer 36 37 .seealso: DMPlexGenerate() 38 @*/ 39 PetscErrorCode DMPlexInvertCell(PetscInt dim, PetscInt numCorners, int cone[]) 40 { 41 int tmpc; 42 43 PetscFunctionBegin; 44 if (dim != 3) PetscFunctionReturn(0); 45 switch (numCorners) { 46 case 4: 47 tmpc = cone[0]; 48 cone[0] = cone[1]; 49 cone[1] = tmpc; 50 break; 51 case 8: 52 tmpc = cone[1]; 53 cone[1] = cone[3]; 54 cone[3] = tmpc; 55 break; 56 default: break; 57 } 58 PetscFunctionReturn(0); 59 } 60 61 62 /*@C 63 DMPlexTriangleSetOptions - Set the options used for the Triangle mesh generator 64 65 Not Collective 66 67 Inputs Parameters: 68 + dm - The DMPlex object 69 - opts - The command line options 70 71 Level: developer 72 73 .seealso: DMPlexTetgenSetOptions(), DMPlexGenerate() 74 @*/ 75 PetscErrorCode DMPlexTriangleSetOptions(DM dm, const char *opts) 76 { 77 DM_Plex *mesh = (DM_Plex*) dm->data; 78 PetscErrorCode ierr; 79 80 PetscFunctionBegin; 81 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 82 PetscValidPointer(opts, 2); 83 ierr = PetscFree(mesh->triangleOpts);CHKERRQ(ierr); 84 ierr = PetscStrallocpy(opts, &mesh->triangleOpts);CHKERRQ(ierr); 85 PetscFunctionReturn(0); 86 } 87 88 /*@C 89 DMPlexTetgenSetOptions - Set the options used for the Tetgen mesh generator 90 91 Not Collective 92 93 Inputs Parameters: 94 + dm - The DMPlex object 95 - opts - The command line options 96 97 Level: developer 98 99 .seealso: DMPlexTriangleSetOptions(), DMPlexGenerate() 100 @*/ 101 PetscErrorCode DMPlexTetgenSetOptions(DM dm, const char *opts) 102 { 103 DM_Plex *mesh = (DM_Plex*) dm->data; 104 PetscErrorCode ierr; 105 106 PetscFunctionBegin; 107 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 108 PetscValidPointer(opts, 2); 109 ierr = PetscFree(mesh->tetgenOpts);CHKERRQ(ierr); 110 ierr = PetscStrallocpy(opts, &mesh->tetgenOpts);CHKERRQ(ierr); 111 PetscFunctionReturn(0); 112 } 113 114 /* 115 Contains the list of registered DMPlexGenerators routines 116 */ 117 PetscFunctionList DMPlexGenerateList = NULL; 118 119 struct _n_PetscFunctionList { 120 PetscErrorCode (*generate)(DM, PetscBool, DM*); 121 PetscErrorCode (*refine)(DM,double*, DM*); 122 char *name; /* string to identify routine */ 123 PetscInt dim; 124 PetscFunctionList next; /* next pointer */ 125 }; 126 127 /*@C 128 DMPlexGenerate - Generates a mesh. 129 130 Not Collective 131 132 Input Parameters: 133 + boundary - The DMPlex boundary object 134 . name - The mesh generation package name 135 - interpolate - Flag to create intermediate mesh elements 136 137 Output Parameter: 138 . mesh - The DMPlex object 139 140 Level: intermediate 141 142 .seealso: DMPlexCreate(), DMRefine() 143 @*/ 144 PetscErrorCode DMPlexGenerate(DM boundary, const char name[], PetscBool interpolate, DM *mesh) 145 { 146 PetscInt dim; 147 char genname[1024]; 148 PetscBool flg; 149 PetscErrorCode ierr; 150 PetscFunctionList fl; 151 152 PetscFunctionBegin; 153 PetscValidHeaderSpecific(boundary, DM_CLASSID, 1); 154 PetscValidLogicalCollectiveBool(boundary, interpolate, 2); 155 ierr = DMGetDimension(boundary, &dim);CHKERRQ(ierr); 156 ierr = PetscOptionsGetString(((PetscObject) boundary)->options,((PetscObject) boundary)->prefix, "-dm_plex_generator", genname, 1024, &flg);CHKERRQ(ierr); 157 if (flg) name = genname; 158 159 fl = DMPlexGenerateList; 160 if (name) { 161 while (fl) { 162 ierr = PetscStrcmp(fl->name,name,&flg);CHKERRQ(ierr); 163 if (flg) { 164 ierr = (*fl->generate)(boundary,interpolate,mesh);CHKERRQ(ierr); 165 PetscFunctionReturn(0); 166 } 167 fl = fl->next; 168 } 169 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Grid generator %g not registered",name); 170 } else { 171 while (fl) { 172 if (boundary->dim == fl->dim) { 173 ierr = (*fl->generate)(boundary,interpolate,mesh);CHKERRQ(ierr); 174 PetscFunctionReturn(0); 175 } 176 fl = fl->next; 177 } 178 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"No grid generator of dimension %D registered",boundary->dim+1); 179 } 180 PetscFunctionReturn(0); 181 } 182 183 /*@C 184 DMPlexGenerateRegister - Adds a grid generator to DMPlex 185 186 Not Collective 187 188 Input Parameters: 189 + name_solver - name of a new user-defined grid generator 190 . fnc - generator function 191 . rfnc - refinement function 192 - dim - dimension of boundary of domain 193 194 Notes: 195 DMPlexGenerateRegister() may be called multiple times to add several user-defined solvers. 196 197 Sample usage: 198 .vb 199 DMPlexGenerateRegister("my_generator",MyGeneratorCreate,MyGeneratorRefiner,dim); 200 .ve 201 202 Then, your generator can be chosen with the procedural interface via 203 $ DMPlexGenerate(dm,"my_generator",...) 204 or at runtime via the option 205 $ -dm_plex_generator my_generator 206 207 Level: advanced 208 209 .seealso: DMPlexGenerateRegisterAll(), DMPlexGenerate(), DMPlexGenerateRegisterDestroy() 210 211 @*/ 212 PetscErrorCode DMPlexGenerateRegister(const char sname[],PetscErrorCode (*fnc)(DM, PetscBool,DM*), PetscErrorCode (*rfnc)(DM, double*,DM*),PetscInt dim) 213 { 214 PetscErrorCode ierr; 215 PetscFunctionList entry; 216 217 PetscFunctionBegin; 218 ierr = PetscNew(&entry);CHKERRQ(ierr); 219 ierr = PetscStrallocpy(sname,&entry->name);CHKERRQ(ierr); 220 entry->generate = fnc; 221 entry->refine = rfnc; 222 entry->dim = dim; 223 entry->next = NULL; 224 if (!DMPlexGenerateList) DMPlexGenerateList = entry; 225 else { 226 PetscFunctionList fl = DMPlexGenerateList; 227 while (fl->next) fl = fl->next; 228 fl->next = entry; 229 } 230 PetscFunctionReturn(0); 231 } 232 233 extern PetscBool DMPlexGenerateRegisterAllCalled; 234 235 PetscErrorCode DMPlexGenerateRegisterDestroy(void) 236 { 237 PetscFunctionList next,fl; 238 PetscErrorCode ierr; 239 240 PetscFunctionBegin; 241 next = fl = DMPlexGenerateList; 242 while (next) { 243 next = fl ? fl->next : NULL; 244 if (fl) {ierr = PetscFree(fl->name);CHKERRQ(ierr);} 245 ierr = PetscFree(fl);CHKERRQ(ierr); 246 fl = next; 247 } 248 DMPlexGenerateList = NULL; 249 DMPlexGenerateRegisterAllCalled = PETSC_FALSE; 250 PetscFunctionReturn(0); 251 } 252