Actual source code: ex23.c
2: static char help[] = "Scatters from a parallel vector to a sequential vector.\n\
3: Using a blocked send and a strided receive.\n\n";
5: /*
6: 0 1 2 3 | 4 5 6 7 || 8 9 10 11
8: Scatter first and third block to first processor and
9: second and third block to second processor
10: */
11: #include petscvec.h
15: int main(int argc,char **argv)
16: {
18: PetscInt i,blocks[2],nlocal;
19: PetscMPIInt size,rank;
20: PetscScalar value;
21: Vec x,y;
22: IS is1,is2;
23: VecScatter ctx = 0;
25: PetscInitialize(&argc,&argv,(char*)0,help);
26: MPI_Comm_size(PETSC_COMM_WORLD,&size);
27: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
29: if (size != 2) SETERRQ(1,"Must run with 2 processors");
31: /* create two vectors */
32: if (!rank) nlocal = 8;
33: else nlocal = 4;
34: VecCreate(PETSC_COMM_WORLD,&x);
35: VecSetSizes(x,nlocal,12);
36: VecSetFromOptions(x);
37: VecCreateSeq(PETSC_COMM_SELF,8,&y);
39: /* create two index sets */
40: if (!rank) {
41: blocks[0] = 0; blocks[1] = 8;
42: } else {
43: blocks[0] = 4; blocks[1] = 8;
44: }
45: ISCreateBlock(PETSC_COMM_SELF,4,2,blocks,&is1);
46: ISCreateStride(PETSC_COMM_SELF,8,0,1,&is2);
48: for (i=0; i<12; i++) {
49: value = i;
50: VecSetValues(x,1,&i,&value,INSERT_VALUES);
51: }
52: VecAssemblyBegin(x);
53: VecAssemblyEnd(x);
55: VecScatterCreate(x,is1,y,is2,&ctx);
56: VecScatterBegin(ctx,x,y,INSERT_VALUES,SCATTER_FORWARD);
57: VecScatterEnd(ctx,x,y,INSERT_VALUES,SCATTER_FORWARD);
58: VecScatterDestroy(ctx);
59:
60: PetscSleep(2.0*rank);
61: VecView(y,PETSC_VIEWER_STDOUT_SELF);
63: VecDestroy(x);
64: VecDestroy(y);
65: ISDestroy(is1);
66: ISDestroy(is2);
68: PetscFinalize();
69: return 0;
70: }
71: