Actual source code: commonmpvec.c
petsc-3.14.6 2021-03-30
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: }
32: /*@
33: VecGhostGetLocalForm - Obtains the local ghosted representation of
34: a parallel vector (obtained with VecCreateGhost(), VecCreateGhostWithArray()
35: or VecCreateSeq()). Returns NULL if the Vec is not ghosted.
37: Logically Collective
39: Input Parameter:
40: . g - the global vector
42: Output Parameter:
43: . l - the local (ghosted) representation, NULL if g is not ghosted
45: Notes:
46: This routine does not actually update the ghost values, but rather it
47: returns a sequential vector that includes the locations for the ghost
48: values and their current values. The returned vector and the original
49: vector passed in share the same array that contains the actual vector data.
51: To update the ghost values from the locations on the other processes one must call
52: VecGhostUpdateBegin() and VecGhostUpdateEnd() before accessing the ghost values. Thus normal
53: usage is
54: $ VecGhostUpdateBegin(x,INSERT_VALUES,SCATTER_FORWARD);
55: $ VecGhostUpdateEnd(x,INSERT_VALUES,SCATTER_FORWARD);
56: $ VecGhostGetLocalForm(x,&xlocal);
57: $ VecGetArray(xlocal,&xvalues);
58: $ // access the non-ghost values in locations xvalues[0:n-1] and ghost values in locations xvalues[n:n+nghost];
59: $ VecRestoreArray(xlocal,&xvalues);
60: $ VecGhostRestoreLocalForm(x,&xlocal);
62: One should call VecGhostRestoreLocalForm() or VecDestroy() once one is
63: finished using the object.
65: Level: advanced
67: .seealso: VecCreateGhost(), VecGhostRestoreLocalForm(), VecCreateGhostWithArray()
69: @*/
70: PetscErrorCode VecGhostGetLocalForm(Vec g,Vec *l)
71: {
73: PetscBool isseq,ismpi;
79: PetscObjectTypeCompare((PetscObject)g,VECSEQ,&isseq);
80: PetscObjectTypeCompare((PetscObject)g,VECMPI,&ismpi);
81: if (ismpi) {
82: Vec_MPI *v = (Vec_MPI*)g->data;
83: *l = v->localrep;
84: } else if (isseq) {
85: *l = g;
86: } else {
87: *l = NULL;
88: }
89: if (*l) {
90: VecGhostStateSync_Private(g,*l);
91: PetscObjectReference((PetscObject)*l);
92: }
93: return(0);
94: }
96: /*@
97: VecGhostIsLocalForm - Checks if a given vector is the local form of a global vector
99: Not Collective
101: Input Parameter:
102: + g - the global vector
103: - l - the local vector
105: Output Parameter:
106: . flg - PETSC_TRUE if local vector is local form
108: Level: advanced
110: .seealso: VecCreateGhost(), VecGhostRestoreLocalForm(), VecCreateGhostWithArray(), VecGhostGetLocalForm()
112: @*/
113: PetscErrorCode VecGhostIsLocalForm(Vec g,Vec l,PetscBool *flg)
114: {
116: PetscBool isseq,ismpi;
122: *flg = PETSC_FALSE;
123: PetscObjectTypeCompare((PetscObject)g,VECSEQ,&isseq);
124: PetscObjectTypeCompare((PetscObject)g,VECMPI,&ismpi);
125: if (ismpi) {
126: Vec_MPI *v = (Vec_MPI*)g->data;
127: if (l == v->localrep) *flg = PETSC_TRUE;
128: } else if (isseq) {
129: if (l == g) *flg = PETSC_TRUE;
130: } else SETERRQ(PetscObjectComm((PetscObject)g),PETSC_ERR_ARG_WRONG,"Global vector is not ghosted");
131: return(0);
132: }
134: /*@
135: VecGhostRestoreLocalForm - Restores the local ghosted representation of
136: a parallel vector obtained with VecGhostGetLocalForm().
138: Logically Collective
140: Input Parameter:
141: + g - the global vector
142: - l - the local (ghosted) representation
144: Notes:
145: This routine does not actually update the ghost values, but rather it
146: returns a sequential vector that includes the locations for the ghost values
147: and their current values.
149: Level: advanced
151: .seealso: VecCreateGhost(), VecGhostGetLocalForm(), VecCreateGhostWithArray()
152: @*/
153: PetscErrorCode VecGhostRestoreLocalForm(Vec g,Vec *l)
154: {
158: if (*l) {
159: VecGhostStateSync_Private(g,*l);
160: PetscObjectDereference((PetscObject)*l);
161: }
162: return(0);
163: }
165: /*@
166: VecGhostUpdateBegin - Begins the vector scatter to update the vector from
167: local representation to global or global representation to local.
169: Neighbor-wise Collective on Vec
171: Input Parameters:
172: + g - the vector (obtained with VecCreateGhost() or VecDuplicate())
173: . insertmode - one of ADD_VALUES, MAX_VALUES, MIN_VALUES or INSERT_VALUES
174: - scattermode - one of SCATTER_FORWARD or SCATTER_REVERSE
176: Notes:
177: Use the following to update the ghost regions with correct values from the owning process
178: .vb
179: VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
180: VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
181: .ve
183: Use the following to accumulate the ghost region values onto the owning processors
184: .vb
185: VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
186: VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
187: .ve
189: To accumulate the ghost region values onto the owning processors and then update
190: the ghost regions correctly, call the latter followed by the former, i.e.,
191: .vb
192: VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
193: VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
194: VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
195: VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
196: .ve
198: Level: advanced
200: .seealso: VecCreateGhost(), VecGhostUpdateEnd(), VecGhostGetLocalForm(),
201: VecGhostRestoreLocalForm(),VecCreateGhostWithArray()
203: @*/
204: PetscErrorCode VecGhostUpdateBegin(Vec g,InsertMode insertmode,ScatterMode scattermode)
205: {
206: Vec_MPI *v;
208: PetscBool ismpi,isseq;
212: PetscObjectTypeCompare((PetscObject)g,VECMPI,&ismpi);
213: PetscObjectTypeCompare((PetscObject)g,VECSEQ,&isseq);
214: if (ismpi) {
215: v = (Vec_MPI*)g->data;
216: if (!v->localrep) SETERRQ(PetscObjectComm((PetscObject)g),PETSC_ERR_ARG_WRONG,"Vector is not ghosted");
217: if (!v->localupdate) return(0);
218: if (scattermode == SCATTER_REVERSE) {
219: VecScatterBegin(v->localupdate,v->localrep,g,insertmode,scattermode);
220: } else {
221: VecScatterBegin(v->localupdate,g,v->localrep,insertmode,scattermode);
222: }
223: } else if (isseq) {
224: /* Do nothing */
225: } else SETERRQ(PetscObjectComm((PetscObject)g),PETSC_ERR_ARG_WRONG,"Vector is not ghosted");
226: return(0);
227: }
229: /*@
230: VecGhostUpdateEnd - End the vector scatter to update the vector from
231: local representation to global or global representation to local.
233: Neighbor-wise Collective on Vec
235: Input Parameters:
236: + g - the vector (obtained with VecCreateGhost() or VecDuplicate())
237: . insertmode - one of ADD_VALUES, MAX_VALUES, MIN_VALUES or INSERT_VALUES
238: - scattermode - one of SCATTER_FORWARD or SCATTER_REVERSE
240: Notes:
242: Use the following to update the ghost regions with correct values from the owning process
243: .vb
244: VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
245: VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
246: .ve
248: Use the following to accumulate the ghost region values onto the owning processors
249: .vb
250: VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
251: VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
252: .ve
254: To accumulate the ghost region values onto the owning processors and then update
255: the ghost regions correctly, call the later followed by the former, i.e.,
256: .vb
257: VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
258: VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
259: VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
260: VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
261: .ve
263: Level: advanced
265: .seealso: VecCreateGhost(), VecGhostUpdateBegin(), VecGhostGetLocalForm(),
266: VecGhostRestoreLocalForm(),VecCreateGhostWithArray()
268: @*/
269: PetscErrorCode VecGhostUpdateEnd(Vec g,InsertMode insertmode,ScatterMode scattermode)
270: {
271: Vec_MPI *v;
273: PetscBool ismpi;
277: PetscObjectTypeCompare((PetscObject)g,VECMPI,&ismpi);
278: if (ismpi) {
279: v = (Vec_MPI*)g->data;
280: if (!v->localrep) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Vector is not ghosted");
281: if (!v->localupdate) return(0);
282: if (scattermode == SCATTER_REVERSE) {
283: VecScatterEnd(v->localupdate,v->localrep,g,insertmode,scattermode);
284: } else {
285: VecScatterEnd(v->localupdate,g,v->localrep,insertmode,scattermode);
286: }
287: }
288: return(0);
289: }