Actual source code: plexgenerate.c

petsc-master 2020-02-23
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: .seealso: DMPlexTetgenSetOptions(), DMPlexGenerate()
 74: @*/
 75: PetscErrorCode DMPlexTriangleSetOptions(DM dm, const char *opts)
 76: {
 77:   DM_Plex       *mesh = (DM_Plex*) dm->data;

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

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

 91:   Not Collective

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

 97:   Level: developer

 99: .seealso: DMPlexTriangleSetOptions(), DMPlexGenerate()
100: @*/
101: PetscErrorCode DMPlexTetgenSetOptions(DM dm, const char *opts)
102: {
103:   DM_Plex       *mesh = (DM_Plex*) dm->data;

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

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

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

127: /*@C
128:   DMPlexGenerate - Generates a mesh.

130:   Not Collective

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

137:   Output Parameter:
138: . mesh - The DMPlex object

140:   Options Database:
141: +  -dm_plex_generate <name> - package to generate mesh, for example, triangle, ctetgen or tetgen
142: -  -dm_plex_generator <name> - package to generate mesh, for example, triangle, ctetgen or tetgen (deprecated)

144:   Level: intermediate

146: .seealso: DMPlexCreate(), DMRefine()
147: @*/
148: PetscErrorCode DMPlexGenerate(DM boundary, const char name[], PetscBool interpolate, DM *mesh)
149: {
150:   PetscInt          dim;
151:   char              genname[1024];
152:   PetscBool         flg;
153:   PetscErrorCode    ierr;
154:   PetscFunctionList fl;
155:   const char*       suggestions;

160:   DMGetDimension(boundary, &dim);
161:   PetscOptionsGetString(((PetscObject) boundary)->options,((PetscObject) boundary)->prefix, "-dm_plex_generator", genname, 1024, &flg);
162:   if (flg) name = genname;
163:   else {
164:     PetscOptionsGetString(((PetscObject) boundary)->options,((PetscObject) boundary)->prefix, "-dm_plex_generate", genname, 1024, &flg);
165:     if (flg) name = genname;
166:   }

168:   fl = DMPlexGenerateList;
169:   if (name) {
170:     while (fl) {
171:       PetscStrcmp(fl->name,name,&flg);
172:       if (flg) {
173:         (*fl->generate)(boundary,interpolate,mesh);
174:         return(0);
175:       }
176:       fl = fl->next;
177:     }
178:     SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Grid generator %s not registered; you may need to add --download-%s to your ./configure options",name,name);
179:   } else {
180:     while (fl) {
181:       if (boundary->dim == fl->dim) {
182:         (*fl->generate)(boundary,interpolate,mesh);
183:         return(0);
184:       }
185:       fl = fl->next;
186:     }
187:     suggestions = "";
188:     if (boundary->dim+1 == 2) suggestions = " You may need to add --download-triangle to your ./configure options";
189:     else if (boundary->dim+1 == 3) suggestions = " You may need to add --download-ctetgen or --download-tetgen in your ./configure options";
190:     SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"No grid generator of dimension %D registered%s",boundary->dim+1,suggestions);
191:   }
192:   return(0);
193: }

195: /*@C
196:   DMPlexGenerateRegister -  Adds a grid generator to DMPlex

198:    Not Collective

200:    Input Parameters:
201: +  name_solver - name of a new user-defined grid generator
202: .  fnc - generator function
203: .  rfnc - refinement function
204: -  dim - dimension of boundary of domain

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

209:    Sample usage:
210: .vb
211:    DMPlexGenerateRegister("my_generator",MyGeneratorCreate,MyGeneratorRefiner,dim);
212: .ve

214:    Then, your generator can be chosen with the procedural interface via
215: $     DMPlexGenerate(dm,"my_generator",...)
216:    or at runtime via the option
217: $     -dm_plex_generator my_generator

219:    Level: advanced

221: .seealso: DMPlexGenerateRegisterAll(), DMPlexGenerate(), DMPlexGenerateRegisterDestroy()

223: @*/
224: PetscErrorCode  DMPlexGenerateRegister(const char sname[],PetscErrorCode (*fnc)(DM, PetscBool,DM*), PetscErrorCode (*rfnc)(DM, double*,DM*),PetscInt dim)
225: {
226:   PetscErrorCode    ierr;
227:   PetscFunctionList entry;

230:   PetscNew(&entry);
231:   PetscStrallocpy(sname,&entry->name);
232:   entry->generate = fnc;
233:   entry->refine   = rfnc;
234:   entry->dim      = dim;
235:   entry->next     = NULL;
236:   if (!DMPlexGenerateList) DMPlexGenerateList = entry;
237:   else {
238:     PetscFunctionList fl = DMPlexGenerateList;
239:     while (fl->next) fl = fl->next;
240:     fl->next = entry;
241:   }
242:   return(0);
243: }

245: extern PetscBool DMPlexGenerateRegisterAllCalled;

247: PetscErrorCode  DMPlexGenerateRegisterDestroy(void)
248: {
249:   PetscFunctionList next,fl;
250:   PetscErrorCode    ierr;

253:   next = fl =  DMPlexGenerateList;
254:     while (next) {
255:     next = fl ? fl->next : NULL;
256:     if (fl) {PetscFree(fl->name);}
257:     PetscFree(fl);
258:     fl = next;
259:   }
260:   DMPlexGenerateList              = NULL;
261:   DMPlexGenerateRegisterAllCalled = PETSC_FALSE;
262:   return(0);
263: }