1
2 static char help[] = "Tests the routines VecScatterCreateToAll(), VecScatterCreateToZero()\n\n";
3
4 #include <petscvec.h>
5
main(int argc,char ** argv)6 int main(int argc,char **argv)
7 {
8 PetscInt n = 3,i,len,start,end;
9 PetscErrorCode ierr;
10 PetscMPIInt size,rank;
11 PetscScalar value,*yy;
12 Vec x,y,z,y_t;
13 VecScatter toall,tozero;
14
15 ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
16 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
17 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
18
19 /* create two vectors */
20 ierr = VecCreateMPI(PETSC_COMM_WORLD,PETSC_DECIDE,size*n,&x);CHKERRQ(ierr);
21
22 /* each processor inserts its values */
23
24 ierr = VecGetOwnershipRange(x,&start,&end);CHKERRQ(ierr);
25 for (i=start; i<end; i++) {
26 value = (PetscScalar) i;
27 ierr = VecSetValues(x,1,&i,&value,INSERT_VALUES);CHKERRQ(ierr);
28 }
29 ierr = VecAssemblyBegin(x);CHKERRQ(ierr);
30 ierr = VecAssemblyEnd(x);CHKERRQ(ierr);
31 ierr = VecView(x,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
32
33 ierr = VecScatterCreateToAll(x,&toall,&y);CHKERRQ(ierr);
34 ierr = VecScatterBegin(toall,x,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
35 ierr = VecScatterEnd(toall,x,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
36 ierr = VecScatterDestroy(&toall);CHKERRQ(ierr);
37
38 /* Cannot view the above vector with VecView(), so place it in an MPI Vec
39 and do a VecView() */
40 ierr = VecGetArray(y,&yy);CHKERRQ(ierr);
41 ierr = VecGetLocalSize(y,&len);CHKERRQ(ierr);
42 ierr = VecCreateMPIWithArray(PETSC_COMM_WORLD,1,len,PETSC_DECIDE,yy,&y_t);CHKERRQ(ierr);
43 ierr = VecView(y_t,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
44 ierr = VecDestroy(&y_t);CHKERRQ(ierr);
45 ierr = VecRestoreArray(y,&yy);CHKERRQ(ierr);
46
47 ierr = VecScatterCreateToAll(x,&tozero,&z);CHKERRQ(ierr);
48 ierr = VecScatterBegin(tozero,x,z,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
49 ierr = VecScatterEnd(tozero,x,z,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
50 ierr = VecScatterDestroy(&tozero);CHKERRQ(ierr);
51 if (!rank) {
52 ierr = VecView(z,PETSC_VIEWER_STDOUT_SELF);CHKERRQ(ierr);
53 }
54 ierr = VecDestroy(&z);CHKERRQ(ierr);
55
56 ierr = VecScatterCreateToZero(x,&tozero,&z);CHKERRQ(ierr);
57 ierr = VecScatterBegin(tozero,x,z,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
58 ierr = VecScatterEnd(tozero,x,z,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
59 ierr = VecScatterDestroy(&tozero);CHKERRQ(ierr);
60 ierr = VecDestroy(&z);CHKERRQ(ierr);
61
62 ierr = VecDestroy(&x);CHKERRQ(ierr);
63 ierr = VecDestroy(&y);CHKERRQ(ierr);
64
65 ierr = PetscFinalize();
66 return ierr;
67 }
68
69
70
71 /*TEST
72
73 test:
74 nsize: 4
75
76 TEST*/
77