Actual source code: ex11.c
2: static char help[]= "Scatters between parallel vectors. \n\
3: uses block index sets\n\n";
5: #include <petscvec.h>
7: int main(int argc,char **argv)
8: {
9: PetscInt bs=1,n=5,N,i,low;
10: PetscInt ix0[3] = {5,7,9},iy0[3] = {1,2,4},ix1[3] = {2,3,1},iy1[3] = {0,3,9};
11: PetscMPIInt size,rank;
12: PetscScalar *array;
13: Vec x,x1,y;
14: IS isx,isy;
15: VecScatter ctx;
16: VecScatterType type;
17: PetscBool flg;
19: PetscInitialize(&argc,&argv,(char*)0,help);
20: MPI_Comm_size(PETSC_COMM_WORLD,&size);
21: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
25: PetscOptionsGetInt(NULL,NULL,"-bs",&bs,NULL);
26: n = bs*n;
28: /* Create vector x over shared memory */
29: VecCreate(PETSC_COMM_WORLD,&x);
30: VecSetSizes(x,n,PETSC_DECIDE);
31: VecSetFromOptions(x);
33: VecGetOwnershipRange(x,&low,NULL);
34: VecGetArray(x,&array);
35: for (i=0; i<n; i++) {
36: array[i] = (PetscScalar)(i + low);
37: }
38: VecRestoreArray(x,&array);
39: /* VecView(x,PETSC_VIEWER_STDOUT_WORLD); */
41: /* Test some vector functions */
42: VecAssemblyBegin(x);
43: VecAssemblyEnd(x);
45: VecGetSize(x,&N);
46: VecGetLocalSize(x,&n);
48: VecDuplicate(x,&x1);
49: VecCopy(x,x1);
50: VecEqual(x,x1,&flg);
53: VecScale(x1,2.0);
54: VecSet(x1,10.0);
55: /* VecView(x1,PETSC_VIEWER_STDOUT_WORLD); */
57: /* Create vector y over shared memory */
58: VecCreate(PETSC_COMM_WORLD,&y);
59: VecSetSizes(y,n,PETSC_DECIDE);
60: VecSetFromOptions(y);
61: VecGetArray(y,&array);
62: for (i=0; i<n; i++) {
63: array[i] = -(PetscScalar) (i + 100*rank);
64: }
65: VecRestoreArray(y,&array);
66: VecAssemblyBegin(y);
67: VecAssemblyEnd(y);
68: /* VecView(y,PETSC_VIEWER_STDOUT_WORLD); */
70: /* Create two index sets */
71: if (rank == 0) {
72: ISCreateBlock(PETSC_COMM_SELF,bs,3,ix0,PETSC_COPY_VALUES,&isx);
73: ISCreateBlock(PETSC_COMM_SELF,bs,3,iy0,PETSC_COPY_VALUES,&isy);
74: } else {
75: ISCreateBlock(PETSC_COMM_SELF,bs,3,ix1,PETSC_COPY_VALUES,&isx);
76: ISCreateBlock(PETSC_COMM_SELF,bs,3,iy1,PETSC_COPY_VALUES,&isy);
77: }
79: if (rank == 10) {
80: PetscPrintf(PETSC_COMM_SELF,"\n[%d] isx:\n",rank);
81: ISView(isx,PETSC_VIEWER_STDOUT_SELF);
82: PetscPrintf(PETSC_COMM_SELF,"\n[%d] isy:\n",rank);
83: ISView(isy,PETSC_VIEWER_STDOUT_SELF);
84: }
86: /* Create Vector scatter */
87: VecScatterCreate(x,isx,y,isy,&ctx);
88: VecScatterSetFromOptions(ctx);
89: VecScatterGetType(ctx,&type);
90: PetscPrintf(PETSC_COMM_WORLD,"scatter type %s\n",type);
92: /* Test forward vecscatter */
93: VecSet(y,0.0);
94: VecScatterBegin(ctx,x,y,ADD_VALUES,SCATTER_FORWARD);
95: VecScatterEnd(ctx,x,y,ADD_VALUES,SCATTER_FORWARD);
96: PetscPrintf(PETSC_COMM_WORLD,"\nSCATTER_FORWARD y:\n");
97: VecView(y,PETSC_VIEWER_STDOUT_WORLD);
99: /* Test reverse vecscatter */
100: VecSet(x,0.0);
101: VecSet(y,0.0);
102: VecGetOwnershipRange(y,&low,NULL);
103: VecGetArray(y,&array);
104: for (i=0; i<n; i++) {
105: array[i] = (PetscScalar)(i+ low);
106: }
107: VecRestoreArray(y,&array);
108: VecScatterBegin(ctx,y,x,ADD_VALUES,SCATTER_REVERSE);
109: VecScatterEnd(ctx,y,x,ADD_VALUES,SCATTER_REVERSE);
110: PetscPrintf(PETSC_COMM_WORLD,"\nSCATTER_REVERSE x:\n");
111: VecView(x,PETSC_VIEWER_STDOUT_WORLD);
113: /* Free objects */
114: VecScatterDestroy(&ctx);
115: ISDestroy(&isx);
116: ISDestroy(&isy);
117: VecDestroy(&x);
118: VecDestroy(&x1);
119: VecDestroy(&y);
120: PetscFinalize();
121: return 0;
122: }
124: /*TEST
126: test:
127: nsize: 2
129: TEST*/