Actual source code: dacorn.c
petsc-3.9.4 2018-09-11
2: /*
3: Code for manipulating distributed regular arrays in parallel.
4: */
6: #include <petsc/private/dmdaimpl.h>
8: PetscErrorCode DMCreateCoordinateDM_DA(DM dm, DM *cdm)
9: {
12: DMDAGetReducedDMDA(dm,dm->dim,cdm);
13: return(0);
14: }
16: /*@C
17: DMDASetFieldName - Sets the names of individual field components in multicomponent
18: vectors associated with a DMDA.
20: Not Collective
22: Input Parameters:
23: + da - the distributed array
24: . nf - field number for the DMDA (0, 1, ... dof-1), where dof indicates the
25: number of degrees of freedom per node within the DMDA
26: - names - the name of the field (component)
28: Notes: It must be called after having called DMSetUp().
30: Level: intermediate
32: .keywords: distributed array, get, component name
34: .seealso: DMDAGetFieldName(), DMDASetCoordinateName(), DMDAGetCoordinateName(), DMDASetFieldNames(), DMSetUp()
35: @*/
36: PetscErrorCode DMDASetFieldName(DM da,PetscInt nf,const char name[])
37: {
39: DM_DA *dd = (DM_DA*)da->data;
43: if (nf < 0 || nf >= dd->w) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Invalid field number: %D",nf);
44: if (!dd->fieldname) SETERRQ(PetscObjectComm((PetscObject)da),PETSC_ERR_ORDER,"You should call DMSetUp() first");
45: PetscFree(dd->fieldname[nf]);
46: PetscStrallocpy(name,&dd->fieldname[nf]);
47: return(0);
48: }
50: /*@C
51: DMDAGetFieldNames - Gets the name of each component in the vector associated with the DMDA
53: Collective on TS
55: Input Parameter:
56: . dm - the DMDA object
58: Output Parameter:
59: . 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
61: Level: intermediate
63: Not supported from Fortran, use DMDAGetFieldName()
65: .keywords: distributed array, get, component name
67: .seealso: DMDAGetFieldName(), DMDASetCoordinateName(), DMDAGetCoordinateName(), DMDASetFieldName(), DMDASetFieldNames()
68: @*/
69: PetscErrorCode DMDAGetFieldNames(DM da,const char * const **names)
70: {
71: DM_DA *dd = (DM_DA*)da->data;
74: *names = (const char * const *) dd->fieldname;
75: return(0);
76: }
78: /*@C
79: DMDASetFieldNames - Sets the name of each component in the vector associated with the DMDA
81: Collective on TS
83: Input Parameters:
84: + dm - the DMDA object
85: - 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
87: Notes: It must be called after having called DMSetUp().
89: Level: intermediate
91: Not supported from Fortran, use DMDASetFieldName()
93: .keywords: distributed array, get, component name
95: .seealso: DMDAGetFieldName(), DMDASetCoordinateName(), DMDAGetCoordinateName(), DMDASetFieldName(), DMSetUp()
96: @*/
97: PetscErrorCode DMDASetFieldNames(DM da,const char * const *names)
98: {
100: DM_DA *dd = (DM_DA*)da->data;
101: char **fieldname;
102: PetscInt nf = 0;
105: if (!dd->fieldname) SETERRQ(PetscObjectComm((PetscObject)da),PETSC_ERR_ORDER,"You should call DMSetUp() first");
106: while (names[nf++]) {};
107: if (nf != dd->w+1) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Invalid number of fields %D",nf-1);
108: PetscStrArrayallocpy(names,&fieldname);
109: PetscStrArrayDestroy(&dd->fieldname);
110: dd->fieldname = fieldname;
111: return(0);
112: }
114: /*@C
115: DMDAGetFieldName - Gets the names of individual field components in multicomponent
116: vectors associated with a DMDA.
118: Not Collective
120: Input Parameter:
121: + da - the distributed array
122: - nf - field number for the DMDA (0, 1, ... dof-1), where dof indicates the
123: number of degrees of freedom per node within the DMDA
125: Output Parameter:
126: . names - the name of the field (component)
128: Notes: It must be called after having called DMSetUp().
130: Level: intermediate
132: .keywords: distributed array, get, component name
134: .seealso: DMDASetFieldName(), DMDASetCoordinateName(), DMDAGetCoordinateName(), DMSetUp()
135: @*/
136: PetscErrorCode DMDAGetFieldName(DM da,PetscInt nf,const char **name)
137: {
138: DM_DA *dd = (DM_DA*)da->data;
143: if (nf < 0 || nf >= dd->w) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Invalid field number: %D",nf);
144: if (!dd->fieldname) SETERRQ(PetscObjectComm((PetscObject)da),PETSC_ERR_ORDER,"You should call DMSetUp() first");
145: *name = dd->fieldname[nf];
146: return(0);
147: }
149: /*@C
150: DMDASetCoordinateName - Sets the name of the coordinate directions associated with a DMDA, for example "x" or "y"
152: Not Collective
154: Input Parameters:
155: + dm - the DM
156: . nf - coordinate number for the DMDA (0, 1, ... dim-1),
157: - name - the name of the coordinate
159: Notes: It must be called after having called DMSetUp().
161: Level: intermediate
163: Not supported from Fortran
165: .keywords: distributed array, get, component name
167: .seealso: DMDAGetCoordinateName(), DMDASetFieldName(), DMDAGetFieldName(), DMSetUp()
168: @*/
169: PetscErrorCode DMDASetCoordinateName(DM dm,PetscInt nf,const char name[])
170: {
172: DM_DA *dd = (DM_DA*)dm->data;
176: if (nf < 0 || nf >= dm->dim) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Invalid coordinate number: %D",nf);
177: if (!dd->coordinatename) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ORDER,"You should call DMSetUp() first");
178: PetscFree(dd->coordinatename[nf]);
179: PetscStrallocpy(name,&dd->coordinatename[nf]);
180: return(0);
181: }
183: /*@C
184: DMDAGetCoordinateName - Gets the name of a coodinate direction associated with a DMDA.
186: Not Collective
188: Input Parameter:
189: + dm - the DM
190: - nf - number for the DMDA (0, 1, ... dim-1)
192: Output Parameter:
193: . names - the name of the coordinate direction
195: Notes: It must be called after having called DMSetUp().
197: Level: intermediate
199: Not supported from Fortran
201: .keywords: distributed array, get, component name
203: .seealso: DMDASetCoordinateName(), DMDASetFieldName(), DMDAGetFieldName(), DMSetUp()
204: @*/
205: PetscErrorCode DMDAGetCoordinateName(DM dm,PetscInt nf,const char **name)
206: {
207: DM_DA *dd = (DM_DA*)dm->data;
212: if (nf < 0 || nf >= dm->dim) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Invalid coordinate number: %D",nf);
213: if (!dd->coordinatename) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ORDER,"You should call DMSetUp() first");
214: *name = dd->coordinatename[nf];
215: return(0);
216: }
218: /*@C
219: DMDAGetCorners - Returns the global (x,y,z) indices of the lower left
220: corner and size of the local region, excluding ghost points.
222: Not Collective
224: Input Parameter:
225: . da - the distributed array
227: Output Parameters:
228: + x,y,z - the corner indices (where y and z are optional; these are used
229: for 2D and 3D problems)
230: - m,n,p - widths in the corresponding directions (where n and p are optional;
231: these are used for 2D and 3D problems)
233: Note:
234: The corner information is independent of the number of degrees of
235: freedom per node set with the DMDACreateXX() routine. Thus the x, y, z, and
236: m, n, p can be thought of as coordinates on a logical grid, where each
237: grid point has (potentially) several degrees of freedom.
238: Any of y, z, n, and p can be passed in as NULL if not needed.
240: Level: beginner
242: .keywords: distributed array, get, corners, nodes, local indices
244: .seealso: DMDAGetGhostCorners(), DMDAGetOwnershipRanges()
245: @*/
246: PetscErrorCode DMDAGetCorners(DM da,PetscInt *x,PetscInt *y,PetscInt *z,PetscInt *m,PetscInt *n,PetscInt *p)
247: {
248: PetscInt w;
249: DM_DA *dd = (DM_DA*)da->data;
253: /* since the xs, xe ... have all been multiplied by the number of degrees
254: of freedom per cell, w = dd->w, we divide that out before returning.*/
255: w = dd->w;
256: if (x) *x = dd->xs/w + dd->xo;
257: /* the y and z have NOT been multiplied by w */
258: if (y) *y = dd->ys + dd->yo;
259: if (z) *z = dd->zs + dd->zo;
260: if (m) *m = (dd->xe - dd->xs)/w;
261: if (n) *n = (dd->ye - dd->ys);
262: if (p) *p = (dd->ze - dd->zs);
263: return(0);
264: }
266: /*@
267: DMDAGetLocalBoundingBox - Returns the local bounding box for the DMDA.
269: Not Collective
271: Input Parameter:
272: . dm - the DM
274: Output Parameters:
275: + lmin - local minimum coordinates (length dim, optional)
276: - lmax - local maximim coordinates (length dim, optional)
278: Level: beginner
280: Not supported from Fortran
282: .keywords: distributed array, get, coordinates
284: .seealso: DMDAGetCoordinateDA(), DMGetCoordinates(), DMDAGetBoundingBox()
285: @*/
286: PetscErrorCode DMDAGetLocalBoundingBox(DM dm,PetscReal lmin[],PetscReal lmax[])
287: {
288: PetscErrorCode ierr;
289: Vec coords = NULL;
290: PetscInt dim,i,j;
291: const PetscScalar *local_coords;
292: PetscReal min[3]={PETSC_MAX_REAL,PETSC_MAX_REAL,PETSC_MAX_REAL},max[3]={PETSC_MIN_REAL,PETSC_MIN_REAL,PETSC_MIN_REAL};
293: PetscInt N,Ni;
297: dim = dm->dim;
298: DMGetCoordinates(dm,&coords);
299: if (coords) {
300: VecGetArrayRead(coords,&local_coords);
301: VecGetLocalSize(coords,&N);
302: Ni = N/dim;
303: for (i=0; i<Ni; i++) {
304: for (j=0; j<3; j++) {
305: min[j] = j < dim ? PetscMin(min[j],PetscRealPart(local_coords[i*dim+j])) : 0;
306: max[j] = j < dim ? PetscMax(max[j],PetscRealPart(local_coords[i*dim+j])) : 0;
307: }
308: }
309: VecRestoreArrayRead(coords,&local_coords);
310: } else { /* Just use grid indices */
311: DMDALocalInfo info;
312: DMDAGetLocalInfo(dm,&info);
313: min[0] = info.xs;
314: min[1] = info.ys;
315: min[2] = info.zs;
316: max[0] = info.xs + info.xm-1;
317: max[1] = info.ys + info.ym-1;
318: max[2] = info.zs + info.zm-1;
319: }
320: if (lmin) {PetscMemcpy(lmin,min,dim*sizeof(PetscReal));}
321: if (lmax) {PetscMemcpy(lmax,max,dim*sizeof(PetscReal));}
322: return(0);
323: }
325: /*@
326: DMDAGetBoundingBox - Returns the global bounding box for the DMDA.
328: Collective on DMDA
330: Input Parameter:
331: . dm - the DM
333: Output Parameters:
334: + gmin - global minimum coordinates (length dim, optional)
335: - gmax - global maximim coordinates (length dim, optional)
337: Level: beginner
339: .keywords: distributed array, get, coordinates
341: .seealso: DMDAGetCoordinateDA(), DMGetCoordinates(), DMDAGetLocalBoundingBox()
342: @*/
343: PetscErrorCode DMDAGetBoundingBox(DM dm,PetscReal gmin[],PetscReal gmax[])
344: {
346: PetscMPIInt count;
347: PetscReal lmin[3],lmax[3];
351: PetscMPIIntCast(dm->dim,&count);
352: DMDAGetLocalBoundingBox(dm,lmin,lmax);
353: if (gmin) {MPIU_Allreduce(lmin,gmin,count,MPIU_REAL,MPIU_MIN,PetscObjectComm((PetscObject)dm));}
354: if (gmax) {MPIU_Allreduce(lmax,gmax,count,MPIU_REAL,MPIU_MAX,PetscObjectComm((PetscObject)dm));}
355: return(0);
356: }
358: /*@
359: DMDAGetReducedDMDA - Gets the DMDA with the same layout but with fewer or more fields
361: Collective on DMDA
363: Input Parameters:
364: + da - the distributed array
365: - nfields - number of fields in new DMDA
367: Output Parameter:
368: . nda - the new DMDA
370: Level: intermediate
372: .keywords: distributed array, get, corners, nodes, local indices, coordinates
374: .seealso: DMDAGetGhostCorners(), DMSetCoordinates(), DMDASetUniformCoordinates(), DMGetCoordinates(), DMDAGetGhostedCoordinates()
375: @*/
376: PetscErrorCode DMDAGetReducedDMDA(DM da,PetscInt nfields,DM *nda)
377: {
378: PetscErrorCode ierr;
379: DM_DA *dd = (DM_DA*)da->data;
380: PetscInt s,m,n,p,M,N,P,dim,Mo,No,Po;
381: const PetscInt *lx,*ly,*lz;
382: DMBoundaryType bx,by,bz;
383: DMDAStencilType stencil_type;
384: PetscInt ox,oy,oz;
385: PetscInt cl,rl;
388: dim = da->dim;
389: M = dd->M;
390: N = dd->N;
391: P = dd->P;
392: m = dd->m;
393: n = dd->n;
394: p = dd->p;
395: s = dd->s;
396: bx = dd->bx;
397: by = dd->by;
398: bz = dd->bz;
400: stencil_type = dd->stencil_type;
402: DMDAGetOwnershipRanges(da,&lx,&ly,&lz);
403: if (dim == 1) {
404: DMDACreate1d(PetscObjectComm((PetscObject)da),bx,M,nfields,s,dd->lx,nda);
405: } else if (dim == 2) {
406: DMDACreate2d(PetscObjectComm((PetscObject)da),bx,by,stencil_type,M,N,m,n,nfields,s,lx,ly,nda);
407: } else if (dim == 3) {
408: DMDACreate3d(PetscObjectComm((PetscObject)da),bx,by,bz,stencil_type,M,N,P,m,n,p,nfields,s,lx,ly,lz,nda);
409: }
410: DMSetUp(*nda);
411: if (da->coordinates) {
412: PetscObjectReference((PetscObject)da->coordinates);
413: (*nda)->coordinates = da->coordinates;
414: }
416: /* allow for getting a reduced DA corresponding to a domain decomposition */
417: DMDAGetOffset(da,&ox,&oy,&oz,&Mo,&No,&Po);
418: DMDASetOffset(*nda,ox,oy,oz,Mo,No,Po);
420: /* allow for getting a reduced DA corresponding to a coarsened DA */
421: DMGetCoarsenLevel(da,&cl);
422: DMGetRefineLevel(da,&rl);
424: (*nda)->levelup = rl;
425: (*nda)->leveldown = cl;
426: return(0);
427: }
429: /*@C
430: DMDAGetCoordinateArray - Gets an array containing the coordinates of the DMDA
432: Not Collective
434: Input Parameter:
435: . dm - the DM
437: Output Parameter:
438: . xc - the coordinates
440: Level: intermediate
442: Not supported from Fortran
444: .keywords: distributed array, get, component name
446: .seealso: DMDASetCoordinateName(), DMDASetFieldName(), DMDAGetFieldName(), DMDARestoreCoordinateArray()
447: @*/
448: PetscErrorCode DMDAGetCoordinateArray(DM dm,void *xc)
449: {
451: DM cdm;
452: Vec x;
456: DMGetCoordinates(dm,&x);
457: DMGetCoordinateDM(dm,&cdm);
458: DMDAVecGetArray(cdm,x,xc);
459: return(0);
460: }
462: /*@C
463: DMDARestoreCoordinateArray - Sets an array containing the coordinates of the DMDA
465: Not Collective
467: Input Parameter:
468: + dm - the DM
469: - xc - the coordinates
471: Level: intermediate
473: Not supported from Fortran
475: .keywords: distributed array, get, component name
477: .seealso: DMDASetCoordinateName(), DMDASetFieldName(), DMDAGetFieldName(), DMDAGetCoordinateArray()
478: @*/
479: PetscErrorCode DMDARestoreCoordinateArray(DM dm,void *xc)
480: {
482: DM cdm;
483: Vec x;
487: DMGetCoordinates(dm,&x);
488: DMGetCoordinateDM(dm,&cdm);
489: DMDAVecRestoreArray(cdm,x,xc);
490: return(0);
491: }