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