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 .keywords: mesh, points 74 .seealso: DMPlexTetgenSetOptions(), DMPlexGenerate() 75 @*/ 76 PetscErrorCode DMPlexTriangleSetOptions(DM dm, const char *opts) 77 { 78 DM_Plex *mesh = (DM_Plex*) dm->data; 79 PetscErrorCode ierr; 80 81 PetscFunctionBegin; 82 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 83 PetscValidPointer(opts, 2); 84 ierr = PetscFree(mesh->triangleOpts);CHKERRQ(ierr); 85 ierr = PetscStrallocpy(opts, &mesh->triangleOpts);CHKERRQ(ierr); 86 PetscFunctionReturn(0); 87 } 88 89 /*@C 90 DMPlexTetgenSetOptions - Set the options used for the Tetgen mesh generator 91 92 Not Collective 93 94 Inputs Parameters: 95 + dm - The DMPlex object 96 - opts - The command line options 97 98 Level: developer 99 100 .keywords: mesh, points 101 .seealso: DMPlexTriangleSetOptions(), DMPlexGenerate() 102 @*/ 103 PetscErrorCode DMPlexTetgenSetOptions(DM dm, const char *opts) 104 { 105 DM_Plex *mesh = (DM_Plex*) dm->data; 106 PetscErrorCode ierr; 107 108 PetscFunctionBegin; 109 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 110 PetscValidPointer(opts, 2); 111 ierr = PetscFree(mesh->tetgenOpts);CHKERRQ(ierr); 112 ierr = PetscStrallocpy(opts, &mesh->tetgenOpts);CHKERRQ(ierr); 113 PetscFunctionReturn(0); 114 } 115 116 /* 117 Contains the list of registered DMPlexGenerators routines 118 */ 119 PetscFunctionList DMPlexGenerateList = NULL; 120 121 struct _n_PetscFunctionList { 122 PetscErrorCode (*generate)(DM, PetscBool, DM*); 123 PetscErrorCode (*refine)(DM,double*, DM*); 124 char *name; /* string to identify routine */ 125 PetscInt dim; 126 PetscFunctionList next; /* next pointer */ 127 }; 128 129 /*@C 130 DMPlexGenerate - Generates a mesh. 131 132 Not Collective 133 134 Input Parameters: 135 + boundary - The DMPlex boundary object 136 . name - The mesh generation package name 137 - interpolate - Flag to create intermediate mesh elements 138 139 Output Parameter: 140 . mesh - The DMPlex object 141 142 Level: intermediate 143 144 .keywords: mesh, elements 145 .seealso: DMPlexCreate(), DMRefine() 146 @*/ 147 PetscErrorCode DMPlexGenerate(DM boundary, const char name[], PetscBool interpolate, DM *mesh) 148 { 149 PetscInt dim; 150 char genname[1024]; 151 PetscBool flg; 152 PetscErrorCode ierr; 153 PetscFunctionList fl; 154 155 PetscFunctionBegin; 156 PetscValidHeaderSpecific(boundary, DM_CLASSID, 1); 157 PetscValidLogicalCollectiveBool(boundary, interpolate, 2); 158 ierr = DMGetDimension(boundary, &dim);CHKERRQ(ierr); 159 ierr = PetscOptionsGetString(((PetscObject) boundary)->options,((PetscObject) boundary)->prefix, "-dm_plex_generator", genname, 1024, &flg);CHKERRQ(ierr); 160 if (flg) name = genname; 161 162 fl = DMPlexGenerateList; 163 if (name) { 164 while (fl) { 165 ierr = PetscStrcmp(fl->name,name,&flg);CHKERRQ(ierr); 166 if (flg) { 167 (*fl->generate)(boundary,interpolate,mesh);CHKERRQ(ierr); 168 PetscFunctionReturn(0); 169 } 170 fl = fl->next; 171 } 172 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Grid generator %g not registered",name);CHKERRQ(ierr); 173 } else { 174 while (fl) { 175 if (boundary->dim == fl->dim) { 176 (*fl->generate)(boundary,interpolate,mesh);CHKERRQ(ierr); 177 PetscFunctionReturn(0); 178 } 179 fl = fl->next; 180 } 181 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"No grid generator of dimension %D registered",boundary->dim);CHKERRQ(ierr); 182 } 183 PetscFunctionReturn(0); 184 } 185 186 /*@C 187 DMPlexGenerateRegister - Adds a grid generator to DMPlex 188 189 Not Collective 190 191 Input Parameters: 192 + name_solver - name of a new user-defined solver 193 . fnc - generator function 194 . rfnc - refinement function 195 - dim - dimension of boundary of domain 196 197 Notes: 198 DMPlexGenerateRegister() may be called multiple times to add several user-defined solvers. 199 200 Sample usage: 201 .vb 202 DMPlexGenerateRegister("my_generator",MyGeneratorCreate); 203 .ve 204 205 Then, your generator can be chosen with the procedural interface via 206 $ DMPlexGenerateSetType(ksp,"my_generator") 207 or at runtime via the option 208 $ -dmplex_generate_type my_generator 209 210 Level: advanced 211 212 .keywords: DMPlexGenerate, register 213 214 .seealso: DMPlexGenerateRegisterAll(), DMPlexGenerateRegisterDestroy() 215 216 @*/ 217 PetscErrorCode DMPlexGenerateRegister(const char sname[],PetscErrorCode (*fnc)(DM, PetscBool,DM*), PetscErrorCode (*rfnc)(DM, double*,DM*),PetscInt dim) 218 { 219 PetscErrorCode ierr; 220 PetscFunctionList entry; 221 222 PetscFunctionBegin; 223 ierr = PetscNew(&entry);CHKERRQ(ierr); 224 ierr = PetscStrallocpy(sname,&entry->name);CHKERRQ(ierr); 225 entry->generate = fnc; 226 entry->refine = rfnc; 227 entry->dim = dim; 228 entry->next = NULL; 229 if (!DMPlexGenerateList) DMPlexGenerateList = entry; 230 else { 231 PetscFunctionList fl = DMPlexGenerateList; 232 while (fl->next) fl = fl->next; 233 fl->next = entry; 234 } 235 PetscFunctionReturn(0); 236 } 237 238 extern PetscBool DMPlexGenerateRegisterAllCalled; 239 240 PetscErrorCode DMPlexGenerateRegisterDestroy(void) 241 { 242 PetscFunctionList next,fl; 243 PetscErrorCode ierr; 244 245 PetscFunctionBegin; 246 next = fl = DMPlexGenerateList; 247 while (next) { 248 next = fl ? fl->next : NULL; 249 if (fl) {ierr = PetscFree(fl->name);CHKERRQ(ierr);} 250 ierr = PetscFree(fl);CHKERRQ(ierr); 251 fl = next; 252 } 253 DMPlexGenerateList = NULL; 254 DMPlexGenerateRegisterAllCalled = PETSC_FALSE; 255 PetscFunctionReturn(0); 256 } 257