Actual source code: ex8.c
petsc-3.12.2 2019-11-22
2: static char help[] = "Demonstrates scattering with strided index sets.\n\n";
4: #include <petscvec.h>
6: int main(int argc,char **argv)
7: {
9: PetscInt n = 6,loc[6] = {0,1,2,3,4,5};
10: PetscScalar two = 2.0,vals[6] = {10,11,12,13,14,15};
11: Vec x,y;
12: IS is1,is2;
13: VecScatter ctx = 0;
15: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
17: /* create two vectors */
18: VecCreateSeq(PETSC_COMM_SELF,n,&x);
19: VecDuplicate(x,&y);
21: /* create two index sets */
22: ISCreateStride(PETSC_COMM_SELF,3,0,2,&is1);
23: ISCreateStride(PETSC_COMM_SELF,3,1,2,&is2);
25: VecSetValues(x,6,loc,vals,INSERT_VALUES);
26: VecView(x,PETSC_VIEWER_STDOUT_SELF);
27: PetscPrintf(PETSC_COMM_SELF,"----\n");
28: VecSet(y,two);
29: VecScatterCreate(x,is1,y,is2,&ctx);
30: VecScatterBegin(ctx,x,y,INSERT_VALUES,SCATTER_FORWARD);
31: VecScatterEnd(ctx,x,y,INSERT_VALUES,SCATTER_FORWARD);
32: VecScatterDestroy(&ctx);
34: VecView(y,PETSC_VIEWER_STDOUT_SELF);
36: ISDestroy(&is1);
37: ISDestroy(&is2);
38: VecDestroy(&x);
39: VecDestroy(&y);
41: PetscFinalize();
42: return ierr;
43: }
47: /*TEST
49: test:
51: TEST*/