Actual source code: sda2.c
1: /*$Id: sda2.c,v 1.24 2001/08/07 03:04:49 balay Exp $*/
2: /*
3: Simplified interface to PETSC DA (distributed array) object.
4: This is for a user who is not using PETSc Vecs (vectors).
5: */
7: #include petscda.h
9: extern int DALocalToLocalCreate(DA);
11: struct _SDA {
12: DA da;
13: Vec gvec,lvec,Gvec;
14: };
18: int SDAArrayView(SDA da,PetscScalar *values,PetscViewer v)
19: {
23: VecPlaceArray(da->lvec,values);
24: if (!da->Gvec) {
25: DACreateGlobalVector(da->da,&da->Gvec);
26: }
27: DALocalToGlobalBegin(da->da,da->lvec,da->Gvec);
28: DALocalToGlobalEnd(da->da,da->lvec,da->Gvec);
29: VecView(da->Gvec,v);
30: return(0);
31: }
35: /*@C
36: SDACreate1d - Creates a one-dimensional regular array that is
37: distributed across some processors. This is the simplified
38: interface, must be used with SDAXXX operations, NOT DAXXX operations.
40: Input Parameters:
41: . comm - MPI communicator
42: . wrap - type of periodicity should the array have, if any
43: $ DA_NONPERIODIC, DA_XPERIODIC
44: . M - global dimension of the array
45: . w - number of degress of freedom per node
46: . s - stencil width
47: . lc - array containing number of nodes in X direction on each processor, or PETSC_NULL
49: Output Parameter:
50: . sda - the resulting array object
52: Level: beginner
54: .keywords: distributed array, create, two-dimensional
56: .seealso: SDADestroy(), SDACreate2d(), SDACreate3d()
57: @*/
58: int SDACreate1d(MPI_Comm comm,DAPeriodicType wrap,int M,int w,int s,int *lc,SDA *sda)
59: {
60: int ierr;
61: DA da;
62: char **args;
63: int argc = 0;
64: Vec tmp;
66: PetscInitialize(&argc,&args,0,0);
69: PetscNew(struct _SDA,sda);
70: DACreate1d(comm,wrap,M,w,s,lc,&da);
71: (*sda)->da = da;
73: /* set up two dummy work vectors for the vector scatter */
74: DACreateLocalVector(da,&(*sda)->gvec);
75: /* we free the actual space in the vectors because it is not
76: needed since the user provides her/his own with SDA */
77: VecReplaceArray((*sda)->gvec,PETSC_NULL);
78: VecDuplicate((*sda)->gvec,&(*sda)->lvec);
79: VecReplaceArray((*sda)->lvec,PETSC_NULL);
81: /* destroy global vector */
82: DACreateGlobalVector(da,&tmp);
83: VecDestroy(tmp);
84: (*sda)->Gvec = 0;
86: /* free scatters in DA never needed by user */
87: DALocalToLocalCreate(da);
88: /* VecScatterDestroy(da->ltog);da->ltog = 0; */
89: /* VecScatterDestroy(da->gtol);da->gtol = 0; */
91: return(0);
92: }
96: /*@C
97: SDACreate2d - Creates a two-dimensional regular array that is
98: distributed across some processors. This is the simplified
99: interface, must be used with SDAXXX operations, NOT DAXXX operations.
101: Input Parameters:
102: . comm - MPI communicator
103: . wrap - type of periodicity should the array have, if any
104: $ DA_NONPERIODIC, DA_XPERIODIC,
105: $ DA_YPERIODIC, DA_XYPERIODIC
106: . stencil_type - stencil type either DA_STENCIL_BOX or DA_STENCIL_STAR
107: . M,N - global dimension in each direction of the array
108: . m,n - corresponding number of processors in each dimension
109: (or PETSC_DECIDE to have calculated)
110: . w - number of degress of freedom per node
111: . s - stencil width
112: . lx, ly - arrays containing the number of nodes in each cell along
113: $ the x and y coordinates, or PETSC_NULL
115: Output Parameter:
116: . inra - the resulting array object
118: Level: beginner
120: .keywords: distributed array, create, two-dimensional
122: .seealso: DADestroy(), DAView(), DACreate1d(), DACreate3d()
123: @*/
124: int SDACreate2d(MPI_Comm comm,DAPeriodicType wrap,DAStencilType stencil_type,
125: int M,int N,int m,int n,int w,int s,int *lx,int *ly,SDA *sda)
126: {
127: int ierr;
128: DA da;
129: char **args;
130: int argc = 0;
131: Vec tmp;
133: PetscInitialize(&argc,&args,0,0);
136: PetscNew(struct _SDA,sda);
137: DACreate2d(comm,wrap,stencil_type,M,N,m,n,w,s,lx,ly,&da);
138: (*sda)->da = da;
140: /* set up two dummy work vectors for the vector scatter */
141: DACreateLocalVector(da,&(*sda)->gvec);
142: /* we free the actual space in the vectors because it is not
143: needed since the user provides her/his own with SDA */
144: VecReplaceArray((*sda)->gvec,PETSC_NULL);
145: VecDuplicate((*sda)->gvec,&(*sda)->lvec);
146: VecReplaceArray((*sda)->lvec,PETSC_NULL);
148: /* destroy global vector */
149: DACreateGlobalVector(da,&tmp);
150: VecDestroy(tmp);
151: (*sda)->Gvec = 0;
153: /* free scatters in DA never needed by user */
154: DALocalToLocalCreate(da);
155: /*VecScatterDestroy(da->ltog);da->ltog = 0; */
156: /*VecScatterDestroy(da->gtol);da->gtol = 0;*/
158: return(0);
159: }
163: /*@C
164: SDACreate3d - Creates a three-dimensional regular array that is
165: distributed across some processors. This is the simplified
166: interface, must be used with SDAXXX operations, NOT DAXXX operations.
168: Input Parameters:
169: . comm - MPI communicator
170: . wrap - type of periodicity should the array have, if any
171: $ DA_NONPERIODIC, DA_XPERIODIC,
172: $ DA_YPERIODIC, DA_XYPERIODIC
173: . stencil_type - stencil type either DA_STENCIL_BOX or DA_STENCIL_STAR
174: . M,N - global dimension in each direction of the array
175: . m,n - corresponding number of processors in each dimension
176: (or PETSC_DECIDE to have calculated)
177: . w - number of degress of freedom per node
178: . s - stencil width
179: . lx, ly, lz - arrays containing the number of nodes in each cell along
180: $ the x, y, and z coordinates, or PETSC_NUL
182: Output Parameter:
183: . inra - the resulting array object
185: Level: beginner
187: .keywords: distributed array, create, two-dimensional
189: .seealso: DADestroy(), DAView(), DACreate1d(), DACreate3d()
190: @*/
191: int SDACreate3d(MPI_Comm comm,DAPeriodicType wrap,DAStencilType stencil_type,int M,
192: int N,int P,int m,int n,int p,int w,int s,int *lx,int *ly,int *lz,SDA *sda)
193: {
194: int ierr;
195: DA da;
196: Vec tmp;
197: char **args;
198: int argc = 0;
200: PetscInitialize(&argc,&args,0,0);
203: PetscNew(struct _SDA,sda);
204: DACreate3d(comm,wrap,stencil_type,M,N,P,m,n,p,w,s,lx,ly,lz,&da);
205: (*sda)->da = da;
207: /* set up two dummy work vectors for the vector scatter */
208: DACreateLocalVector(da,&(*sda)->gvec);
209: /* we free the actual space in the vectors because it is not
210: needed since the user provides her/his own with SDA */
211: VecReplaceArray((*sda)->gvec,PETSC_NULL);
212: VecDuplicate((*sda)->gvec,&(*sda)->lvec);
213: VecReplaceArray((*sda)->lvec,PETSC_NULL);
215: /* destroy global vector */
216: DACreateGlobalVector(da,&tmp);
217: VecDestroy(tmp);
218: (*sda)->Gvec = 0;
219: /* free scatters in DA never needed by user */
220: DALocalToLocalCreate(da);
221: /*VecScatterDestroy(da->ltog);da->ltog = 0;*/
222: /*VecScatterDestroy(da->gtol);da->gtol = 0;*/
224: return(0);
225: }
229: /*@C
230: SDADestroy - Destroys simple distributed array.
232: Input parameters:
233: sda - distributed array
235: Level: beginner
237: .keywords: distributed array
239: .seealso: SDACreate1d(), SDACreate2d(), SDACreate3d()
240: @*/
241: int SDADestroy(SDA sda)
242: {
246: VecDestroy(sda->gvec);
247: VecDestroy(sda->lvec);
248: if (sda->Gvec) {VecDestroy(sda->Gvec);}
249: DADestroy(sda->da);
250: PetscFree(sda);
251: return(0);
252: }
256: /*@C
257: SDALocalToLocalBegin - Maps from a local representation (including
258: ghostpoints) to another where the ghostpoints in the second are
259: set correctly. Must be followed by SDALocalToLocalEnd().
261: Input Parameters:
262: . da - the distributed array context
263: . g - the original vector
264: . mode - one of INSERT_VALUES or ADD_VALUES
266: Output Parameter:
267: . l - the vector with correct ghost values
269: Level: beginner
271: .keywords: distributed array, global to local, begin
273: .seealso: SDALocalToLocalEnd(), SDACreate2d()
274: @*/
275: int SDALocalToLocalBegin(SDA sda,PetscScalar *g,InsertMode mode,PetscScalar *l)
276: {
278: DA da = sda->da;
279: Vec gvec = sda->gvec,lvec = sda->lvec;
282: VecPlaceArray(gvec,g);
283: VecPlaceArray(lvec,l);
284: DALocalToLocalBegin(da,gvec,mode,lvec);
285: return(0);
286: }
290: /*@C
291: SDALocalToLocalEnd - Maps from a local representation (including
292: ghostpoints) to another where the ghostpoints in the second are
293: set correctly. Must be preceeded by SDALocalToLocalBegin().
295: Input Parameters:
296: . da - the distributed array context
297: . g - the original vector
298: . mode - one of INSERT_VALUES or ADD_VALUES
300: Output Parameter:
301: . l - the vector with correct ghost values
303: Level: beginner
305: .keywords: distributed array, global to local, end
307: .seealso: SDALocalToLocalBegin(), SDACreate2d()
308: @*/
309: int SDALocalToLocalEnd(SDA sda,PetscScalar *g,InsertMode mode,PetscScalar *l)
310: {
312: DA da = sda->da;
313: Vec gvec = sda->gvec,lvec = sda->lvec;
316: VecPlaceArray(gvec,g);
317: VecPlaceArray(lvec,l);
318: DALocalToLocalEnd(da,gvec,mode,lvec);
319: return(0);
320: }
321:
324: /*@C
325: SDAGetCorners - Returns the global (x,y,z) indices of the lower left
326: corner of the local region, excluding ghost points.
328: Input Parameter:
329: . da - the distributed array
331: Output Parameters:
332: . x,y,z - the corner indices
333: $ y and z are optional (used for 2D and 3D problems)
334: . m,n,p - widths in the corresponding directions
335: $ n and p are optional (used for 2D and 3D problems)
337: Note:
338: Any of y, z, n, and p should be set to PETSC_NULL if not needed.
340: Level: beginner
342: .keywords: distributed array, get, corners, nodes, local indices
344: .seealso: SDAGetGhostCorners()
345: @*/
346: int SDAGetCorners(SDA da,int *x,int *y,int *z,int *m,int *n,int *p)
347: {
351: DAGetCorners(da->da,x,y,z,m,n,p);
352: return(0);
353: }
357: /*@C
358: SDAGetGhostCorners - Returns the global (x,y,z) indices of the lower left
359: corner of the local region, including ghost points.
361: Input Parameter:
362: . da - the distributed array
364: Output Parameters:
365: . x,y,z - the corner indices
366: $ y and z are optional (used for 2D and 3D problems)
367: . m,n,p - widths in the corresponding directions
368: $ n and p are optional (used for 2D and 3D problems)
370: Note:
371: Any of y, z, n, and p should be set to PETSC_NULL if not needed.
374: Level: beginner
376: .keywords: distributed array, get, ghost, corners, nodes, local indices
378: .seealso: SDAGetCorners()
379: @*/
380: int SDAGetGhostCorners(SDA da,int *x,int *y,int *z,int *m,int *n,int *p)
381: {
385: DAGetGhostCorners(da->da,x,y,z,m,n,p);
386: return(0);
387: }