Actual source code: sda2.c
1: #define PETSCDM_DLL
3: /*
4: Simplified interface to PETSC DA (distributed array) object.
5: This is for a user who is not using PETSc Vecs (vectors).
6: */
8: #include petscda.h
10: EXTERN PetscErrorCode DALocalToLocalCreate(DA);
12: struct _n_SDA {
13: DA da;
14: Vec gvec,lvec,Gvec;
15: };
19: PetscErrorCode SDAArrayView(SDA da,PetscScalar *values,PetscViewer v)
20: {
24: VecPlaceArray(da->lvec,values);
25: if (!da->Gvec) {
26: DACreateGlobalVector(da->da,&da->Gvec);
27: }
28: DALocalToGlobalBegin(da->da,da->lvec,da->Gvec);
29: DALocalToGlobalEnd(da->da,da->lvec,da->Gvec);
30: VecView(da->Gvec,v);
31: return(0);
32: }
36: /*@C
37: SDACreate1d - Creates a one-dimensional regular array that is
38: distributed across some processors. This is the simplified
39: interface, must be used with SDAXXX operations, NOT DAXXX operations.
41: Input Parameters:
42: . comm - MPI communicator
43: . wrap - type of periodicity should the array have, if any
44: $ DA_NONPERIODIC, DA_XPERIODIC
45: . M - global dimension of the array
46: . w - number of degress of freedom per node
47: . s - stencil width
48: . lc - array containing number of nodes in X direction on each processor, or PETSC_NULL
50: Output Parameter:
51: . sda - the resulting array object
53: Level: beginner
55: .keywords: distributed array, create, two-dimensional
57: .seealso: SDADestroy(), SDACreate2d(), SDACreate3d()
58: @*/
59: PetscErrorCode SDACreate1d(MPI_Comm comm,DAPeriodicType wrap,PetscInt M,PetscInt w,PetscInt s,PetscInt *lc,SDA *sda)
60: {
62: DA da;
63: char **args;
64: int argc = 0;
65: Vec tmp;
67: PetscInitialize(&argc,&args,0,0);
70: PetscNew(struct _n_SDA,sda);
71: DACreate1d(comm,wrap,M,w,s,lc,&da);
72: (*sda)->da = da;
74: /* set up two dummy work vectors for the vector scatter */
75: DACreateLocalVector(da,&(*sda)->gvec);
76: /* we free the actual space in the vectors because it is not
77: needed since the user provides her/his own with SDA */
78: VecReplaceArray((*sda)->gvec,PETSC_NULL);
79: VecDuplicate((*sda)->gvec,&(*sda)->lvec);
80: VecReplaceArray((*sda)->lvec,PETSC_NULL);
82: /* destroy global vector */
83: DACreateGlobalVector(da,&tmp);
84: VecDestroy(tmp);
85: (*sda)->Gvec = 0;
87: /* free scatters in DA never needed by user */
88: DALocalToLocalCreate(da);
89: /* VecScatterDestroy(da->ltog);da->ltog = 0; */
90: /* VecScatterDestroy(da->gtol);da->gtol = 0; */
92: return(0);
93: }
97: /*@C
98: SDACreate2d - Creates a two-dimensional regular array that is
99: distributed across some processors. This is the simplified
100: interface, must be used with SDAXXX operations, NOT DAXXX operations.
102: Input Parameters:
103: . comm - MPI communicator
104: . wrap - type of periodicity should the array have, if any
105: $ DA_NONPERIODIC, DA_XPERIODIC,
106: $ DA_YPERIODIC, DA_XYPERIODIC
107: . stencil_type - stencil type either DA_STENCIL_BOX or DA_STENCIL_STAR
108: . M,N - global dimension in each direction of the array
109: . m,n - corresponding number of processors in each dimension
110: (or PETSC_DECIDE to have calculated)
111: . w - number of degress of freedom per node
112: . s - stencil width
113: . lx, ly - arrays containing the number of nodes in each cell along
114: $ the x and y coordinates, or PETSC_NULL
116: Output Parameter:
117: . inra - the resulting array object
119: Level: beginner
121: .keywords: distributed array, create, two-dimensional
123: .seealso: DADestroy(), DAView(), DACreate1d(), DACreate3d()
124: @*/
125: PetscErrorCode SDACreate2d(MPI_Comm comm,DAPeriodicType wrap,DAStencilType stencil_type,
126: PetscInt M,PetscInt N,PetscInt m,PetscInt n,PetscInt w,PetscInt s,PetscInt *lx,PetscInt *ly,SDA *sda)
127: {
129: DA da;
130: char **args;
131: int argc = 0;
132: Vec tmp;
134: PetscInitialize(&argc,&args,0,0);
137: PetscNew(struct _n_SDA,sda);
138: DACreate2d(comm,wrap,stencil_type,M,N,m,n,w,s,lx,ly,&da);
139: (*sda)->da = da;
141: /* set up two dummy work vectors for the vector scatter */
142: DACreateLocalVector(da,&(*sda)->gvec);
143: /* we free the actual space in the vectors because it is not
144: needed since the user provides her/his own with SDA */
145: VecReplaceArray((*sda)->gvec,PETSC_NULL);
146: VecDuplicate((*sda)->gvec,&(*sda)->lvec);
147: VecReplaceArray((*sda)->lvec,PETSC_NULL);
149: /* destroy global vector */
150: DACreateGlobalVector(da,&tmp);
151: VecDestroy(tmp);
152: (*sda)->Gvec = 0;
154: /* free scatters in DA never needed by user */
155: DALocalToLocalCreate(da);
156: /*VecScatterDestroy(da->ltog);da->ltog = 0; */
157: /*VecScatterDestroy(da->gtol);da->gtol = 0;*/
159: return(0);
160: }
164: /*@C
165: SDACreate3d - Creates a three-dimensional regular array that is
166: distributed across some processors. This is the simplified
167: interface, must be used with SDAXXX operations, NOT DAXXX operations.
169: Input Parameters:
170: . comm - MPI communicator
171: . wrap - type of periodicity should the array have, if any
172: $ DA_NONPERIODIC, DA_XPERIODIC,
173: $ DA_YPERIODIC, DA_XYPERIODIC
174: . stencil_type - stencil type either DA_STENCIL_BOX or DA_STENCIL_STAR
175: . M,N - global dimension in each direction of the array
176: . m,n - corresponding number of processors in each dimension
177: (or PETSC_DECIDE to have calculated)
178: . w - number of degress of freedom per node
179: . s - stencil width
180: . lx, ly, lz - arrays containing the number of nodes in each cell along
181: $ the x, y, and z coordinates, or PETSC_NUL
183: Output Parameter:
184: . inra - the resulting array object
186: Level: beginner
188: .keywords: distributed array, create, two-dimensional
190: .seealso: DADestroy(), DAView(), DACreate1d(), DACreate3d()
191: @*/
192: PetscErrorCode SDACreate3d(MPI_Comm comm,DAPeriodicType wrap,DAStencilType stencil_type,PetscInt M,
193: PetscInt N,PetscInt P,PetscInt m,PetscInt n,PetscInt p,PetscInt w,PetscInt s,PetscInt *lx,PetscInt *ly,PetscInt *lz,SDA *sda)
194: {
196: DA da;
197: Vec tmp;
198: char **args;
199: int argc = 0;
201: PetscInitialize(&argc,&args,0,0);
204: PetscNew(struct _n_SDA,sda);
205: DACreate3d(comm,wrap,stencil_type,M,N,P,m,n,p,w,s,lx,ly,lz,&da);
206: (*sda)->da = da;
208: /* set up two dummy work vectors for the vector scatter */
209: DACreateLocalVector(da,&(*sda)->gvec);
210: /* we free the actual space in the vectors because it is not
211: needed since the user provides her/his own with SDA */
212: VecReplaceArray((*sda)->gvec,PETSC_NULL);
213: VecDuplicate((*sda)->gvec,&(*sda)->lvec);
214: VecReplaceArray((*sda)->lvec,PETSC_NULL);
216: /* destroy global vector */
217: DACreateGlobalVector(da,&tmp);
218: VecDestroy(tmp);
219: (*sda)->Gvec = 0;
220: /* free scatters in DA never needed by user */
221: DALocalToLocalCreate(da);
222: /*VecScatterDestroy(da->ltog);da->ltog = 0;*/
223: /*VecScatterDestroy(da->gtol);da->gtol = 0;*/
225: return(0);
226: }
230: /*@
231: SDADestroy - Destroys simple distributed array.
233: Input parameters:
234: sda - distributed array
236: Level: beginner
238: .keywords: distributed array
240: .seealso: SDACreate1d(), SDACreate2d(), SDACreate3d()
241: @*/
242: PetscErrorCode SDADestroy(SDA sda)
243: {
247: VecDestroy(sda->gvec);
248: VecDestroy(sda->lvec);
249: if (sda->Gvec) {VecDestroy(sda->Gvec);}
250: DADestroy(sda->da);
251: PetscFree(sda);
252: return(0);
253: }
257: /*@
258: SDALocalToLocalBegin - Maps from a local representation (including
259: ghostpoints) to another where the ghostpoints in the second are
260: set correctly. Must be followed by SDALocalToLocalEnd().
262: Input Parameters:
263: . da - the distributed array context
264: . g - the original vector
265: . mode - one of INSERT_VALUES or ADD_VALUES
267: Output Parameter:
268: . l - the vector with correct ghost values
270: Level: beginner
272: .keywords: distributed array, global to local, begin
274: .seealso: SDALocalToLocalEnd(), SDACreate2d()
275: @*/
276: PetscErrorCode SDALocalToLocalBegin(SDA sda,PetscScalar *g,InsertMode mode,PetscScalar *l)
277: {
279: DA da = sda->da;
280: Vec gvec = sda->gvec,lvec = sda->lvec;
283: VecPlaceArray(gvec,g);
284: VecPlaceArray(lvec,l);
285: DALocalToLocalBegin(da,gvec,mode,lvec);
286: return(0);
287: }
291: /*@
292: SDALocalToLocalEnd - Maps from a local representation (including
293: ghostpoints) to another where the ghostpoints in the second are
294: set correctly. Must be preceeded by SDALocalToLocalBegin().
296: Input Parameters:
297: . da - the distributed array context
298: . g - the original vector
299: . mode - one of INSERT_VALUES or ADD_VALUES
301: Output Parameter:
302: . l - the vector with correct ghost values
304: Level: beginner
306: .keywords: distributed array, global to local, end
308: .seealso: SDALocalToLocalBegin(), SDACreate2d()
309: @*/
310: PetscErrorCode SDALocalToLocalEnd(SDA sda,PetscScalar *g,InsertMode mode,PetscScalar *l)
311: {
313: DA da = sda->da;
314: Vec gvec = sda->gvec,lvec = sda->lvec;
317: VecPlaceArray(gvec,g);
318: VecPlaceArray(lvec,l);
319: DALocalToLocalEnd(da,gvec,mode,lvec);
320: return(0);
321: }
322:
325: /*@C
326: SDAGetCorners - Returns the global (x,y,z) indices of the lower left
327: corner of the local region, excluding ghost points.
329: Input Parameter:
330: . da - the distributed array
332: Output Parameters:
333: . x,y,z - the corner indices
334: $ y and z are optional (used for 2D and 3D problems)
335: . m,n,p - widths in the corresponding directions
336: $ n and p are optional (used for 2D and 3D problems)
338: Note:
339: Any of y, z, n, and p should be set to PETSC_NULL if not needed.
341: Level: beginner
343: .keywords: distributed array, get, corners, nodes, local indices
345: .seealso: SDAGetGhostCorners()
346: @*/
347: PetscErrorCode SDAGetCorners(SDA da,PetscInt *x,PetscInt *y,PetscInt *z,PetscInt *m,PetscInt *n,PetscInt *p)
348: {
352: DAGetCorners(da->da,x,y,z,m,n,p);
353: return(0);
354: }
358: /*@C
359: SDAGetGhostCorners - Returns the global (x,y,z) indices of the lower left
360: corner of the local region, including ghost points.
362: Input Parameter:
363: . da - the distributed array
365: Output Parameters:
366: . x,y,z - the corner indices
367: $ y and z are optional (used for 2D and 3D problems)
368: . m,n,p - widths in the corresponding directions
369: $ n and p are optional (used for 2D and 3D problems)
371: Note:
372: Any of y, z, n, and p should be set to PETSC_NULL if not needed.
375: Level: beginner
377: .keywords: distributed array, get, ghost, corners, nodes, local indices
379: .seealso: SDAGetCorners()
380: @*/
381: PetscErrorCode SDAGetGhostCorners(SDA da,PetscInt *x,PetscInt *y,PetscInt *z,PetscInt *m,PetscInt *n,PetscInt *p)
382: {
386: DAGetGhostCorners(da->da,x,y,z,m,n,p);
387: return(0);
388: }