Actual source code: dacorn.c
1: #define PETSCDM_DLL
3: /*
4: Code for manipulating distributed regular arrays in parallel.
5: */
7: #include src/dm/da/daimpl.h
11: /*@
12: DASetCoordinates - Sets into the DA a vector that indicates the
13: coordinates of the local nodes (NOT including ghost nodes).
15: Not Collective
17: Input Parameter:
18: + da - the distributed array
19: - c - coordinate vector
21: Note:
22: The coordinates should NOT include those for all ghost points
24: Does NOT increase the reference count of this vector, so caller should NOT
25: destroy the vector.
27: Level: intermediate
29: .keywords: distributed array, get, corners, nodes, local indices, coordinates
31: .seealso: DAGetGhostCorners(), DAGetCoordinates(), DASetUniformCoordinates(). DAGetGhostCoordinates(), DAGetCoordinateDA()
32: @*/
33: PetscErrorCode DASetCoordinates(DA da,Vec c)
34: {
40: da->coordinates = c;
41: VecSetBlockSize(c,da->dim);
42: return(0);
43: }
47: /*@
48: DAGetCoordinates - Gets the node coordinates associated with a DA.
50: Not Collective
52: Input Parameter:
53: . da - the distributed array
55: Output Parameter:
56: . c - coordinate vector
58: Note:
59: Each process has only the coordinates for its local nodes (does NOT have the
60: coordinates for the ghost nodes).
62: For two and three dimensions coordinates are interlaced (x_0,y_0,x_1,y_1,...)
63: and (x_0,y_0,z_0,x_1,y_1,z_1...)
65: You should not destroy or keep around this vector after the DA is destroyed.
67: Level: intermediate
69: .keywords: distributed array, get, corners, nodes, local indices, coordinates
71: .seealso: DAGetGhostCorners(), DASetCoordinates(), DASetUniformCoordinates(), DAGetCoordinates(), DAGetCoordinateDA()
72: @*/
73: PetscErrorCode DAGetCoordinates(DA da,Vec *c)
74: {
76:
79: *c = da->coordinates;
80: return(0);
81: }
85: /*@
86: DAGetCoordinateDA - Gets the DA that scatters between global and local DA coordinates
88: Collective on DA
90: Input Parameter:
91: . da - the distributed array
93: Output Parameter:
94: . dac - coordinate DA
96: Note: You should not destroy or keep around this vector after the DA is destroyed.
98: Level: intermediate
100: .keywords: distributed array, get, corners, nodes, local indices, coordinates
102: .seealso: DAGetGhostCorners(), DASetCoordinates(), DASetUniformCoordinates(), DAGetCoordinates(), DAGetGhostedCoordinates()
103: @*/
104: PetscErrorCode DAGetCoordinateDA(DA da,DA *cda)
105: {
106: DAStencilType st;
108: PetscMPIInt size;
111: if (da->da_coordinates) {
112: *cda = da->da_coordinates;
113: return(0);
114: }
115: MPI_Comm_size(da->comm,&size);
116: if (da->dim == 1) {
117: if (da->w == 1) {
118: da->da_coordinates = da;
119: } else {
120: PetscInt s,m,*lc,l;
121: DAPeriodicType pt;
122: DAGetInfo(da,0,&m,0,0,0,0,0,0,&s,&pt,0);
123: DAGetCorners(da,0,0,0,&l,0,0);
124: PetscMalloc(size*sizeof(PetscInt),&lc);
125: MPI_Allgather(&l,1,MPIU_INT,lc,1,MPIU_INT,da->comm);
126: DACreate1d(da->comm,pt,m,1,s,lc,&da->da_coordinates);
127: PetscFree(lc);
128: }
129: } else if (da->dim == 2) {
130: DAGetInfo(da,0,0,0,0,0,0,0,0,0,0,&st);
131: if (da->w == 2 && st == DA_STENCIL_BOX) {
132: da->da_coordinates = da;
133: } else {
134: PetscInt i,s,m,*lc,*ld,l,k,n,M,N;
135: DAPeriodicType pt;
136: DAGetInfo(da,0,&m,&n,0,&M,&N,0,0,&s,&pt,0);
137: DAGetCorners(da,0,0,0,&l,&k,0);
138: PetscMalloc(size*sizeof(PetscInt),&lc);
139: PetscMalloc(size*sizeof(PetscInt),&ld);
140: /* only first M values in lc matter */
141: MPI_Allgather(&l,1,MPIU_INT,lc,1,MPIU_INT,da->comm);
142: /* every Mth value in ld matters */
143: MPI_Allgather(&k,1,MPIU_INT,ld,1,MPIU_INT,da->comm);
144: for ( i=0; i<N; i++) {
145: ld[i] = ld[M*i];
146: }
147: DACreate2d(da->comm,pt,DA_STENCIL_BOX,m,n,M,N,2,s,lc,ld,&da->da_coordinates);
148: PetscFree(lc);
149: PetscFree(ld);
150: }
151: } else if (da->dim == 3) {
152: DAGetInfo(da,0,0,0,0,0,0,0,0,0,0,&st);
153: if (da->w == 3 && st == DA_STENCIL_BOX) {
154: da->da_coordinates = da;
155: } else {
156: PetscInt i,s,m,*lc,*ld,*le,l,k,q,n,M,N,P,p;
157: DAPeriodicType pt;
158: DAGetInfo(da,0,&m,&n,&p,&M,&N,&P,0,&s,&pt,0);
159: DAGetCorners(da,0,0,0,&l,&k,&q);
160: PetscMalloc(size*sizeof(PetscInt),&lc);
161: PetscMalloc(size*sizeof(PetscInt),&ld);
162: PetscMalloc(size*sizeof(PetscInt),&le);
163: /* only first M values in lc matter */
164: MPI_Allgather(&l,1,MPIU_INT,lc,1,MPIU_INT,da->comm);
165: /* every Mth value in ld matters */
166: MPI_Allgather(&k,1,MPIU_INT,ld,1,MPIU_INT,da->comm);
167: for ( i=0; i<N; i++) {
168: ld[i] = ld[M*i];
169: }
170: MPI_Allgather(&q,1,MPIU_INT,le,1,MPIU_INT,da->comm);
171: for ( i=0; i<P; i++) {
172: le[i] = le[M*N*i];
173: }
174: DACreate3d(da->comm,pt,DA_STENCIL_BOX,m,n,p,M,N,P,3,s,lc,ld,le,&da->da_coordinates);
175: PetscFree(lc);
176: PetscFree(ld);
177: PetscFree(le);
178: }
179: }
180: *cda = da->da_coordinates;
181: return(0);
182: }
187: /*@
188: DAGetGhostedCoordinates - Gets the node coordinates associated with a DA.
190: Collective on DA the first time it is called
192: Input Parameter:
193: . da - the distributed array
195: Output Parameter:
196: . c - coordinate vector
198: Note:
199: Each process has only the coordinates for its local AND ghost nodes
201: For two and three dimensions coordinates are interlaced (x_0,y_0,x_1,y_1,...)
202: and (x_0,y_0,z_0,x_1,y_1,z_1...)
204: You should not destroy or keep around this vector after the DA is destroyed.
206: Level: intermediate
208: .keywords: distributed array, get, corners, nodes, local indices, coordinates
210: .seealso: DAGetGhostCorners(), DASetCoordinates(), DASetUniformCoordinates(), DAGetCoordinates(), DAGetCoordinateDA()
211: @*/
212: PetscErrorCode DAGetGhostedCoordinates(DA da,Vec *c)
213: {
215:
218: if (!da->coordinates) SETERRQ(PETSC_ERR_ORDER,"You must call DASetCoordinates() before this call");
219: if (!da->ghosted_coordinates) {
220: DA dac;
222: DAGetCoordinateDA(da,&dac);
223: DACreateLocalVector(dac,&da->ghosted_coordinates);
224: if (dac == da) {PetscObjectDereference((PetscObject)dac);}
225: DAGlobalToLocalBegin(dac,da->coordinates,INSERT_VALUES,da->ghosted_coordinates);
226: DAGlobalToLocalEnd(dac,da->coordinates,INSERT_VALUES,da->ghosted_coordinates);
227: }
228: *c = da->ghosted_coordinates;
229: return(0);
230: }
234: /*@C
235: DASetFieldName - Sets the names of individual field components in multicomponent
236: vectors associated with a DA.
238: Not Collective
240: Input Parameters:
241: + da - the distributed array
242: . nf - field number for the DA (0, 1, ... dof-1), where dof indicates the
243: number of degrees of freedom per node within the DA
244: - names - the name of the field (component)
246: Level: intermediate
248: .keywords: distributed array, get, component name
250: .seealso: DAGetFieldName()
251: @*/
252: PetscErrorCode DASetFieldName(DA da,PetscInt nf,const char name[])
253: {
257:
259: if (nf < 0 || nf >= da->w) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Invalid field number: %D",nf);
260: if (da->fieldname[nf]) {PetscFree(da->fieldname[nf]);}
261:
262: PetscStrallocpy(name,&da->fieldname[nf]);
263: return(0);
264: }
268: /*@C
269: DAGetFieldName - Gets the names of individual field components in multicomponent
270: vectors associated with a DA.
272: Not Collective
274: Input Parameter:
275: + da - the distributed array
276: - nf - field number for the DA (0, 1, ... dof-1), where dof indicates the
277: number of degrees of freedom per node within the DA
279: Output Parameter:
280: . names - the name of the field (component)
282: Level: intermediate
284: .keywords: distributed array, get, component name
286: .seealso: DASetFieldName()
287: @*/
288: PetscErrorCode DAGetFieldName(DA da,PetscInt nf,char **name)
289: {
291:
294: if (nf < 0 || nf >= da->w) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Invalid field number: %D",nf);
295: *name = da->fieldname[nf];
296: return(0);
297: }
301: /*@
302: DAGetCorners - Returns the global (x,y,z) indices of the lower left
303: corner of the local region, excluding ghost points.
305: Not Collective
307: Input Parameter:
308: . da - the distributed array
310: Output Parameters:
311: + x,y,z - the corner indices (where y and z are optional; these are used
312: for 2D and 3D problems)
313: - m,n,p - widths in the corresponding directions (where n and p are optional;
314: these are used for 2D and 3D problems)
316: Note:
317: The corner information is independent of the number of degrees of
318: freedom per node set with the DACreateXX() routine. Thus the x, y, z, and
319: m, n, p can be thought of as coordinates on a logical grid, where each
320: grid point has (potentially) several degrees of freedom.
321: Any of y, z, n, and p can be passed in as PETSC_NULL if not needed.
323: Level: beginner
325: .keywords: distributed array, get, corners, nodes, local indices
327: .seealso: DAGetGhostCorners()
328: @*/
329: PetscErrorCode DAGetCorners(DA da,PetscInt *x,PetscInt *y,PetscInt *z,PetscInt *m,PetscInt *n,PetscInt *p)
330: {
331: PetscInt w;
335: /* since the xs, xe ... have all been multiplied by the number of degrees
336: of freedom per cell, w = da->w, we divide that out before returning.*/
337: w = da->w;
338: if (x) *x = da->xs/w; if(m) *m = (da->xe - da->xs)/w;
339: /* the y and z have NOT been multiplied by w */
340: if (y) *y = da->ys; if (n) *n = (da->ye - da->ys);
341: if (z) *z = da->zs; if (p) *p = (da->ze - da->zs);
342: return(0);
343: }