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:   PetscObjectState gstate, lstate;

 22:   PetscFunctionBegin;
 23:   PetscCall(PetscObjectStateGet((PetscObject)g, &gstate));
 24:   PetscCall(PetscObjectStateGet((PetscObject)l, &lstate));
 25:   PetscCall(PetscObjectStateSet((PetscObject)g, PetscMax(gstate, lstate)));
 26:   PetscCall(PetscObjectStateSet((PetscObject)l, PetscMax(gstate, lstate)));
 27:   PetscFunctionReturn(PETSC_SUCCESS);
 28: }

 30: /*@
 31:     VecGhostGetLocalForm - Obtains the local ghosted representation of
 32:     a parallel vector (obtained with VecCreateGhost(), VecCreateGhostWithArray()
 33:     or VecCreateSeq()). Returns NULL if the Vec is not ghosted.

 35:     Logically Collective

 37:     Input Parameter:
 38: .   g - the global vector

 40:     Output Parameter:
 41: .   l - the local (ghosted) representation, NULL if g is not ghosted

 43:     Notes:
 44:     This routine does not actually update the ghost values, but rather it
 45:     returns a sequential vector that includes the locations for the ghost
 46:     values and their current values. The returned vector and the original
 47:     vector passed in share the same array that contains the actual vector data.

 49:     To update the ghost values from the locations on the other processes one must call
 50:     VecGhostUpdateBegin() and VecGhostUpdateEnd() before accessing the ghost values. Thus normal
 51:     usage is
 52: $     VecGhostUpdateBegin(x,INSERT_VALUES,SCATTER_FORWARD);
 53: $     VecGhostUpdateEnd(x,INSERT_VALUES,SCATTER_FORWARD);
 54: $     VecGhostGetLocalForm(x,&xlocal);
 55: $     VecGetArray(xlocal,&xvalues);
 56: $        // access the non-ghost values in locations xvalues[0:n-1] and ghost values in locations xvalues[n:n+nghost];
 57: $     VecRestoreArray(xlocal,&xvalues);
 58: $     VecGhostRestoreLocalForm(x,&xlocal);

 60:     One should call VecGhostRestoreLocalForm() or VecDestroy() once one is
 61:     finished using the object.

 63:     Level: advanced

 65: .seealso: `VecCreateGhost()`, `VecGhostRestoreLocalForm()`, `VecCreateGhostWithArray()`

 67: @*/
 68: PetscErrorCode VecGhostGetLocalForm(Vec g, Vec *l)
 69: {
 70:   PetscBool isseq, ismpi;

 72:   PetscFunctionBegin;

 76:   PetscCall(PetscObjectTypeCompare((PetscObject)g, VECSEQ, &isseq));
 77:   PetscCall(PetscObjectTypeCompare((PetscObject)g, VECMPI, &ismpi));
 78:   if (ismpi) {
 79:     Vec_MPI *v = (Vec_MPI *)g->data;
 80:     *l         = v->localrep;
 81:   } else if (isseq) {
 82:     *l = g;
 83:   } else {
 84:     *l = NULL;
 85:   }
 86:   if (*l) {
 87:     PetscCall(VecGhostStateSync_Private(g, *l));
 88:     PetscCall(PetscObjectReference((PetscObject)*l));
 89:   }
 90:   PetscFunctionReturn(PETSC_SUCCESS);
 91: }

 93: /*@
 94:     VecGhostIsLocalForm - Checks if a given vector is the local form of a global vector

 96:     Not Collective

 98:     Input Parameters:
 99: +   g - the global vector
100: -   l - the local vector

102:     Output Parameter:
103: .   flg - PETSC_TRUE if local vector is local form

105:     Level: advanced

107: .seealso: `VecCreateGhost()`, `VecGhostRestoreLocalForm()`, `VecCreateGhostWithArray()`, `VecGhostGetLocalForm()`

109: @*/
110: PetscErrorCode VecGhostIsLocalForm(Vec g, Vec l, PetscBool *flg)
111: {
112:   PetscBool isseq, ismpi;

114:   PetscFunctionBegin;

118:   *flg = PETSC_FALSE;
119:   PetscCall(PetscObjectTypeCompare((PetscObject)g, VECSEQ, &isseq));
120:   PetscCall(PetscObjectTypeCompare((PetscObject)g, VECMPI, &ismpi));
121:   if (ismpi) {
122:     Vec_MPI *v = (Vec_MPI *)g->data;
123:     if (l == v->localrep) *flg = PETSC_TRUE;
124:   } else if (isseq) {
125:     if (l == g) *flg = PETSC_TRUE;
126:   } else SETERRQ(PetscObjectComm((PetscObject)g), PETSC_ERR_ARG_WRONG, "Global vector is not ghosted");
127:   PetscFunctionReturn(PETSC_SUCCESS);
128: }

130: /*@
131:     VecGhostRestoreLocalForm - Restores the local ghosted representation of
132:     a parallel vector obtained with VecGhostGetLocalForm().

134:     Logically Collective

136:     Input Parameters:
137: +   g - the global vector
138: -   l - the local (ghosted) representation

140:     Notes:
141:     This routine does not actually update the ghost values, but rather it
142:     returns a sequential vector that includes the locations for the ghost values
143:     and their current values.

145:     Level: advanced

147: .seealso: `VecCreateGhost()`, `VecGhostGetLocalForm()`, `VecCreateGhostWithArray()`
148: @*/
149: PetscErrorCode VecGhostRestoreLocalForm(Vec g, Vec *l)
150: {
151:   PetscFunctionBegin;
152:   if (*l) {
153:     PetscCall(VecGhostStateSync_Private(g, *l));
154:     PetscCall(PetscObjectDereference((PetscObject)*l));
155:   }
156:   PetscFunctionReturn(PETSC_SUCCESS);
157: }

159: /*@
160:    VecGhostUpdateBegin - Begins the vector scatter to update the vector from
161:    local representation to global or global representation to local.

163:    Neighbor-wise Collective

165:    Input Parameters:
166: +  g - the vector (obtained with `VecCreateGhost()` or `VecDuplicate()`)
167: .  insertmode - one of ADD_VALUES, MAX_VALUES, MIN_VALUES or INSERT_VALUES
168: -  scattermode - one of SCATTER_FORWARD or SCATTER_REVERSE

170:    Notes:
171:    Use the following to update the ghost regions with correct values from the owning process
172: .vb
173:        VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
174:        VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
175: .ve

177:    Use the following to accumulate the ghost region values onto the owning processors
178: .vb
179:        VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
180:        VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
181: .ve

183:    To accumulate the ghost region values onto the owning processors and then update
184:    the ghost regions correctly, call the latter followed by the former, i.e.,
185: .vb
186:        VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
187:        VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
188:        VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
189:        VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
190: .ve

192:    Level: advanced

194: .seealso: `VecCreateGhost()`, `VecGhostUpdateEnd()`, `VecGhostGetLocalForm()`,
195:           `VecGhostRestoreLocalForm()`, `VecCreateGhostWithArray()`

197: @*/
198: PetscErrorCode VecGhostUpdateBegin(Vec g, InsertMode insertmode, ScatterMode scattermode)
199: {
200:   Vec_MPI  *v;
201:   PetscBool ismpi, isseq;

203:   PetscFunctionBegin;
205:   PetscCall(PetscObjectTypeCompare((PetscObject)g, VECMPI, &ismpi));
206:   PetscCall(PetscObjectTypeCompare((PetscObject)g, VECSEQ, &isseq));
207:   if (ismpi) {
208:     v = (Vec_MPI *)g->data;
209:     PetscCheck(v->localrep, PetscObjectComm((PetscObject)g), PETSC_ERR_ARG_WRONG, "Vector is not ghosted");
210:     if (!v->localupdate) PetscFunctionReturn(PETSC_SUCCESS);
211:     if (scattermode == SCATTER_REVERSE) {
212:       PetscCall(VecScatterBegin(v->localupdate, v->localrep, g, insertmode, scattermode));
213:     } else {
214:       PetscCall(VecScatterBegin(v->localupdate, g, v->localrep, insertmode, scattermode));
215:     }
216:   } else if (isseq) {
217:     /* Do nothing */
218:   } else SETERRQ(PetscObjectComm((PetscObject)g), PETSC_ERR_ARG_WRONG, "Vector is not ghosted");
219:   PetscFunctionReturn(PETSC_SUCCESS);
220: }

222: /*@
223:    VecGhostUpdateEnd - End the vector scatter to update the vector from
224:    local representation to global or global representation to local.

226:    Neighbor-wise Collective

228:    Input Parameters:
229: +  g - the vector (obtained with `VecCreateGhost()` or `VecDuplicate()`)
230: .  insertmode - one of ADD_VALUES, MAX_VALUES, MIN_VALUES or INSERT_VALUES
231: -  scattermode - one of SCATTER_FORWARD or SCATTER_REVERSE

233:    Notes:

235:    Use the following to update the ghost regions with correct values from the owning process
236: .vb
237:        VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
238:        VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
239: .ve

241:    Use the following to accumulate the ghost region values onto the owning processors
242: .vb
243:        VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
244:        VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
245: .ve

247:    To accumulate the ghost region values onto the owning processors and then update
248:    the ghost regions correctly, call the later followed by the former, i.e.,
249: .vb
250:        VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
251:        VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
252:        VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
253:        VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
254: .ve

256:    Level: advanced

258: .seealso: `VecCreateGhost()`, `VecGhostUpdateBegin()`, `VecGhostGetLocalForm()`,
259:           `VecGhostRestoreLocalForm()`, `VecCreateGhostWithArray()`

261: @*/
262: PetscErrorCode VecGhostUpdateEnd(Vec g, InsertMode insertmode, ScatterMode scattermode)
263: {
264:   Vec_MPI  *v;
265:   PetscBool ismpi;

267:   PetscFunctionBegin;
269:   PetscCall(PetscObjectTypeCompare((PetscObject)g, VECMPI, &ismpi));
270:   if (ismpi) {
271:     v = (Vec_MPI *)g->data;
272:     PetscCheck(v->localrep, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Vector is not ghosted");
273:     if (!v->localupdate) PetscFunctionReturn(PETSC_SUCCESS);
274:     if (scattermode == SCATTER_REVERSE) {
275:       PetscCall(VecScatterEnd(v->localupdate, v->localrep, g, insertmode, scattermode));
276:     } else {
277:       PetscCall(VecScatterEnd(v->localupdate, g, v->localrep, insertmode, scattermode));
278:     }
279:   }
280:   PetscFunctionReturn(PETSC_SUCCESS);
281: }