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