Actual source code: ex33.c
1: /*$Id: ex33.c,v 1.1 2001/09/24 21:04:45 balay Exp $*/
3: static char help[] = "Tests the routines VecScatterCreateToAll(), VecScatterCreateToZero()\n\n";
5: #include petscvec.h
9: int main(int argc,char **argv)
10: {
11: PetscInt n = 3,i,len,start,end;
13: PetscMPIInt size,rank;
14: PetscScalar value,*yy;
15: Vec x,y,z,y_t;
16: VecScatter toall,tozero;
18: PetscInitialize(&argc,&argv,(char*)0,help);
19: MPI_Comm_size(PETSC_COMM_WORLD,&size);
20: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
22: /* create two vectors */
23: VecCreateMPI(PETSC_COMM_WORLD,PETSC_DECIDE,size*n,&x);
25: /* each processor inserts its values */
27: VecGetOwnershipRange(x,&start,&end);
28: for (i=start; i<end; i++) {
29: value = (PetscScalar) i;
30: VecSetValues(x,1,&i,&value,INSERT_VALUES);
31: }
32: VecAssemblyBegin(x);
33: VecAssemblyEnd(x);
34: VecView(x,PETSC_VIEWER_STDOUT_WORLD);
36: VecScatterCreateToAll(x,&toall,&y);
37: VecScatterBegin(toall,x,y,INSERT_VALUES,SCATTER_FORWARD);
38: VecScatterEnd(toall,x,y,INSERT_VALUES,SCATTER_FORWARD);
39: VecScatterDestroy(toall);
41: /* Cannot view the above vector with VecView(), so place it in an MPI Vec
42: and do a VecView() */
43: VecGetArray(y,&yy);
44: VecGetLocalSize(y,&len);
45: VecCreateMPIWithArray(PETSC_COMM_WORLD,len,PETSC_DECIDE,yy,&y_t);
46: VecView(y_t,PETSC_VIEWER_STDOUT_WORLD);
47: VecDestroy(y_t);
48: VecRestoreArray(y,&yy);
50: VecScatterCreateToAll(x,&tozero,&z);
51: VecScatterBegin(tozero,x,z,INSERT_VALUES,SCATTER_FORWARD);
52: VecScatterEnd(tozero,x,z,INSERT_VALUES,SCATTER_FORWARD);
53: VecScatterDestroy(tozero);
54: if (!rank) {
55: VecView(z,PETSC_VIEWER_STDOUT_SELF);
56: }
57: VecDestroy(z);
59: VecScatterCreateToZero(x,&tozero,&z);
60: VecScatterBegin(tozero,x,z,INSERT_VALUES,SCATTER_FORWARD);
61: VecScatterEnd(tozero,x,z,INSERT_VALUES,SCATTER_FORWARD);
62: VecScatterDestroy(tozero);
63: VecDestroy(z);
65: VecDestroy(x);
66: VecDestroy(y);
68: PetscFinalize();
69: return 0;
70: }
71: