Actual source code: dagtol.c
1: /*$Id: dagtol.c,v 1.29 2001/03/23 23:25:00 balay Exp $*/
2:
3: /*
4: Code for manipulating distributed regular arrays in parallel.
5: */
7: #include "src/dm/da/daimpl.h" /*I "petscda.h" I*/
9: /*@
10: DAGlobalToLocalBegin - Maps values from the global vector to the local
11: patch; the ghost points are included. Must be followed by
12: DAGlobalToLocalEnd() to complete the exchange.
14: Collective on DA
16: Input Parameters:
17: + da - the distributed array context
18: . g - the global vector
19: - mode - one of INSERT_VALUES or ADD_VALUES
21: Output Parameter:
22: . l - the local values
24: Level: beginner
26: Notes:
27: The global and local vectors used here need not be the same as those
28: obtained from DACreateGlobalVector() and DACreateLocalVector(), BUT they
29: must have the same parallel data layout; they could, for example, be
30: obtained with VecDuplicate() from the DA originating vectors.
32: .keywords: distributed array, global to local, begin
34: .seealso: DAGlobalToLocalEnd(), DALocalToGlobal(), DACreate2d(),
35: DALocalToLocalBegin(), DALocalToLocalEnd()
37: @*/
38: int DAGlobalToLocalBegin(DA da,Vec g,InsertMode mode,Vec l)
39: {
46: VecScatterBegin(g,l,mode,SCATTER_FORWARD,da->gtol);
47: return(0);
48: }
50: /*@
51: DAGlobalToLocalEnd - Maps values from the global vector to the local
52: patch; the ghost points are included. Must be preceeded by
53: DAGlobalToLocalBegin().
55: Collective on DA
57: Input Parameters:
58: + da - the distributed array context
59: . g - the global vector
60: - mode - one of INSERT_VALUES or ADD_VALUES
62: Output Parameter:
63: . l - the local values
65: Level: beginner
67: Notes:
68: The global and local vectors used here need not be the same as those
69: obtained from DACreateGlobalVector() and DACreateLocalVector(), BUT they
70: must have the same parallel data layout; they could, for example, be
71: obtained with VecDuplicate() from the DA originating vectors.
73: .keywords: distributed array, global to local, end
75: .seealso: DAGlobalToLocalBegin(), DALocalToGlobal(), DACreate2d(),
76: DALocalToLocalBegin(), DALocalToLocalEnd()
77: @*/
78: int DAGlobalToLocalEnd(DA da,Vec g,InsertMode mode,Vec l)
79: {
86: VecScatterEnd(g,l,mode,SCATTER_FORWARD,da->gtol);
87: return(0);
88: }
90: /*
91: DAGlobalToNatural_Create - Create the global to natural scatter object
93: Collective on DA
95: Input Parameter:
96: . da - the distributed array context
98: Level: developer
100: Notes: This is an internal routine called by DAGlobalToNatural() to
101: create the scatter context.
103: .keywords: distributed array, global to local, begin
105: .seealso: DAGlobalToNaturalEnd(), DALocalToGlobal(), DACreate2d(),
106: DAGlobalToLocalBegin(), DAGlobalToLocalEnd(), DACreateNaturalVector()
107: */
108: int DAGlobalToNatural_Create(DA da)
109: {
110: int ierr,m,start;
111: IS from,to;
115: if (!da->natural) {
116: SETERRQ(1,"Natural layout vector not yet created; cannot scatter into it");
117: }
118: if (!da->ao) {
119: SETERRQ(1,"Cannot use -da_noao with this function");
120: }
122: /* create the scatter context */
123: VecGetLocalSize(da->natural,&m);
124: VecGetOwnershipRange(da->natural,&start,PETSC_NULL);
125: ISCreateStride(da->comm,m,start,1,&to);
126: AOPetscToApplicationIS(da->ao,to);
127: ISCreateStride(da->comm,m,start,1,&from);
128: VecScatterCreate(da->global,from,da->natural,to,&da->gton);
129: ISDestroy(from);
130: ISDestroy(to);
131: return(0);
132: }
134: /*@
135: DAGlobalToNaturalBegin - Maps values from the global vector to a global vector
136: in the "natural" grid ordering. Must be followed by
137: DAGlobalToNaturalEnd() to complete the exchange.
139: Collective on DA
141: Input Parameters:
142: + da - the distributed array context
143: . g - the global vector
144: - mode - one of INSERT_VALUES or ADD_VALUES
146: Output Parameter:
147: . l - the natural ordering values
149: Level: advanced
151: Notes:
152: The global and natrual vectors used here need not be the same as those
153: obtained from DACreateGlobalVector() and DACreateNaturalVector(), BUT they
154: must have the same parallel data layout; they could, for example, be
155: obtained with VecDuplicate() from the DA originating vectors.
157: .keywords: distributed array, global to local, begin
159: .seealso: DAGlobalToNaturalEnd(), DALocalToGlobal(), DACreate2d(),
160: DAGlobalToLocalBegin(), DAGlobalToLocalEnd(), DACreateNaturalVector()
161: @*/
162: int DAGlobalToNaturalBegin(DA da,Vec g,InsertMode mode,Vec l)
163: {
170: if (!da->gton) {
171: /* create the scatter context */
172: DAGlobalToNatural_Create(da);
173: }
174: VecScatterBegin(g,l,mode,SCATTER_FORWARD,da->gton);
175: return(0);
176: }
178: /*@
179: DAGlobalToNaturalEnd - Maps values from the global vector to a global vector
180: in the natural ordering. Must be preceeded by DAGlobalToNaturalBegin().
182: Collective on DA
184: Input Parameters:
185: + da - the distributed array context
186: . g - the global vector
187: - mode - one of INSERT_VALUES or ADD_VALUES
189: Output Parameter:
190: . l - the global values in the natural ordering
192: Level: advanced
194: Notes:
195: The global and local vectors used here need not be the same as those
196: obtained from DACreateGlobalVector() and DACreateNaturalVector(), BUT they
197: must have the same parallel data layout; they could, for example, be
198: obtained with VecDuplicate() from the DA originating vectors.
200: .keywords: distributed array, global to local, end
202: .seealso: DAGlobalToNaturalBegin(), DALocalToGlobal(), DACreate2d(),
203: DAGlobalToLocalBegin(), DAGlobalToLocalEnd(), DACreateNaturalVector()
204: @*/
205: int DAGlobalToNaturalEnd(DA da,Vec g,InsertMode mode,Vec l)
206: {
213: VecScatterEnd(g,l,mode,SCATTER_FORWARD,da->gton);
214: return(0);
215: }
217: /*@
218: DANaturalToGlobalBegin - Maps values from a global vector in the "natural" ordering
219: to a global vector in the PETSc DA grid ordering. Must be followed by
220: DANaturalToGlobalEnd() to complete the exchange.
222: Collective on DA
224: Input Parameters:
225: + da - the distributed array context
226: . g - the global vector in a natural ordering
227: - mode - one of INSERT_VALUES or ADD_VALUES
229: Output Parameter:
230: . l - the values in the DA ordering
232: Level: advanced
234: Notes:
235: The global and natural vectors used here need not be the same as those
236: obtained from DACreateGlobalVector() and DACreateNaturalVector(), BUT they
237: must have the same parallel data layout; they could, for example, be
238: obtained with VecDuplicate() from the DA originating vectors.
240: .keywords: distributed array, global to local, begin
242: .seealso: DAGlobalToNaturalEnd(), DAGlobalToNaturalBegin(), DALocalToGlobal(), DACreate2d(),
243: DAGlobalToLocalBegin(), DAGlobalToLocalEnd(), DACreateNaturalVector()
244: @*/
245: int DANaturalToGlobalBegin(DA da,Vec g,InsertMode mode,Vec l)
246: {
253: if (!da->gton) {
254: /* create the scatter context */
255: DAGlobalToNatural_Create(da);
256: }
257: VecScatterBegin(g,l,mode,SCATTER_REVERSE,da->gton);
258: return(0);
259: }
261: /*@
262: DANaturalToGlobalEnd - Maps values from the natural ordering global vector
263: to a global vector in the PETSc DA ordering. Must be preceeded by DANaturalToGlobalBegin().
265: Collective on DA
267: Input Parameters:
268: + da - the distributed array context
269: . g - the global vector in a natural ordering
270: - mode - one of INSERT_VALUES or ADD_VALUES
272: Output Parameter:
273: . l - the global values in the PETSc DA ordering
275: Level: intermediate
277: Notes:
278: The global and local vectors used here need not be the same as those
279: obtained from DACreateGlobalVector() and DACreateNaturalVector(), BUT they
280: must have the same parallel data layout; they could, for example, be
281: obtained with VecDuplicate() from the DA originating vectors.
283: .keywords: distributed array, global to local, end
285: .seealso: DAGlobalToNaturalBegin(), DAGlobalToNaturalEnd(), DALocalToGlobal(), DACreate2d(),
286: DAGlobalToLocalBegin(), DAGlobalToLocalEnd(), DACreateNaturalVector()
287: @*/
288: int DANaturalToGlobalEnd(DA da,Vec g,InsertMode mode,Vec l)
289: {
296: VecScatterEnd(g,l,mode,SCATTER_REVERSE,da->gton);
297: return(0);
298: }