Actual source code: ex42.c
2: static char help[] = "Scatters from a parallel vector to a parallel vector.\n\n";
4: #include <petscvec.h>
6: int main(int argc,char **argv)
7: {
8: PetscInt n = 5,N,i;
9: PetscMPIInt size,rank;
10: PetscScalar value,zero = 0.0;
11: Vec x,y;
12: IS is1,is2;
13: VecScatter ctx = 0;
15: PetscInitialize(&argc,&argv,(char*)0,help);
16: MPI_Comm_size(PETSC_COMM_WORLD,&size);
17: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
19: /* create two vectors */
20: N = size*n;
21: VecCreate(PETSC_COMM_WORLD,&y);
22: VecSetSizes(y,n,PETSC_DECIDE);
23: VecSetFromOptions(y);
25: VecCreate(PETSC_COMM_WORLD,&x);
26: VecSetSizes(x,n,PETSC_DECIDE);
27: VecSetFromOptions(x);
29: /* create two index sets */
30: ISCreateStride(PETSC_COMM_WORLD,n,n*rank,1,&is1);
31: ISCreateStride(PETSC_COMM_WORLD,n,(n*(rank+1))%N,1,&is2);
33: /* fill local part of parallel vector x */
34: value = (PetscScalar)(rank+1);
35: for (i=n*rank; i<n*(rank+1); i++) {
36: VecSetValues(x,1,&i,&value,INSERT_VALUES);
37: }
38: VecAssemblyBegin(x);
39: VecAssemblyEnd(x);
41: VecSet(y,zero);
43: VecScatterCreate(x,is1,y,is2,&ctx);
44: for (i=0; i<100; i++) {
45: PetscReal ynorm;
46: PetscInt j;
47: VecNormBegin(y,NORM_2,&ynorm);
48: PetscCommSplitReductionBegin(PetscObjectComm((PetscObject)y));
49: for (j=0; j<3; j++) {
50: VecScatterBegin(ctx,x,y,ADD_VALUES,SCATTER_FORWARD);
51: VecScatterEnd(ctx,x,y,ADD_VALUES,SCATTER_FORWARD);
52: }
53: VecNormEnd(y,NORM_2,&ynorm);
54: /* PetscPrintf(PETSC_COMM_WORLD,"ynorm = %8.2G\n",ynorm); */
55: }
56: VecScatterDestroy(&ctx);
57: VecView(y,PETSC_VIEWER_STDOUT_WORLD);
59: VecDestroy(&x);
60: VecDestroy(&y);
61: ISDestroy(&is1);
62: ISDestroy(&is2);
64: PetscFinalize();
65: return 0;
66: }