Actual source code: ex7.c
1: static char help[]= "Test vecscatter of different block sizes across processes\n\n";
3: #include <petscvec.h>
4: int main(int argc,char **argv)
5: {
6: PetscInt i,bs,n,low,high;
7: PetscMPIInt nproc,rank;
8: Vec x,y,z;
9: IS ix,iy;
10: VecScatter vscat;
11: const PetscScalar *yv;
13: PetscInitialize(&argc,&argv,(char*)0,help);
14: MPI_Comm_size(PETSC_COMM_WORLD,&nproc);
15: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
18: /* Create an MPI vector x of size 12 on two processes, and set x = {0, 1, 2, .., 11} */
19: VecCreateMPI(PETSC_COMM_WORLD,6,PETSC_DECIDE,&x);
20: VecGetOwnershipRange(x,&low,&high);
21: for (i=low; i<high; i++) VecSetValue(x,i,(PetscScalar)i,INSERT_VALUES);
22: VecAssemblyBegin(x);
23: VecAssemblyEnd(x);
25: /* Create a seq vector y, and a parallel to sequential (PtoS) vecscatter to scatter x to y */
26: if (rank == 0) {
27: /* On rank 0, seq y is of size 6. We will scatter x[0,1,2,6,7,8] to y[0,1,2,3,4,5] using IS with bs=3 */
28: PetscInt idx[2]={0,2};
29: PetscInt idy[2]={0,1};
30: n = 6;
31: bs = 3;
32: VecCreateSeq(PETSC_COMM_SELF,n,&y);
33: ISCreateBlock(PETSC_COMM_SELF,bs,2,idx,PETSC_COPY_VALUES,&ix);
34: ISCreateBlock(PETSC_COMM_SELF,bs,2,idy,PETSC_COPY_VALUES,&iy);
35: } else {
36: /* On rank 1, seq y is of size 4. We will scatter x[4,5,10,11] to y[0,1,2,3] using IS with bs=2 */
37: PetscInt idx[2]= {2,5};
38: PetscInt idy[2]= {0,1};
39: n = 4;
40: bs = 2;
41: VecCreateSeq(PETSC_COMM_SELF,n,&y);
42: ISCreateBlock(PETSC_COMM_SELF,bs,2,idx,PETSC_COPY_VALUES,&ix);
43: ISCreateBlock(PETSC_COMM_SELF,bs,2,idy,PETSC_COPY_VALUES,&iy);
44: }
45: VecScatterCreate(x,ix,y,iy,&vscat);
47: /* Do the vecscatter */
48: VecScatterBegin(vscat,x,y,INSERT_VALUES,SCATTER_FORWARD);
49: VecScatterEnd(vscat,x,y,INSERT_VALUES,SCATTER_FORWARD);
51: /* Print y. Since y is sequential, we put y in a parallel z to print its value on both ranks */
52: VecGetArrayRead(y,&yv);
53: VecCreateMPIWithArray(PETSC_COMM_WORLD,1,n,PETSC_DECIDE,yv,&z);
54: VecView(z,PETSC_VIEWER_STDOUT_WORLD);
55: VecRestoreArrayRead(y,&yv);
57: ISDestroy(&ix);
58: ISDestroy(&iy);
59: VecDestroy(&x);
60: VecDestroy(&y);
61: VecDestroy(&z);
62: VecScatterDestroy(&vscat);
64: PetscFinalize();
65: return 0;
66: }
68: /*TEST
70: test:
71: nsize: 2
72: args:
73: requires:
74: TEST*/