Actual source code: dacorn.c


  2: /*
  3:   Code for manipulating distributed regular arrays in parallel.
  4: */

  6: #include <petsc/private/dmdaimpl.h>
  7: #include <petscdmfield.h>

  9: PetscErrorCode DMCreateCoordinateDM_DA(DM dm, DM *cdm)
 10: {
 11:   DMDACreateCompatibleDMDA(dm, dm->dim, cdm);
 12:   return 0;
 13: }

 15: PetscErrorCode DMCreateCoordinateField_DA(DM dm, DMField *field)
 16: {
 17:   PetscReal   gmin[3], gmax[3];
 18:   PetscScalar corners[24];
 19:   PetscInt    dim;
 20:   PetscInt    i, j;
 21:   DM          cdm;

 23:   DMGetDimension(dm, &dim);
 24:   /* TODO: this is wrong if coordinates are not rectilinear */
 25:   DMGetBoundingBox(dm, gmin, gmax);
 26:   for (i = 0; i < (1 << dim); i++) {
 27:     for (j = 0; j < dim; j++) corners[i * dim + j] = (i & (1 << j)) ? gmax[j] : gmin[j];
 28:   }
 29:   DMClone(dm, &cdm);
 30:   DMFieldCreateDA(cdm, dim, corners, field);
 31:   DMDestroy(&cdm);
 32:   return 0;
 33: }

 35: /*@C
 36:    DMDASetFieldName - Sets the names of individual field components in multicomponent
 37:    vectors associated with a `DMDA`.

 39:    Logically collective; name must contain a common value

 41:    Input Parameters:
 42: +  da - the distributed array
 43: .  nf - field number for the `DMDA` (0, 1, ... dof-1), where dof indicates the
 44:         number of degrees of freedom per node within the `DMDA`
 45: -  names - the name of the field (component)

 47:   Level: intermediate

 49:   Note:
 50:     It must be called after having called `DMSetUp()`.

 52: .seealso: `DM`, `DMDA`, `DMDAGetFieldName()`, `DMDASetCoordinateName()`, `DMDAGetCoordinateName()`, `DMDASetFieldNames()`, `DMSetUp()`
 53: @*/
 54: PetscErrorCode DMDASetFieldName(DM da, PetscInt nf, const char name[])
 55: {
 56:   DM_DA *dd = (DM_DA *)da->data;

 61:   PetscFree(dd->fieldname[nf]);
 62:   PetscStrallocpy(name, &dd->fieldname[nf]);
 63:   return 0;
 64: }

 66: /*@C
 67:    DMDAGetFieldNames - Gets the name of each component in the vector associated with the `DMDA`

 69:    Not collective; names will contain a common value

 71:    Input Parameter:
 72: .  dm - the `DMDA` object

 74:    Output Parameter:
 75: .  names - the names of the components, final string is NULL, will have the same number of entries as the dof used in creating the `DMDA`

 77:    Level: intermediate

 79:    Fortran Note:
 80:    Not supported from Fortran, use `DMDAGetFieldName()`

 82: .seealso: `DM`, `DMDA`, `DMDAGetFieldName()`, `DMDASetCoordinateName()`, `DMDAGetCoordinateName()`, `DMDASetFieldName()`, `DMDASetFieldNames()`
 83: @*/
 84: PetscErrorCode DMDAGetFieldNames(DM da, const char *const **names)
 85: {
 86:   DM_DA *dd = (DM_DA *)da->data;

 88:   *names = (const char *const *)dd->fieldname;
 89:   return 0;
 90: }

 92: /*@C
 93:    DMDASetFieldNames - Sets the name of each component in the vector associated with the DMDA

 95:    Logically collective; names must contain a common value

 97:    Input Parameters:
 98: +  dm - the `DMDA` object
 99: -  names - the names of the components, final string must be NULL, must have the same number of entries as the dof used in creating the `DMDA`

101:    Level: intermediate

103:    Note:
104:     It must be called after having called `DMSetUp()`.

106:    Fortran Note:
107:    Not supported from Fortran, use `DMDASetFieldName()`

109: .seealso: `DM`, `DMDA`, `DMDAGetFieldName()`, `DMDASetCoordinateName()`, `DMDAGetCoordinateName()`, `DMDASetFieldName()`, `DMSetUp()`
110: @*/
111: PetscErrorCode DMDASetFieldNames(DM da, const char *const *names)
112: {
113:   DM_DA   *dd = (DM_DA *)da->data;
114:   char   **fieldname;
115:   PetscInt nf = 0;

118:   while (names[nf++]) { };
120:   PetscStrArrayallocpy(names, &fieldname);
121:   PetscStrArrayDestroy(&dd->fieldname);
122:   dd->fieldname = fieldname;
123:   return 0;
124: }

126: /*@C
127:    DMDAGetFieldName - Gets the names of individual field components in multicomponent
128:    vectors associated with a `DMDA`.

130:    Not collective; name will contain a common value

132:    Input Parameters:
133: +  da - the distributed array
134: -  nf - field number for the `DMDA` (0, 1, ... dof-1), where dof indicates the
135:         number of degrees of freedom per node within the `DMDA`

137:    Output Parameter:
138: .  names - the name of the field (component)

140:   Level: intermediate

142:   Note:
143:     It must be called after having called `DMSetUp()`.

145: .seealso: `DM`, `DMDA`, `DMDASetFieldName()`, `DMDASetCoordinateName()`, `DMDAGetCoordinateName()`, `DMSetUp()`
146: @*/
147: PetscErrorCode DMDAGetFieldName(DM da, PetscInt nf, const char **name)
148: {
149:   DM_DA *dd = (DM_DA *)da->data;

155:   *name = dd->fieldname[nf];
156:   return 0;
157: }

159: /*@C
160:    DMDASetCoordinateName - Sets the name of the coordinate directions associated with a `DMDA`, for example "x" or "y"

162:    Logically collective; name must contain a common value

164:    Input Parameters:
165: +  dm - the `DMDA`
166: .  nf - coordinate number for the DMDA (0, 1, ... dim-1),
167: -  name - the name of the coordinate

169:   Level: intermediate

171:   Note:
172:     It must be called after having called `DMSetUp()`.

174:   Fortran Note:
175:   Not supported from Fortran

177: .seealso: `DM`, `DMDA`, `DMDAGetCoordinateName()`, `DMDASetFieldName()`, `DMDAGetFieldName()`, `DMSetUp()`
178: @*/
179: PetscErrorCode DMDASetCoordinateName(DM dm, PetscInt nf, const char name[])
180: {
181:   DM_DA *dd = (DM_DA *)dm->data;

186:   PetscFree(dd->coordinatename[nf]);
187:   PetscStrallocpy(name, &dd->coordinatename[nf]);
188:   return 0;
189: }

191: /*@C
192:    DMDAGetCoordinateName - Gets the name of a coodinate direction associated with a `DMDA`.

194:    Not collective; name will contain a common value

196:    Input Parameters:
197: +  dm - the `DMDA`
198: -  nf -  number for the DMDA (0, 1, ... dim-1)

200:    Output Parameter:
201: .  names - the name of the coordinate direction

203:   Level: intermediate

205:   Note:
206:     It must be called after having called `DMSetUp()`.

208:   Fortran Note:
209:   Not supported from Fortran

211: .seealso: `DM`, `DMDA`, `DMDASetCoordinateName()`, `DMDASetFieldName()`, `DMDAGetFieldName()`, `DMSetUp()`
212: @*/
213: PetscErrorCode DMDAGetCoordinateName(DM dm, PetscInt nf, const char **name)
214: {
215:   DM_DA *dd = (DM_DA *)dm->data;

221:   *name = dd->coordinatename[nf];
222:   return 0;
223: }

225: /*@C
226:    DMDAGetCorners - Returns the global (x,y,z) indices of the lower left
227:    corner and size of the local region, excluding ghost points.

229:    Not collective

231:    Input Parameter:
232: .  da - the distributed array

234:    Output Parameters:
235: +  x - the corner index for the first dimension
236: .  y - the corner index for the second dimension (only used in 2D and 3D problems)
237: .  z - the corner index for the third dimension (only used in 3D problems)
238: .  m - the width in the first dimension
239: .  n - the width in the second dimension (only used in 2D and 3D problems)
240: -  p - the width in the third dimension (only used in 3D problems)

242:   Level: beginner

244:    Note:
245:    The corner information is independent of the number of degrees of
246:    freedom per node set with the `DMDACreateXX()` routine. Thus the x, y, z, and
247:    m, n, p can be thought of as coordinates on a logical grid, where each
248:    grid point has (potentially) several degrees of freedom.
249:    Any of y, z, n, and p can be passed in as NULL if not needed.

251: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetOwnershipRanges()`, `DMStagGetCorners()`
252: @*/
253: PetscErrorCode DMDAGetCorners(DM da, PetscInt *x, PetscInt *y, PetscInt *z, PetscInt *m, PetscInt *n, PetscInt *p)
254: {
255:   PetscInt w;
256:   DM_DA   *dd = (DM_DA *)da->data;

259:   /* since the xs, xe ... have all been multiplied by the number of degrees
260:      of freedom per cell, w = dd->w, we divide that out before returning.*/
261:   w = dd->w;
262:   if (x) *x = dd->xs / w + dd->xo;
263:   /* the y and z have NOT been multiplied by w */
264:   if (y) *y = dd->ys + dd->yo;
265:   if (z) *z = dd->zs + dd->zo;
266:   if (m) *m = (dd->xe - dd->xs) / w;
267:   if (n) *n = (dd->ye - dd->ys);
268:   if (p) *p = (dd->ze - dd->zs);
269:   return 0;
270: }

272: PetscErrorCode DMGetLocalBoundingIndices_DMDA(DM dm, PetscReal lmin[], PetscReal lmax[])
273: {
274:   DMDALocalInfo info;

276:   DMDAGetLocalInfo(dm, &info);
277:   lmin[0] = info.xs;
278:   lmin[1] = info.ys;
279:   lmin[2] = info.zs;
280:   lmax[0] = info.xs + info.xm - 1;
281:   lmax[1] = info.ys + info.ym - 1;
282:   lmax[2] = info.zs + info.zm - 1;
283:   return 0;
284: }

286: /*@
287:    DMDAGetReducedDMDA - Deprecated; use DMDACreateCompatibleDMDA()

289:    Level: deprecated
290: @*/
291: PetscErrorCode DMDAGetReducedDMDA(DM da, PetscInt nfields, DM *nda)
292: {
293:   DMDACreateCompatibleDMDA(da, nfields, nda);
294:   return 0;
295: }

297: /*@
298:    DMDACreateCompatibleDMDA - Creates a `DMDA` with the same layout but with fewer or more fields

300:    Collective

302:    Input Parameters:
303: +  da - the distributed array
304: -  nfields - number of fields in new `DMDA`

306:    Output Parameter:
307: .  nda - the new `DMDA`

309:   Level: intermediate

311: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMSetCoordinates()`, `DMDASetUniformCoordinates()`, `DMGetCoordinates()`, `DMDAGetGhostedCoordinates()`, `DMStagCreateCompatibleDMStag()`
312: @*/
313: PetscErrorCode DMDACreateCompatibleDMDA(DM da, PetscInt nfields, DM *nda)
314: {
315:   DM_DA          *dd = (DM_DA *)da->data;
316:   PetscInt        s, m, n, p, M, N, P, dim, Mo, No, Po;
317:   const PetscInt *lx, *ly, *lz;
318:   DMBoundaryType  bx, by, bz;
319:   DMDAStencilType stencil_type;
320:   Vec             coords;
321:   PetscInt        ox, oy, oz;
322:   PetscInt        cl, rl;

324:   dim = da->dim;
325:   M   = dd->M;
326:   N   = dd->N;
327:   P   = dd->P;
328:   m   = dd->m;
329:   n   = dd->n;
330:   p   = dd->p;
331:   s   = dd->s;
332:   bx  = dd->bx;
333:   by  = dd->by;
334:   bz  = dd->bz;

336:   stencil_type = dd->stencil_type;

338:   DMDAGetOwnershipRanges(da, &lx, &ly, &lz);
339:   if (dim == 1) {
340:     DMDACreate1d(PetscObjectComm((PetscObject)da), bx, M, nfields, s, dd->lx, nda);
341:   } else if (dim == 2) {
342:     DMDACreate2d(PetscObjectComm((PetscObject)da), bx, by, stencil_type, M, N, m, n, nfields, s, lx, ly, nda);
343:   } else if (dim == 3) {
344:     DMDACreate3d(PetscObjectComm((PetscObject)da), bx, by, bz, stencil_type, M, N, P, m, n, p, nfields, s, lx, ly, lz, nda);
345:   }
346:   DMSetUp(*nda);
347:   DMGetCoordinates(da, &coords);
348:   DMSetCoordinates(*nda, coords);

350:   /* allow for getting a reduced DA corresponding to a domain decomposition */
351:   DMDAGetOffset(da, &ox, &oy, &oz, &Mo, &No, &Po);
352:   DMDASetOffset(*nda, ox, oy, oz, Mo, No, Po);

354:   /* allow for getting a reduced DA corresponding to a coarsened DA */
355:   DMGetCoarsenLevel(da, &cl);
356:   DMGetRefineLevel(da, &rl);

358:   (*nda)->levelup   = rl;
359:   (*nda)->leveldown = cl;
360:   return 0;
361: }

363: /*@C
364:    DMDAGetCoordinateArray - Gets an array containing the coordinates of the `DMDA`

366:    Not collective

368:    Input Parameter:
369: .  dm - the `DMDA`

371:    Output Parameter:
372: .  xc - the coordinates

374:   Level: intermediate

376:   Fortran Note:
377:   Not supported from Fortran

379: .seealso: `DM`, `DMDA`, `DMDASetCoordinateName()`, `DMDASetFieldName()`, `DMDAGetFieldName()`, `DMDARestoreCoordinateArray()`
380: @*/
381: PetscErrorCode DMDAGetCoordinateArray(DM dm, void *xc)
382: {
383:   DM  cdm;
384:   Vec x;

387:   DMGetCoordinates(dm, &x);
388:   DMGetCoordinateDM(dm, &cdm);
389:   DMDAVecGetArray(cdm, x, xc);
390:   return 0;
391: }

393: /*@C
394:    DMDARestoreCoordinateArray - Sets an array containing the coordinates of the `DMDA`

396:    Not collective

398:    Input Parameters:
399: +  dm - the `DMDA`
400: -  xc - the coordinates

402:   Level: intermediate

404:   Fortran Note:
405:   Not supported from Fortran

407: .seealso: `DM`, `DMDA`, `DMDASetCoordinateName()`, `DMDASetFieldName()`, `DMDAGetFieldName()`, `DMDAGetCoordinateArray()`
408: @*/
409: PetscErrorCode DMDARestoreCoordinateArray(DM dm, void *xc)
410: {
411:   DM  cdm;
412:   Vec x;

415:   DMGetCoordinates(dm, &x);
416:   DMGetCoordinateDM(dm, &cdm);
417:   DMDAVecRestoreArray(cdm, x, xc);
418:   return 0;
419: }