Actual source code: ex4.c
petsc-3.8.3 2017-12-09
2: static char help[] = "Scatters from a parallel vector into seqential vectors.\n\n";
4: #include <petscvec.h>
6: int main(int argc,char **argv)
7: {
9: PetscMPIInt rank;
10: PetscInt n = 5,idx1[2] = {0,3},idx2[2] = {1,4};
11: PetscScalar one = 1.0,two = 2.0;
12: Vec x,y;
13: IS is1,is2;
14: VecScatter ctx = 0;
16: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
17: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
18: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
20: /* create two vectors */
21: VecCreate(PETSC_COMM_WORLD,&x);
22: VecSetSizes(x,n,PETSC_DECIDE);
23: VecSetFromOptions(x);
24: VecCreateSeq(PETSC_COMM_SELF,n,&y);
26: /* create two index sets */
27: ISCreateGeneral(PETSC_COMM_SELF,2,idx1,PETSC_COPY_VALUES,&is1);
28: ISCreateGeneral(PETSC_COMM_SELF,2,idx2,PETSC_COPY_VALUES,&is2);
30: VecSet(x,one);
31: VecSet(y,two);
32: VecScatterCreate(x,is1,y,is2,&ctx);
33: VecScatterBegin(ctx,x,y,INSERT_VALUES,SCATTER_FORWARD);
34: VecScatterEnd(ctx,x,y,INSERT_VALUES,SCATTER_FORWARD);
35: VecScatterDestroy(&ctx);
37: if (!rank) {VecView(y,PETSC_VIEWER_STDOUT_SELF);}
39: ISDestroy(&is1);
40: ISDestroy(&is2);
42: VecDestroy(&x);
43: VecDestroy(&y);
44: PetscFinalize();
45: return ierr;
46: }