Actual source code: sda2.c
1: /*
2: Simplified interface to PETSC DA (distributed array) object.
3: This is for a user who is not using PETSc Vecs (vectors).
4: */
6: #include petscda.h
8: EXTERN PetscErrorCode DALocalToLocalCreate(DA);
10: struct _SDA {
11: DA da;
12: Vec gvec,lvec,Gvec;
13: };
17: PetscErrorCode SDAArrayView(SDA da,PetscScalar *values,PetscViewer v)
18: {
22: VecPlaceArray(da->lvec,values);
23: if (!da->Gvec) {
24: DACreateGlobalVector(da->da,&da->Gvec);
25: }
26: DALocalToGlobalBegin(da->da,da->lvec,da->Gvec);
27: DALocalToGlobalEnd(da->da,da->lvec,da->Gvec);
28: VecView(da->Gvec,v);
29: return(0);
30: }
34: /*@C
35: SDACreate1d - Creates a one-dimensional regular array that is
36: distributed across some processors. This is the simplified
37: interface, must be used with SDAXXX operations, NOT DAXXX operations.
39: Input Parameters:
40: . comm - MPI communicator
41: . wrap - type of periodicity should the array have, if any
42: $ DA_NONPERIODIC, DA_XPERIODIC
43: . M - global dimension of the array
44: . w - number of degress of freedom per node
45: . s - stencil width
46: . lc - array containing number of nodes in X direction on each processor, or PETSC_NULL
48: Output Parameter:
49: . sda - the resulting array object
51: Level: beginner
53: .keywords: distributed array, create, two-dimensional
55: .seealso: SDADestroy(), SDACreate2d(), SDACreate3d()
56: @*/
57: PetscErrorCode SDACreate1d(MPI_Comm comm,DAPeriodicType wrap,PetscInt M,PetscInt w,PetscInt s,PetscInt *lc,SDA *sda)
58: {
60: DA da;
61: char **args;
62: int argc = 0;
63: Vec tmp;
65: PetscInitialize(&argc,&args,0,0);
68: PetscNew(struct _SDA,sda);
69: DACreate1d(comm,wrap,M,w,s,lc,&da);
70: (*sda)->da = da;
72: /* set up two dummy work vectors for the vector scatter */
73: DACreateLocalVector(da,&(*sda)->gvec);
74: /* we free the actual space in the vectors because it is not
75: needed since the user provides her/his own with SDA */
76: VecReplaceArray((*sda)->gvec,PETSC_NULL);
77: VecDuplicate((*sda)->gvec,&(*sda)->lvec);
78: VecReplaceArray((*sda)->lvec,PETSC_NULL);
80: /* destroy global vector */
81: DACreateGlobalVector(da,&tmp);
82: VecDestroy(tmp);
83: (*sda)->Gvec = 0;
85: /* free scatters in DA never needed by user */
86: DALocalToLocalCreate(da);
87: /* VecScatterDestroy(da->ltog);da->ltog = 0; */
88: /* VecScatterDestroy(da->gtol);da->gtol = 0; */
90: return(0);
91: }
95: /*@C
96: SDACreate2d - Creates a two-dimensional regular array that is
97: distributed across some processors. This is the simplified
98: interface, must be used with SDAXXX operations, NOT DAXXX operations.
100: Input Parameters:
101: . comm - MPI communicator
102: . wrap - type of periodicity should the array have, if any
103: $ DA_NONPERIODIC, DA_XPERIODIC,
104: $ DA_YPERIODIC, DA_XYPERIODIC
105: . stencil_type - stencil type either DA_STENCIL_BOX or DA_STENCIL_STAR
106: . M,N - global dimension in each direction of the array
107: . m,n - corresponding number of processors in each dimension
108: (or PETSC_DECIDE to have calculated)
109: . w - number of degress of freedom per node
110: . s - stencil width
111: . lx, ly - arrays containing the number of nodes in each cell along
112: $ the x and y coordinates, or PETSC_NULL
114: Output Parameter:
115: . inra - the resulting array object
117: Level: beginner
119: .keywords: distributed array, create, two-dimensional
121: .seealso: DADestroy(), DAView(), DACreate1d(), DACreate3d()
122: @*/
123: PetscErrorCode SDACreate2d(MPI_Comm comm,DAPeriodicType wrap,DAStencilType stencil_type,
124: PetscInt M,PetscInt N,PetscInt m,PetscInt n,PetscInt w,PetscInt s,PetscInt *lx,PetscInt *ly,SDA *sda)
125: {
127: DA da;
128: char **args;
129: int argc = 0;
130: Vec tmp;
132: PetscInitialize(&argc,&args,0,0);
135: PetscNew(struct _SDA,sda);
136: DACreate2d(comm,wrap,stencil_type,M,N,m,n,w,s,lx,ly,&da);
137: (*sda)->da = da;
139: /* set up two dummy work vectors for the vector scatter */
140: DACreateLocalVector(da,&(*sda)->gvec);
141: /* we free the actual space in the vectors because it is not
142: needed since the user provides her/his own with SDA */
143: VecReplaceArray((*sda)->gvec,PETSC_NULL);
144: VecDuplicate((*sda)->gvec,&(*sda)->lvec);
145: VecReplaceArray((*sda)->lvec,PETSC_NULL);
147: /* destroy global vector */
148: DACreateGlobalVector(da,&tmp);
149: VecDestroy(tmp);
150: (*sda)->Gvec = 0;
152: /* free scatters in DA never needed by user */
153: DALocalToLocalCreate(da);
154: /*VecScatterDestroy(da->ltog);da->ltog = 0; */
155: /*VecScatterDestroy(da->gtol);da->gtol = 0;*/
157: return(0);
158: }
162: /*@C
163: SDACreate3d - Creates a three-dimensional regular array that is
164: distributed across some processors. This is the simplified
165: interface, must be used with SDAXXX operations, NOT DAXXX operations.
167: Input Parameters:
168: . comm - MPI communicator
169: . wrap - type of periodicity should the array have, if any
170: $ DA_NONPERIODIC, DA_XPERIODIC,
171: $ DA_YPERIODIC, DA_XYPERIODIC
172: . stencil_type - stencil type either DA_STENCIL_BOX or DA_STENCIL_STAR
173: . M,N - global dimension in each direction of the array
174: . m,n - corresponding number of processors in each dimension
175: (or PETSC_DECIDE to have calculated)
176: . w - number of degress of freedom per node
177: . s - stencil width
178: . lx, ly, lz - arrays containing the number of nodes in each cell along
179: $ the x, y, and z coordinates, or PETSC_NUL
181: Output Parameter:
182: . inra - the resulting array object
184: Level: beginner
186: .keywords: distributed array, create, two-dimensional
188: .seealso: DADestroy(), DAView(), DACreate1d(), DACreate3d()
189: @*/
190: PetscErrorCode SDACreate3d(MPI_Comm comm,DAPeriodicType wrap,DAStencilType stencil_type,PetscInt M,
191: PetscInt N,PetscInt P,PetscInt m,PetscInt n,PetscInt p,PetscInt w,PetscInt s,PetscInt *lx,PetscInt *ly,PetscInt *lz,SDA *sda)
192: {
194: DA da;
195: Vec tmp;
196: char **args;
197: int argc = 0;
199: PetscInitialize(&argc,&args,0,0);
202: PetscNew(struct _SDA,sda);
203: DACreate3d(comm,wrap,stencil_type,M,N,P,m,n,p,w,s,lx,ly,lz,&da);
204: (*sda)->da = da;
206: /* set up two dummy work vectors for the vector scatter */
207: DACreateLocalVector(da,&(*sda)->gvec);
208: /* we free the actual space in the vectors because it is not
209: needed since the user provides her/his own with SDA */
210: VecReplaceArray((*sda)->gvec,PETSC_NULL);
211: VecDuplicate((*sda)->gvec,&(*sda)->lvec);
212: VecReplaceArray((*sda)->lvec,PETSC_NULL);
214: /* destroy global vector */
215: DACreateGlobalVector(da,&tmp);
216: VecDestroy(tmp);
217: (*sda)->Gvec = 0;
218: /* free scatters in DA never needed by user */
219: DALocalToLocalCreate(da);
220: /*VecScatterDestroy(da->ltog);da->ltog = 0;*/
221: /*VecScatterDestroy(da->gtol);da->gtol = 0;*/
223: return(0);
224: }
228: /*@C
229: SDADestroy - Destroys simple distributed array.
231: Input parameters:
232: sda - distributed array
234: Level: beginner
236: .keywords: distributed array
238: .seealso: SDACreate1d(), SDACreate2d(), SDACreate3d()
239: @*/
240: PetscErrorCode SDADestroy(SDA sda)
241: {
245: VecDestroy(sda->gvec);
246: VecDestroy(sda->lvec);
247: if (sda->Gvec) {VecDestroy(sda->Gvec);}
248: DADestroy(sda->da);
249: PetscFree(sda);
250: return(0);
251: }
255: /*@C
256: SDALocalToLocalBegin - Maps from a local representation (including
257: ghostpoints) to another where the ghostpoints in the second are
258: set correctly. Must be followed by SDALocalToLocalEnd().
260: Input Parameters:
261: . da - the distributed array context
262: . g - the original vector
263: . mode - one of INSERT_VALUES or ADD_VALUES
265: Output Parameter:
266: . l - the vector with correct ghost values
268: Level: beginner
270: .keywords: distributed array, global to local, begin
272: .seealso: SDALocalToLocalEnd(), SDACreate2d()
273: @*/
274: PetscErrorCode SDALocalToLocalBegin(SDA sda,PetscScalar *g,InsertMode mode,PetscScalar *l)
275: {
277: DA da = sda->da;
278: Vec gvec = sda->gvec,lvec = sda->lvec;
281: VecPlaceArray(gvec,g);
282: VecPlaceArray(lvec,l);
283: DALocalToLocalBegin(da,gvec,mode,lvec);
284: return(0);
285: }
289: /*@C
290: SDALocalToLocalEnd - Maps from a local representation (including
291: ghostpoints) to another where the ghostpoints in the second are
292: set correctly. Must be preceeded by SDALocalToLocalBegin().
294: Input Parameters:
295: . da - the distributed array context
296: . g - the original vector
297: . mode - one of INSERT_VALUES or ADD_VALUES
299: Output Parameter:
300: . l - the vector with correct ghost values
302: Level: beginner
304: .keywords: distributed array, global to local, end
306: .seealso: SDALocalToLocalBegin(), SDACreate2d()
307: @*/
308: PetscErrorCode SDALocalToLocalEnd(SDA sda,PetscScalar *g,InsertMode mode,PetscScalar *l)
309: {
311: DA da = sda->da;
312: Vec gvec = sda->gvec,lvec = sda->lvec;
315: VecPlaceArray(gvec,g);
316: VecPlaceArray(lvec,l);
317: DALocalToLocalEnd(da,gvec,mode,lvec);
318: return(0);
319: }
320:
323: /*@C
324: SDAGetCorners - Returns the global (x,y,z) indices of the lower left
325: corner of the local region, excluding ghost points.
327: Input Parameter:
328: . da - the distributed array
330: Output Parameters:
331: . x,y,z - the corner indices
332: $ y and z are optional (used for 2D and 3D problems)
333: . m,n,p - widths in the corresponding directions
334: $ n and p are optional (used for 2D and 3D problems)
336: Note:
337: Any of y, z, n, and p should be set to PETSC_NULL if not needed.
339: Level: beginner
341: .keywords: distributed array, get, corners, nodes, local indices
343: .seealso: SDAGetGhostCorners()
344: @*/
345: PetscErrorCode SDAGetCorners(SDA da,PetscInt *x,PetscInt *y,PetscInt *z,PetscInt *m,PetscInt *n,PetscInt *p)
346: {
350: DAGetCorners(da->da,x,y,z,m,n,p);
351: return(0);
352: }
356: /*@C
357: SDAGetGhostCorners - Returns the global (x,y,z) indices of the lower left
358: corner of the local region, including ghost points.
360: Input Parameter:
361: . da - the distributed array
363: Output Parameters:
364: . x,y,z - the corner indices
365: $ y and z are optional (used for 2D and 3D problems)
366: . m,n,p - widths in the corresponding directions
367: $ n and p are optional (used for 2D and 3D problems)
369: Note:
370: Any of y, z, n, and p should be set to PETSC_NULL if not needed.
373: Level: beginner
375: .keywords: distributed array, get, ghost, corners, nodes, local indices
377: .seealso: SDAGetCorners()
378: @*/
379: PetscErrorCode SDAGetGhostCorners(SDA da,PetscInt *x,PetscInt *y,PetscInt *z,PetscInt *m,PetscInt *n,PetscInt *p)
380: {
384: DAGetGhostCorners(da->da,x,y,z,m,n,p);
385: return(0);
386: }