Actual source code: commonmpvec.c
2: #include <../src/vec/vec/impls/mpi/pvecimpl.h> /*I "petscvec.h" I*/
6: /*
7: This is used in VecGhostGetLocalForm and VecGhostRestoreLocalForm to ensure
8: that the state is updated if either vector has changed since the last time
9: one of these functions was called. It could apply to any PetscObject, but
10: VecGhost is quite different from other objects in that two separate vectors
11: look at the same memory.
13: In principle, we could only propagate state to the local vector on
14: GetLocalForm and to the global vector on RestoreLocalForm, but this version is
15: more conservative (i.e. robust against misuse) and simpler.
17: Note that this function is correct and changes nothing if both arguments are the
18: same, which is the case in serial.
19: */
20: static PetscErrorCode VecGhostStateSync_Private(Vec g,Vec l)
21: {
23: PetscInt gstate,lstate;
26: PetscObjectStateQuery((PetscObject)g,&gstate);
27: PetscObjectStateQuery((PetscObject)l,&lstate);
28: PetscObjectSetState((PetscObject)g,PetscMax(gstate,lstate));
29: PetscObjectSetState((PetscObject)l,PetscMax(gstate,lstate));
30: return(0);
31: }
35: /*@
36: VecGhostGetLocalForm - Obtains the local ghosted representation of
37: a parallel vector created with VecCreateGhost().
39: Not Collective
41: Input Parameter:
42: . g - the global vector. Vector must be have been obtained with either
43: VecCreateGhost(), VecCreateGhostWithArray() or VecCreateSeq().
45: Output Parameter:
46: . l - the local (ghosted) representation
48: Notes:
49: This routine does not actually update the ghost values, but rather it
50: returns a sequential vector that includes the locations for the ghost
51: values and their current values. The returned vector and the original
52: vector passed in share the same array that contains the actual vector data.
54: One should call VecGhostRestoreLocalForm() or VecDestroy() once one is
55: finished using the object.
57: Level: advanced
59: Concepts: vectors^ghost point access
61: .seealso: VecCreateGhost(), VecGhostRestoreLocalForm(), VecCreateGhostWithArray()
63: @*/
64: PetscErrorCode VecGhostGetLocalForm(Vec g,Vec *l)
65: {
67: PetscBool isseq,ismpi;
73: PetscTypeCompare((PetscObject)g,VECSEQ,&isseq);
74: PetscTypeCompare((PetscObject)g,VECMPI,&ismpi);
75: if (ismpi) {
76: Vec_MPI *v = (Vec_MPI*)g->data;
77: if (!v->localrep) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Vector is not ghosted");
78: *l = v->localrep;
79: } else if (isseq) {
80: *l = g;
81: } else {
82: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Vector type %s does not have local representation",((PetscObject)g)->type_name);
83: }
84: VecGhostStateSync_Private(g,*l);
85: PetscObjectReference((PetscObject)*l);
86: return(0);
87: }
91: /*@
92: VecGhostRestoreLocalForm - Restores the local ghosted representation of
93: a parallel vector obtained with VecGhostGetLocalForm().
95: Not Collective
97: Input Parameter:
98: + g - the global vector
99: - l - the local (ghosted) representation
101: Notes:
102: This routine does not actually update the ghost values, but rather it
103: returns a sequential vector that includes the locations for the ghost values
104: and their current values.
106: Level: advanced
108: .seealso: VecCreateGhost(), VecGhostGetLocalForm(), VecCreateGhostWithArray()
109: @*/
110: PetscErrorCode VecGhostRestoreLocalForm(Vec g,Vec *l)
111: {
115: VecGhostStateSync_Private(g,*l);
116: PetscObjectDereference((PetscObject)*l);
117: return(0);
118: }
122: /*@
123: VecGhostUpdateBegin - Begins the vector scatter to update the vector from
124: local representation to global or global representation to local.
126: Neighbor-wise Collective on Vec
128: Input Parameters:
129: + g - the vector (obtained with VecCreateGhost() or VecDuplicate())
130: . insertmode - one of ADD_VALUES or INSERT_VALUES
131: - scattermode - one of SCATTER_FORWARD or SCATTER_REVERSE
133: Notes:
134: Use the following to update the ghost regions with correct values from the owning process
135: .vb
136: VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
137: VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
138: .ve
140: Use the following to accumulate the ghost region values onto the owning processors
141: .vb
142: VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
143: VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
144: .ve
146: To accumulate the ghost region values onto the owning processors and then update
147: the ghost regions correctly, call the later followed by the former, i.e.,
148: .vb
149: VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
150: VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
151: VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
152: VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
153: .ve
155: Level: advanced
157: .seealso: VecCreateGhost(), VecGhostUpdateEnd(), VecGhostGetLocalForm(),
158: VecGhostRestoreLocalForm(),VecCreateGhostWithArray()
160: @*/
161: PetscErrorCode VecGhostUpdateBegin(Vec g,InsertMode insertmode,ScatterMode scattermode)
162: {
163: Vec_MPI *v;
169: v = (Vec_MPI*)g->data;
170: if (!v->localrep) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Vector is not ghosted");
171: if (!v->localupdate) return(0);
172:
173: if (scattermode == SCATTER_REVERSE) {
174: VecScatterBegin(v->localupdate,v->localrep,g,insertmode,scattermode);
175: } else {
176: VecScatterBegin(v->localupdate,g,v->localrep,insertmode,scattermode);
177: }
178: return(0);
179: }
183: /*@
184: VecGhostUpdateEnd - End the vector scatter to update the vector from
185: local representation to global or global representation to local.
187: Neighbor-wise Collective on Vec
189: Input Parameters:
190: + g - the vector (obtained with VecCreateGhost() or VecDuplicate())
191: . insertmode - one of ADD_VALUES or INSERT_VALUES
192: - scattermode - one of SCATTER_FORWARD or SCATTER_REVERSE
194: Notes:
196: Use the following to update the ghost regions with correct values from the owning process
197: .vb
198: VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
199: VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
200: .ve
202: Use the following to accumulate the ghost region values onto the owning processors
203: .vb
204: VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
205: VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
206: .ve
208: To accumulate the ghost region values onto the owning processors and then update
209: the ghost regions correctly, call the later followed by the former, i.e.,
210: .vb
211: VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
212: VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
213: VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
214: VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
215: .ve
217: Level: advanced
219: .seealso: VecCreateGhost(), VecGhostUpdateBegin(), VecGhostGetLocalForm(),
220: VecGhostRestoreLocalForm(),VecCreateGhostWithArray()
222: @*/
223: PetscErrorCode VecGhostUpdateEnd(Vec g,InsertMode insertmode,ScatterMode scattermode)
224: {
225: Vec_MPI *v;
231: v = (Vec_MPI*)g->data;
232: if (!v->localrep) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Vector is not ghosted");
233: if (!v->localupdate) return(0);
235: if (scattermode == SCATTER_REVERSE) {
236: VecScatterEnd(v->localupdate,v->localrep,g,insertmode,scattermode);
237: } else {
238: VecScatterEnd(v->localupdate,g,v->localrep,insertmode,scattermode);
239: }
240: return(0);
241: }
245: PetscErrorCode VecSetOption_MPI(Vec v,VecOption op,PetscBool flag)
246: {
248: if (op == VEC_IGNORE_OFF_PROC_ENTRIES) {
249: v->stash.donotstash = flag;
250: } else if (op == VEC_IGNORE_NEGATIVE_INDICES) {
251: v->stash.ignorenegidx = flag;
252: }
253: return(0);
254: }
259: PetscErrorCode VecResetArray_MPI(Vec vin)
260: {
261: Vec_MPI *v = (Vec_MPI *)vin->data;
265: v->array = v->unplacedarray;
266: v->unplacedarray = 0;
267: if (v->localrep) {
268: VecResetArray(v->localrep);
269: }
270: return(0);
271: }