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