Actual source code: vecmpitoseq.c

  1: #define PETSCVEC_DLL

 3:  #include private/vecimpl.h

  7: /*@
  8:       VecScatterCreateToAll - Creates a vector and a scatter context that copies all 
  9:           vector values to each processor

 11:   Collective

 13:   Input Parameter: 
 14: .  vin  - input MPIVEC

 16:   Output Parameter:
 17: +  ctx - scatter context
 18: -  vout - output SEQVEC that is large enough to scatter into

 20:   Level: intermediate

 22:    Note: vout may be PETSC_NULL [PETSC_NULL_OBJECT from fortran] if you do not 
 23:    need to have it created

 25:    Usage:
 26: $        VecScatterCreateToAll(vin,&ctx,&vout);
 27: $
 28: $        // scatter as many times as you need 
 29: $        VecScatterBegin(ctx,vin,vout,INSERT_VALUES,SCATTER_FORWARD);
 30: $        VecScatterEnd(ctx,vin,vout,INSERT_VALUES,SCATTER_FORWARD);
 31: $
 32: $        // destroy scatter context and local vector when no longer needed
 33: $        VecScatterDestroy(ctx);
 34: $        VecDestroy(vout);

 36:     Do NOT create a vector and then pass it in as the final argument vout! vout is created by this routine
 37:   automatically (unless you pass PETSC_NULL in for that argument if you do not need it).

 39: .seealso VecScatterCreate(), VecScatterCreateToZero(), VecScatterBegin(), VecScatterEnd()

 41: @*/
 42: PetscErrorCode  VecScatterCreateToAll(Vec vin,VecScatter *ctx,Vec *vout)
 43: {

 46:   PetscInt       N;
 47:   IS             is;
 48:   Vec            tmp;
 49:   Vec           *tmpv;
 50:   PetscTruth     tmpvout = PETSC_FALSE;

 56:   if (vout) {
 58:     tmpv = vout;
 59:   } else {
 60:     tmpvout = PETSC_TRUE;
 61:     tmpv = &tmp;
 62:   }

 64:   /* Create seq vec on each proc, with the same size of the original mpi vec */
 65:   VecGetSize(vin,&N);
 66:   VecCreateSeq(PETSC_COMM_SELF,N,tmpv);
 67:   /* Create the VecScatter ctx with the communication info */
 68:   ISCreateStride(PETSC_COMM_SELF,N,0,1,&is);
 69:   VecScatterCreate(vin,is,*tmpv,is,ctx);
 70:   ISDestroy(is);
 71:   if (tmpvout) {VecDestroy(*tmpv);}
 72:   return(0);
 73: }


 78: /*@
 79:       VecScatterCreateToZero - Creates an output vector and a scatter context used to 
 80:               copy all vector values into the output vector on the zeroth processor

 82:   Collective

 84:   Input Parameter: 
 85: .  vin  - input MPIVEC

 87:   Output Parameter:
 88: +  ctx - scatter context
 89: -  vout - output SEQVEC that is large enough to scatter into on processor 0 and
 90:           of length zero on all other processors

 92:   Level: intermediate

 94:    Note: vout may be PETSC_NULL [PETSC_NULL_OBJECT from fortran] if you do not 
 95:    need to have it created

 97:    Usage:
 98: $        VecScatterCreateToZero(vin,&ctx,&vout);
 99: $
100: $        // scatter as many times as you need 
101: $        VecScatterBegin(ctx,vin,vout,INSERT_VALUES,SCATTER_FORWARD);
102: $        VecScatterEnd(ctx,vin,vout,INSERT_VALUES,SCATTER_FORWARD);
103: $
104: $        // destroy scatter context and local vector when no longer needed
105: $        VecScatterDestroy(ctx);
106: $        VecDestroy(vout);

108: .seealso VecScatterCreate(), VecScatterCreateToAll(), VecScatterBegin(), VecScatterEnd()

110:     Do NOT create a vector and then pass it in as the final argument vout! vout is created by this routine
111:   automatically (unless you pass PETSC_NULL in for that argument if you do not need it).

113: @*/
114: PetscErrorCode  VecScatterCreateToZero(Vec vin,VecScatter *ctx,Vec *vout)
115: {

118:   PetscInt       N;
119:   PetscMPIInt    rank;
120:   IS             is;
121:   Vec            tmp;
122:   Vec           *tmpv;
123:   PetscTruth     tmpvout = PETSC_FALSE;

129:   if (vout) {
131:     tmpv = vout;
132:   } else {
133:     tmpvout = PETSC_TRUE;
134:     tmpv = &tmp;
135:   }

137:   /* Create vec on each proc, with the same size of the original mpi vec (all on process 0)*/
138:   VecGetSize(vin,&N);
139:   MPI_Comm_rank(((PetscObject)vin)->comm,&rank);
140:   if (rank) N = 0;
141:   VecCreateSeq(PETSC_COMM_SELF,N,tmpv);
142:   /* Create the VecScatter ctx with the communication info */
143:   ISCreateStride(PETSC_COMM_SELF,N,0,1,&is);
144:   VecScatterCreate(vin,is,*tmpv,is,ctx);
145:   ISDestroy(is);
146:   if (tmpvout) {VecDestroy(*tmpv);}
147:   return(0);
148: }