Actual source code: plexgenerate.c

petsc-3.9.4 2018-09-11
Report Typos and Errors
  1:  #include <petsc/private/dmpleximpl.h>

  3: PetscErrorCode DMPlexInvertCell_Internal(PetscInt dim, PetscInt numCorners, PetscInt cone[])
  4: {
  5:   int tmpc;

  8:   if (dim != 3) return(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:   return(0);
 23: }

 25: /*@C
 26:   DMPlexInvertCell - This flips tetrahedron and hexahedron orientation since Plex stores them internally with outward normals. Other cells are left untouched.

 28:   Input Parameters:
 29: + numCorners - The number of vertices in a cell
 30: - cone - The incoming cone

 32:   Output Parameter:
 33: . cone - The inverted cone (in-place)

 35:   Level: developer

 37: .seealso: DMPlexGenerate()
 38: @*/
 39: PetscErrorCode DMPlexInvertCell(PetscInt dim, PetscInt numCorners, int cone[])
 40: {
 41:   int tmpc;

 44:   if (dim != 3) return(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:   return(0);
 59: }


 62: /*@C
 63:   DMPlexTriangleSetOptions - Set the options used for the Triangle mesh generator

 65:   Not Collective

 67:   Inputs Parameters:
 68: + dm - The DMPlex object
 69: - opts - The command line options

 71:   Level: developer

 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;

 84:   PetscFree(mesh->triangleOpts);
 85:   PetscStrallocpy(opts, &mesh->triangleOpts);
 86:   return(0);
 87: }

 89: /*@C
 90:   DMPlexTetgenSetOptions - Set the options used for the Tetgen mesh generator

 92:   Not Collective

 94:   Inputs Parameters:
 95: + dm - The DMPlex object
 96: - opts - The command line options

 98:   Level: developer

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;

111:   PetscFree(mesh->tetgenOpts);
112:   PetscStrallocpy(opts, &mesh->tetgenOpts);
113:   return(0);
114: }

116: /*
117:    Contains the list of registered DMPlexGenerators routines
118: */
119: PetscFunctionList DMPlexGenerateList = NULL;

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: };

129: /*@C
130:   DMPlexGenerate - Generates a mesh.

132:   Not Collective

134:   Input Parameters:
135: + boundary - The DMPlex boundary object
136: . name - The mesh generation package name
137: - interpolate - Flag to create intermediate mesh elements

139:   Output Parameter:
140: . mesh - The DMPlex object

142:   Level: intermediate

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;

158:   DMGetDimension(boundary, &dim);
159:   PetscOptionsGetString(((PetscObject) boundary)->options,((PetscObject) boundary)->prefix, "-dm_plex_generator", genname, 1024, &flg);
160:   if (flg) name = genname;

162:   fl = DMPlexGenerateList;
163:   if (name) {
164:     while (fl) {
165:       PetscStrcmp(fl->name,name,&flg);
166:       if (flg) {
167:         (*fl->generate)(boundary,interpolate,mesh);
168:         return(0);
169:       }
170:       fl = fl->next;
171:     }
172:     SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Grid generator %g not registered",name);
173:   } else {
174:     while (fl) {
175:       if (boundary->dim == fl->dim) {
176:         (*fl->generate)(boundary,interpolate,mesh);
177:         return(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);
182:   }
183:   return(0);
184: }

186: /*@C
187:   DMPlexGenerateRegister -  Adds a grid generator to DMPlex

189:    Not Collective

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

197:    Notes:
198:    DMPlexGenerateRegister() may be called multiple times to add several user-defined solvers.

200:    Sample usage:
201: .vb
202:    DMPlexGenerateRegister("my_generator",MyGeneratorCreate);
203: .ve

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

210:    Level: advanced

212: .keywords: DMPlexGenerate, register

214: .seealso: DMPlexGenerateRegisterAll(), DMPlexGenerateRegisterDestroy()

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;

223:   PetscNew(&entry);
224:   PetscStrallocpy(sname,&entry->name);
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:   return(0);
236: }

238: extern PetscBool DMPlexGenerateRegisterAllCalled;

240: PetscErrorCode  DMPlexGenerateRegisterDestroy(void)
241: {
242:   PetscFunctionList next,fl;
243:   PetscErrorCode    ierr;

246:   next = fl =  DMPlexGenerateList;
247:     while (next) {
248:     next = fl ? fl->next : NULL;
249:     if (fl) {PetscFree(fl->name);}
250:     PetscFree(fl);
251:     fl = next;
252:   }
253:   DMPlexGenerateList              = NULL;
254:   DMPlexGenerateRegisterAllCalled = PETSC_FALSE;
255:   return(0);
256: }