Actual source code: dmgenerate.c

  1: #include <petsc/private/dmimpl.h>

  3: PETSC_EXTERN PetscErrorCode DMIsForest(DM, PetscBool *);

  5: DMGeneratorFunctionList DMGenerateList              = NULL;
  6: PetscBool               DMGenerateRegisterAllCalled = PETSC_FALSE;

  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 *);

 32: /*@C
 33:   DMGenerateRegisterAll - Registers all of the mesh generation methods in the DM package.

 35:   Not Collective

 37:   Level: advanced

 39: .seealso: `DMGenerateRegisterDestroy()`
 40: @*/
 41: PetscErrorCode DMGenerateRegisterAll(void)
 42: {
 43:   if (DMGenerateRegisterAllCalled) return 0;
 44:   DMGenerateRegisterAllCalled = PETSC_TRUE;
 45: #if defined(PETSC_HAVE_TRIANGLE)
 46:   DMGenerateRegister("triangle", DMPlexGenerate_Triangle, DMPlexRefine_Triangle, NULL, 1);
 47: #endif
 48: #if defined(PETSC_HAVE_CTETGEN)
 49:   DMGenerateRegister("ctetgen", DMPlexGenerate_CTetgen, DMPlexRefine_CTetgen, NULL, 2);
 50: #endif
 51: #if defined(PETSC_HAVE_TETGEN)
 52:   DMGenerateRegister("tetgen", DMPlexGenerate_Tetgen, DMPlexRefine_Tetgen, NULL, 2);
 53: #endif
 54: #if defined(PETSC_HAVE_PRAGMATIC)
 55:   DMGenerateRegister("pragmatic", NULL, NULL, DMAdaptMetric_Pragmatic_Plex, -1);
 56: #endif
 57: #if defined(PETSC_HAVE_MMG)
 58:   DMGenerateRegister("mmg", NULL, NULL, DMAdaptMetric_Mmg_Plex, -1);
 59: #endif
 60: #if defined(PETSC_HAVE_PARMMG)
 61:   DMGenerateRegister("parmmg", NULL, NULL, DMAdaptMetric_ParMmg_Plex, -1);
 62: #endif
 63:   DMGenerateRegister("cellrefiner", NULL, NULL, DMPlexTransformAdaptLabel, -1);
 64:   DMGenerateRegister("forest", NULL, NULL, DMAdaptLabel_Forest, -1);
 65:   return 0;
 66: }

 68: /*@C
 69:   DMGenerateRegister -  Adds a grid generator to DM

 71:    Not Collective

 73:    Input Parameters:
 74: +  name_solver - name of a new user-defined grid generator
 75: .  fnc - generator function
 76: .  rfnc - refinement function
 77: .  alfnc - adapt by label function
 78: -  dim - dimension of boundary of domain

 80:    Notes:
 81:    DMGenerateRegister() may be called multiple times to add several user-defined solvers.

 83:    Sample usage:
 84: .vb
 85:    DMGenerateRegister("my_generator",MyGeneratorCreate,MyGeneratorRefiner,MyGeneratorAdaptor,dim);
 86: .ve

 88:    Then, your generator can be chosen with the procedural interface via
 89: $     DMGenerate(dm,"my_generator",...)
 90:    or at runtime via the option
 91: $     -dm_generator my_generator

 93:    Level: advanced

 95: .seealso: `DMGenerateRegisterAll()`, `DMPlexGenerate()`, `DMGenerateRegisterDestroy()`

 97: @*/
 98: PetscErrorCode DMGenerateRegister(const char sname[], PetscErrorCode (*fnc)(DM, PetscBool, DM *), PetscErrorCode (*rfnc)(DM, PetscReal *, DM *), PetscErrorCode (*alfnc)(DM, Vec, DMLabel, DMLabel, DM *), PetscInt dim)
 99: {
100:   DMGeneratorFunctionList entry;

102:   PetscNew(&entry);
103:   PetscStrallocpy(sname, &entry->name);
104:   entry->generate = fnc;
105:   entry->refine   = rfnc;
106:   entry->adapt    = alfnc;
107:   entry->dim      = dim;
108:   entry->next     = NULL;
109:   if (!DMGenerateList) DMGenerateList = entry;
110:   else {
111:     DMGeneratorFunctionList fl = DMGenerateList;
112:     while (fl->next) fl = fl->next;
113:     fl->next = entry;
114:   }
115:   return 0;
116: }

118: extern PetscBool DMGenerateRegisterAllCalled;

120: PetscErrorCode DMGenerateRegisterDestroy(void)
121: {
122:   DMGeneratorFunctionList next, fl;

124:   next = fl = DMGenerateList;
125:   while (next) {
126:     next = fl ? fl->next : NULL;
127:     if (fl) PetscFree(fl->name);
128:     PetscFree(fl);
129:     fl = next;
130:   }
131:   DMGenerateList              = NULL;
132:   DMGenerateRegisterAllCalled = PETSC_FALSE;
133:   return 0;
134: }

136: /*@C
137:   DMAdaptLabel - Adapt a dm based on a label with values interpreted as coarsening and refining flags.  Specific implementations of DM maybe have
138:                  specialized flags, but all implementations should accept flag values DM_ADAPT_DETERMINE, DM_ADAPT_KEEP, DM_ADAPT_REFINE, and DM_ADAPT_COARSEN.

140:   Collective on dm

142:   Input parameters:
143: + dm - the pre-adaptation DM object
144: - label - label with the flags

146:   Output parameters:
147: . dmAdapt - the adapted DM object: may be NULL if an adapted DM could not be produced.

149:   Level: intermediate

151: .seealso: `DMAdaptMetric()`, `DMCoarsen()`, `DMRefine()`
152: @*/
153: PetscErrorCode DMAdaptLabel(DM dm, DMLabel label, DM *dmAdapt)
154: {
155:   DMGeneratorFunctionList fl;
156:   char                    adaptname[PETSC_MAX_PATH_LEN];
157:   const char             *name;
158:   PetscInt                dim;
159:   PetscBool               flg, isForest, found = PETSC_FALSE;

164:   *dmAdapt = NULL;
165:   DMGetDimension(dm, &dim);
166:   DMIsForest(dm, &isForest);
167:   name = isForest ? "forest" : "cellrefiner";
168:   PetscOptionsGetString(((PetscObject)dm)->options, ((PetscObject)dm)->prefix, "-dm_adaptor", adaptname, sizeof(adaptname), &flg);
169:   if (flg) name = adaptname;

171:   fl = DMGenerateList;
172:   while (fl) {
173:     PetscStrcmp(fl->name, name, &flg);
174:     if (flg) {
175:       (*fl->adapt)(dm, NULL, label, NULL, dmAdapt);
176:       found = PETSC_TRUE;
177:     }
178:     fl = fl->next;
179:   }
181:   if (*dmAdapt) {
182:     (*dmAdapt)->prealloc_only = dm->prealloc_only; /* maybe this should go .... */
183:     PetscFree((*dmAdapt)->vectype);
184:     PetscStrallocpy(dm->vectype, (char **)&(*dmAdapt)->vectype);
185:     PetscFree((*dmAdapt)->mattype);
186:     PetscStrallocpy(dm->mattype, (char **)&(*dmAdapt)->mattype);
187:   }
188:   return 0;
189: }

191: /*@C
192:   DMAdaptMetric - Generates a mesh adapted to the specified metric field.

194:   Input Parameters:
195: + dm - The DM object
196: . metric - The metric to which the mesh is adapted, defined vertex-wise.
197: . 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_".
198: - 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_".

200:   Output Parameter:
201: . dmAdapt  - Pointer to the DM object containing the adapted mesh

203:   Note: The label in the adapted mesh will be registered under the name of the input DMLabel object

205:   Level: advanced

207: .seealso: `DMAdaptLabel()`, `DMCoarsen()`, `DMRefine()`
208: @*/
209: PetscErrorCode DMAdaptMetric(DM dm, Vec metric, DMLabel bdLabel, DMLabel rgLabel, DM *dmAdapt)
210: {
211:   DMGeneratorFunctionList fl;
212:   char                    adaptname[PETSC_MAX_PATH_LEN];
213:   const char             *name;
214:   const char *const       adaptors[3] = {"pragmatic", "mmg", "parmmg"};
215:   PetscInt                dim;
216:   PetscBool               flg, found = PETSC_FALSE;

223:   *dmAdapt = NULL;
224:   DMGetDimension(dm, &dim);
225:   PetscOptionsGetString(((PetscObject)dm)->options, ((PetscObject)dm)->prefix, "-dm_adaptor", adaptname, sizeof(adaptname), &flg);

227:   /* Default to Mmg in serial and ParMmg in parallel */
228:   if (flg) name = adaptname;
229:   else {
230:     MPI_Comm    comm;
231:     PetscMPIInt size;

233:     PetscObjectGetComm((PetscObject)dm, &comm);
234:     MPI_Comm_size(comm, &size);
235:     if (size == 1) name = adaptors[1];
236:     else name = adaptors[2];
237:   }

239:   fl = DMGenerateList;
240:   while (fl) {
241:     PetscStrcmp(fl->name, name, &flg);
242:     if (flg) {
243:       (*fl->adapt)(dm, metric, bdLabel, rgLabel, dmAdapt);
244:       found = PETSC_TRUE;
245:     }
246:     fl = fl->next;
247:   }
249:   if (*dmAdapt) {
250:     (*dmAdapt)->prealloc_only = dm->prealloc_only; /* maybe this should go .... */
251:     PetscFree((*dmAdapt)->vectype);
252:     PetscStrallocpy(dm->vectype, (char **)&(*dmAdapt)->vectype);
253:     PetscFree((*dmAdapt)->mattype);
254:     PetscStrallocpy(dm->mattype, (char **)&(*dmAdapt)->mattype);
255:   }
256:   return 0;
257: }