Actual source code: ex17f.F
1: !
2: !
3: ! "Scatters from a parallel vector to a sequential vector. In
4: ! this case each local vector is as long as the entire parallel vector.
5: !
6: implicit none
7: #include "finclude/petsc.h"
8: #include "finclude/petscis.h"
9: #include "finclude/petscvec.h"
10: #include "finclude/petscsys.h"
12: PetscErrorCode ierr
13: PetscMPIInt size,rank
14: PetscInt n,NN,low,high
15: PetscInt iglobal,i,ione
16: PetscInt first,stride
17: PetscScalar value,zero
18: Vec x,y
19: IS is1,is2
20: VecScatter ctx
22: n = 5
23: zero = 0.d0
24: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
26: call MPI_Comm_size(MPI_COMM_WORLD,size,ierr)
27: call MPI_Comm_rank(MPI_COMM_WORLD,rank,ierr)
29: ! create two vectors
30: ! one parallel and one sequential. The sequential one on each processor
31: ! is as long as the entire parallel one.
33: NN = size*n
35: call VecCreateMPI(MPI_COMM_WORLD,PETSC_DECIDE,NN,y,ierr)
36: call VecCreateSeq(MPI_COMM_SELF,NN,x,ierr)
38: call VecSet(x,zero,ierr)
39: call VecGetOwnershipRange(y,low,high,ierr)
40: ione = 1
41: do 10, i=0,n-1
42: iglobal = i + low
43: value = i + 10*rank
44: call VecSetValues(y,ione,iglobal,value,INSERT_VALUES,ierr)
45: 10 continue
47: call VecAssemblyBegin(y,ierr)
48: call VecAssemblyEnd(y,ierr)
49: !
50: ! View the parallel vector
51: !
52: call VecView(y,PETSC_VIEWER_STDOUT_WORLD,ierr)
54: ! create two index sets and the scatter context to move the contents of
55: ! of the parallel vector to each sequential vector. If you want the
56: ! parallel vector delivered to only one processor then create a is2
57: ! of length zero on all processors except the one to receive the parallel vector
59: first = 0
60: stride = 1
61: call ISCreateStride(PETSC_COMM_SELF,NN,first,stride,is1,ierr)
62: call ISCreateStride(PETSC_COMM_SELF,NN,first,stride,is2,ierr)
63: call VecScatterCreate(y,is2,x,is1,ctx,ierr)
64: call VecScatterBegin(ctx,y,x,ADD_VALUES,SCATTER_FORWARD,ierr)
65: call VecScatterEnd(ctx,y,x,ADD_VALUES,SCATTER_FORWARD,ierr)
66: call VecScatterDestroy(ctx,ierr)
67: !
68: ! View the sequential vector on the 0th processor
69: !
70: if (rank .eq. 0) then
71: call VecView(x,PETSC_VIEWER_STDOUT_SELF,ierr)
72: endif
74: call VecDestroy(x,ierr)
75: call VecDestroy(y,ierr)
76: call ISDestroy(is1,ierr)
77: call ISDestroy(is2,ierr)
79: call PetscFinalize(ierr)
80: end
81: