Actual source code: commonmpvec.c
2: #include <../src/vec/vec/impls/mpi/pvecimpl.h>
4: /*
5: This is used in VecGhostGetLocalForm and VecGhostRestoreLocalForm to ensure
6: that the state is updated if either vector has changed since the last time
7: one of these functions was called. It could apply to any PetscObject, but
8: VecGhost is quite different from other objects in that two separate vectors
9: look at the same memory.
11: In principle, we could only propagate state to the local vector on
12: GetLocalForm and to the global vector on RestoreLocalForm, but this version is
13: more conservative (i.e. robust against misuse) and simpler.
15: Note that this function is correct and changes nothing if both arguments are the
16: same, which is the case in serial.
17: */
18: static PetscErrorCode VecGhostStateSync_Private(Vec g,Vec l)
19: {
20: PetscErrorCode ierr;
21: PetscObjectState gstate,lstate;
24: PetscObjectStateGet((PetscObject)g,&gstate);
25: PetscObjectStateGet((PetscObject)l,&lstate);
26: PetscObjectStateSet((PetscObject)g,PetscMax(gstate,lstate));
27: PetscObjectStateSet((PetscObject)l,PetscMax(gstate,lstate));
28: return(0);
29: }
31: /*@
32: VecGhostGetLocalForm - Obtains the local ghosted representation of
33: a parallel vector (obtained with VecCreateGhost(), VecCreateGhostWithArray()
34: or VecCreateSeq()). Returns NULL if the Vec is not ghosted.
36: Logically Collective
38: Input Parameter:
39: . g - the global vector
41: Output Parameter:
42: . l - the local (ghosted) representation, NULL if g is not ghosted
44: Notes:
45: This routine does not actually update the ghost values, but rather it
46: returns a sequential vector that includes the locations for the ghost
47: values and their current values. The returned vector and the original
48: vector passed in share the same array that contains the actual vector data.
50: To update the ghost values from the locations on the other processes one must call
51: VecGhostUpdateBegin() and VecGhostUpdateEnd() before accessing the ghost values. Thus normal
52: usage is
53: $ VecGhostUpdateBegin(x,INSERT_VALUES,SCATTER_FORWARD);
54: $ VecGhostUpdateEnd(x,INSERT_VALUES,SCATTER_FORWARD);
55: $ VecGhostGetLocalForm(x,&xlocal);
56: $ VecGetArray(xlocal,&xvalues);
57: $ // access the non-ghost values in locations xvalues[0:n-1] and ghost values in locations xvalues[n:n+nghost];
58: $ VecRestoreArray(xlocal,&xvalues);
59: $ VecGhostRestoreLocalForm(x,&xlocal);
61: One should call VecGhostRestoreLocalForm() or VecDestroy() once one is
62: finished using the object.
64: Level: advanced
66: .seealso: VecCreateGhost(), VecGhostRestoreLocalForm(), VecCreateGhostWithArray()
68: @*/
69: PetscErrorCode VecGhostGetLocalForm(Vec g,Vec *l)
70: {
72: PetscBool isseq,ismpi;
78: PetscObjectTypeCompare((PetscObject)g,VECSEQ,&isseq);
79: PetscObjectTypeCompare((PetscObject)g,VECMPI,&ismpi);
80: if (ismpi) {
81: Vec_MPI *v = (Vec_MPI*)g->data;
82: *l = v->localrep;
83: } else if (isseq) {
84: *l = g;
85: } else {
86: *l = NULL;
87: }
88: if (*l) {
89: VecGhostStateSync_Private(g,*l);
90: PetscObjectReference((PetscObject)*l);
91: }
92: return(0);
93: }
95: /*@
96: VecGhostIsLocalForm - Checks if a given vector is the local form of a global vector
98: Not Collective
100: Input Parameters:
101: + g - the global vector
102: - l - the local vector
104: Output Parameter:
105: . flg - PETSC_TRUE if local vector is local form
107: Level: advanced
109: .seealso: VecCreateGhost(), VecGhostRestoreLocalForm(), VecCreateGhostWithArray(), VecGhostGetLocalForm()
111: @*/
112: PetscErrorCode VecGhostIsLocalForm(Vec g,Vec l,PetscBool *flg)
113: {
115: PetscBool isseq,ismpi;
121: *flg = PETSC_FALSE;
122: PetscObjectTypeCompare((PetscObject)g,VECSEQ,&isseq);
123: PetscObjectTypeCompare((PetscObject)g,VECMPI,&ismpi);
124: if (ismpi) {
125: Vec_MPI *v = (Vec_MPI*)g->data;
126: if (l == v->localrep) *flg = PETSC_TRUE;
127: } else if (isseq) {
128: if (l == g) *flg = PETSC_TRUE;
129: } else SETERRQ(PetscObjectComm((PetscObject)g),PETSC_ERR_ARG_WRONG,"Global vector is not ghosted");
130: return(0);
131: }
133: /*@
134: VecGhostRestoreLocalForm - Restores the local ghosted representation of
135: a parallel vector obtained with VecGhostGetLocalForm().
137: Logically Collective
139: Input Parameters:
140: + g - the global vector
141: - l - the local (ghosted) representation
143: Notes:
144: This routine does not actually update the ghost values, but rather it
145: returns a sequential vector that includes the locations for the ghost values
146: and their current values.
148: Level: advanced
150: .seealso: VecCreateGhost(), VecGhostGetLocalForm(), VecCreateGhostWithArray()
151: @*/
152: PetscErrorCode VecGhostRestoreLocalForm(Vec g,Vec *l)
153: {
157: if (*l) {
158: VecGhostStateSync_Private(g,*l);
159: PetscObjectDereference((PetscObject)*l);
160: }
161: return(0);
162: }
164: /*@
165: VecGhostUpdateBegin - Begins the vector scatter to update the vector from
166: local representation to global or global representation to local.
168: Neighbor-wise Collective on Vec
170: Input Parameters:
171: + g - the vector (obtained with VecCreateGhost() or VecDuplicate())
172: . insertmode - one of ADD_VALUES, MAX_VALUES, MIN_VALUES or INSERT_VALUES
173: - scattermode - one of SCATTER_FORWARD or SCATTER_REVERSE
175: Notes:
176: Use the following to update the ghost regions with correct values from the owning process
177: .vb
178: VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
179: VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
180: .ve
182: Use the following to accumulate the ghost region values onto the owning processors
183: .vb
184: VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
185: VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
186: .ve
188: To accumulate the ghost region values onto the owning processors and then update
189: the ghost regions correctly, call the latter followed by the former, i.e.,
190: .vb
191: VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
192: VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
193: VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
194: VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
195: .ve
197: Level: advanced
199: .seealso: VecCreateGhost(), VecGhostUpdateEnd(), VecGhostGetLocalForm(),
200: VecGhostRestoreLocalForm(),VecCreateGhostWithArray()
202: @*/
203: PetscErrorCode VecGhostUpdateBegin(Vec g,InsertMode insertmode,ScatterMode scattermode)
204: {
205: Vec_MPI *v;
207: PetscBool ismpi,isseq;
211: PetscObjectTypeCompare((PetscObject)g,VECMPI,&ismpi);
212: PetscObjectTypeCompare((PetscObject)g,VECSEQ,&isseq);
213: if (ismpi) {
214: v = (Vec_MPI*)g->data;
215: if (!v->localrep) SETERRQ(PetscObjectComm((PetscObject)g),PETSC_ERR_ARG_WRONG,"Vector is not ghosted");
216: if (!v->localupdate) return(0);
217: if (scattermode == SCATTER_REVERSE) {
218: VecScatterBegin(v->localupdate,v->localrep,g,insertmode,scattermode);
219: } else {
220: VecScatterBegin(v->localupdate,g,v->localrep,insertmode,scattermode);
221: }
222: } else if (isseq) {
223: /* Do nothing */
224: } else SETERRQ(PetscObjectComm((PetscObject)g),PETSC_ERR_ARG_WRONG,"Vector is not ghosted");
225: return(0);
226: }
228: /*@
229: VecGhostUpdateEnd - End the vector scatter to update the vector from
230: local representation to global or global representation to local.
232: Neighbor-wise Collective on Vec
234: Input Parameters:
235: + g - the vector (obtained with VecCreateGhost() or VecDuplicate())
236: . insertmode - one of ADD_VALUES, MAX_VALUES, MIN_VALUES or INSERT_VALUES
237: - scattermode - one of SCATTER_FORWARD or SCATTER_REVERSE
239: Notes:
241: Use the following to update the ghost regions with correct values from the owning process
242: .vb
243: VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
244: VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
245: .ve
247: Use the following to accumulate the ghost region values onto the owning processors
248: .vb
249: VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
250: VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
251: .ve
253: To accumulate the ghost region values onto the owning processors and then update
254: the ghost regions correctly, call the later followed by the former, i.e.,
255: .vb
256: VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
257: VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
258: VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
259: VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
260: .ve
262: Level: advanced
264: .seealso: VecCreateGhost(), VecGhostUpdateBegin(), VecGhostGetLocalForm(),
265: VecGhostRestoreLocalForm(),VecCreateGhostWithArray()
267: @*/
268: PetscErrorCode VecGhostUpdateEnd(Vec g,InsertMode insertmode,ScatterMode scattermode)
269: {
270: Vec_MPI *v;
272: PetscBool ismpi;
276: PetscObjectTypeCompare((PetscObject)g,VECMPI,&ismpi);
277: if (ismpi) {
278: v = (Vec_MPI*)g->data;
279: if (!v->localrep) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Vector is not ghosted");
280: if (!v->localupdate) return(0);
281: if (scattermode == SCATTER_REVERSE) {
282: VecScatterEnd(v->localupdate,v->localrep,g,insertmode,scattermode);
283: } else {
284: VecScatterEnd(v->localupdate,g,v->localrep,insertmode,scattermode);
285: }
286: }
287: return(0);
288: }