Actual source code: ex14.c
petsc-3.10.5 2019-03-28
2: static char help[] = "Scatters from a sequential vector to a parallel vector.\n\
3: This does the tricky case.\n\n";
5: /*T
6: requires: x
7: T*/
9: #include <petscvec.h>
11: int main(int argc,char **argv)
12: {
14: PetscInt n = 5,N;
15: PetscMPIInt size,rank;
16: PetscScalar value,zero = 0.0;
17: Vec x,y;
18: IS is1,is2;
19: VecScatter ctx = 0;
21: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
22: MPI_Comm_size(PETSC_COMM_WORLD,&size);
23: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
25: /* create two vectors */
26: N = size*n;
27: VecCreate(PETSC_COMM_WORLD,&y);
28: VecSetSizes(y,PETSC_DECIDE,N);
29: VecSetFromOptions(y);
30: VecCreateSeq(PETSC_COMM_SELF,N,&x);
32: /* create two index sets */
33: ISCreateStride(PETSC_COMM_SELF,n,0,1,&is1);
34: ISCreateStride(PETSC_COMM_SELF,n,rank,1,&is2);
36: value = rank+1;
37: VecSet(x,value);
38: VecSet(y,zero);
40: VecScatterCreate(x,is1,y,is2,&ctx);
41: VecScatterBegin(ctx,x,y,ADD_VALUES,SCATTER_FORWARD);
42: VecScatterEnd(ctx,x,y,ADD_VALUES,SCATTER_FORWARD);
43: VecScatterDestroy(&ctx);
45: VecView(y,PETSC_VIEWER_STDOUT_WORLD);
47: VecDestroy(&x);
48: VecDestroy(&y);
49: ISDestroy(&is1);
50: ISDestroy(&is2);
52: PetscFinalize();
53: return ierr;
54: }
58: /*TEST
60: test:
61: nsize: 2
63: TEST*/