Actual source code: commonmpvec.c

  1: #include <../src/vec/vec/impls/mpi/pvecimpl.h>

  3: /*
  4:   This is used in VecGhostGetLocalForm and VecGhostRestoreLocalForm to ensure
  5:   that the state is updated if either vector has changed since the last time
  6:   one of these functions was called.  It could apply to any PetscObject, but
  7:   VecGhost is quite different from other objects in that two separate vectors
  8:   look at the same memory.

 10:   In principle, we could only propagate state to the local vector on
 11:   GetLocalForm and to the global vector on RestoreLocalForm, but this version is
 12:   more conservative (i.e. robust against misuse) and simpler.

 14:   Note that this function is correct and changes nothing if both arguments are the
 15:   same, which is the case in serial.
 16: */
 17: static PetscErrorCode VecGhostStateSync_Private(Vec g, Vec l)
 18: {
 19:   PetscObjectState gstate, lstate;

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

 29: /*@
 30:   VecGhostGetLocalForm - Obtains the local ghosted representation of
 31:   a parallel vector (obtained with `VecCreateGhost()`, `VecCreateGhostWithArray()` or `VecCreateSeq()`).

 33:   Logically Collective

 35:   Input Parameter:
 36: . g - the global vector

 38:   Output Parameter:
 39: . l - the local (ghosted) representation,`NULL` if `g` is not ghosted

 41:   Level: advanced

 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: .vb
 53:      VecGhostUpdateBegin(x, INSERT_VALUES, SCATTER_FORWARD);
 54:      VecGhostUpdateEnd(x, INSERT_VALUES, SCATTER_FORWARD);
 55:      VecGhostGetLocalForm(x, &xlocal);
 56:      VecGetArrayRead(xlocal, &xvalues);
 57:         // access the non-ghost values in locations xvalues[0:n-1] and ghost values in locations xvalues[n:n+nghost];
 58:      VecRestoreArrayRead(xlocal, &xvalues);
 59:      VecGhostRestoreLocalForm(x, &xlocal);
 60: .ve

 62:   One should call `VecGhostRestoreLocalForm()` or `VecDestroy()` once one is
 63:   finished using the object.

 65: .seealso: [](ch_vectors), `VecGhostUpdateBegin()`, `VecGhostUpdateEnd()`, `Vec`, `VecType`, `VecCreateGhost()`, `VecGhostRestoreLocalForm()`, `VecCreateGhostWithArray()`
 66: @*/
 67: PetscErrorCode VecGhostGetLocalForm(Vec g, Vec *l)
 68: {
 69:   PetscBool isseq, ismpi;

 71:   PetscFunctionBegin;
 73:   PetscAssertPointer(l, 2);

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

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

 95:   Not Collective

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

101:   Output Parameter:
102: . flg - `PETSC_TRUE` if `l` is the local form

104:   Level: advanced

106: .seealso: [](ch_vectors), `Vec`, `VecType`, `VecCreateGhost()`, `VecGhostRestoreLocalForm()`, `VecCreateGhostWithArray()`, `VecGhostGetLocalForm()`
107: @*/
108: PetscErrorCode VecGhostIsLocalForm(Vec g, Vec l, PetscBool *flg)
109: {
110:   PetscBool isseq, ismpi;

112:   PetscFunctionBegin;

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

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

132:   Logically Collective

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

138:   Level: advanced

140: .seealso: [](ch_vectors), `VecGhostUpdateBegin()`, `VecGhostUpdateEnd()`, `Vec`, `VecType`, `VecCreateGhost()`, `VecGhostGetLocalForm()`, `VecCreateGhostWithArray()`
141: @*/
142: PetscErrorCode VecGhostRestoreLocalForm(Vec g, Vec *l)
143: {
144:   PetscFunctionBegin;
145:   if (*l) {
146:     PetscCall(VecGhostStateSync_Private(g, *l));
147:     PetscCall(PetscObjectDereference((PetscObject)*l));
148:   }
149:   PetscFunctionReturn(PETSC_SUCCESS);
150: }

152: /*@
153:   VecGhostUpdateBegin - Begins the vector scatter to update the vector from
154:   local representation to global or global representation to local.

156:   Neighbor-wise Collective

158:   Input Parameters:
159: + g           - the vector (obtained with `VecCreateGhost()` or `VecDuplicate()`)
160: . insertmode  - one of `ADD_VALUES`, `MAX_VALUES`, `MIN_VALUES` or `INSERT_VALUES`
161: - scattermode - one of `SCATTER_FORWARD` (update ghosts) or `SCATTER_REVERSE` (update local values from ghosts)

163:   Level: advanced

165:   Notes:
166:   Use the following to update the ghost regions with correct values from the owning process
167: .vb
168:        VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
169:        VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
170: .ve

172:   Use the following to accumulate the ghost region values onto the owning processors
173: .vb
174:        VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
175:        VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
176: .ve

178:   To accumulate the ghost region values onto the owning processors and then update
179:   the ghost regions correctly, call the latter followed by the former, i.e.,
180: .vb
181:        VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
182:        VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
183:        VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
184:        VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
185: .ve

187: .seealso: [](ch_vectors), `Vec`, `VecType`, `VecCreateGhost()`, `VecGhostUpdateEnd()`, `VecGhostGetLocalForm()`,
188:           `VecGhostRestoreLocalForm()`, `VecCreateGhostWithArray()`
189: @*/
190: PetscErrorCode VecGhostUpdateBegin(Vec g, InsertMode insertmode, ScatterMode scattermode)
191: {
192:   Vec_MPI  *v;
193:   PetscBool ismpi, isseq;

195:   PetscFunctionBegin;
197:   PetscCall(PetscObjectTypeCompare((PetscObject)g, VECMPI, &ismpi));
198:   PetscCall(PetscObjectTypeCompare((PetscObject)g, VECSEQ, &isseq));
199:   if (ismpi) {
200:     v = (Vec_MPI *)g->data;
201:     PetscCheck(v->localrep, PetscObjectComm((PetscObject)g), PETSC_ERR_ARG_WRONG, "Vector is not ghosted");
202:     if (!v->localupdate) PetscFunctionReturn(PETSC_SUCCESS);
203:     if (scattermode == SCATTER_REVERSE) {
204:       PetscCall(VecScatterBegin(v->localupdate, v->localrep, g, insertmode, scattermode));
205:     } else {
206:       PetscCall(VecScatterBegin(v->localupdate, g, v->localrep, insertmode, scattermode));
207:     }
208:   } else if (isseq) {
209:     /* Do nothing */
210:   } else SETERRQ(PetscObjectComm((PetscObject)g), PETSC_ERR_ARG_WRONG, "Vector is not ghosted");
211:   PetscFunctionReturn(PETSC_SUCCESS);
212: }

214: /*@
215:   VecGhostUpdateEnd - End the vector scatter to update the vector from
216:   local representation to global or global representation to local.

218:   Neighbor-wise Collective

220:   Input Parameters:
221: + g           - the vector (obtained with `VecCreateGhost()` or `VecDuplicate()`)
222: . insertmode  - one of `ADD_VALUES`, `MAX_VALUES`, `MIN_VALUES` or `INSERT_VALUES`
223: - scattermode - one of `SCATTER_FORWARD` (update ghosts) or `SCATTER_REVERSE` (update local values from ghosts)

225:   Level: advanced

227:   Notes:
228:   Use the following to update the ghost regions with correct values from the owning process
229: .vb
230:        VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
231:        VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
232: .ve

234:   Use the following to accumulate the ghost region values onto the owning processors
235: .vb
236:        VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
237:        VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
238: .ve

240:   To accumulate the ghost region values onto the owning processors and then update
241:   the ghost regions correctly, call the later followed by the former, i.e.,
242: .vb
243:        VecGhostUpdateBegin(v,ADD_VALUES,SCATTER_REVERSE);
244:        VecGhostUpdateEnd(v,ADD_VALUES,SCATTER_REVERSE);
245:        VecGhostUpdateBegin(v,INSERT_VALUES,SCATTER_FORWARD);
246:        VecGhostUpdateEnd(v,INSERT_VALUES,SCATTER_FORWARD);
247: .ve

249: .seealso: [](ch_vectors), `Vec`, `VecType`, `VecCreateGhost()`, `VecGhostUpdateBegin()`, `VecGhostGetLocalForm()`,
250:           `VecGhostRestoreLocalForm()`, `VecCreateGhostWithArray()`
251: @*/
252: PetscErrorCode VecGhostUpdateEnd(Vec g, InsertMode insertmode, ScatterMode scattermode)
253: {
254:   Vec_MPI  *v;
255:   PetscBool ismpi;

257:   PetscFunctionBegin;
259:   PetscCall(PetscObjectTypeCompare((PetscObject)g, VECMPI, &ismpi));
260:   if (ismpi) {
261:     v = (Vec_MPI *)g->data;
262:     PetscCheck(v->localrep, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Vector is not ghosted");
263:     if (!v->localupdate) PetscFunctionReturn(PETSC_SUCCESS);
264:     if (scattermode == SCATTER_REVERSE) {
265:       PetscCall(VecScatterEnd(v->localupdate, v->localrep, g, insertmode, scattermode));
266:     } else {
267:       PetscCall(VecScatterEnd(v->localupdate, g, v->localrep, insertmode, scattermode));
268:     }
269:   }
270:   PetscFunctionReturn(PETSC_SUCCESS);
271: }