Actual source code: ex30f.F90
1: !
2: !
3: ! Tests parallel to parallel scatter where a to from index are
4: ! duplicated
5: program main
6: #include <petsc/finclude/petscvec.h>
7: use petscvec
8: implicit none
10: PetscErrorCode ierr
11: PetscInt nlocal, n, row
12: PetscInt nlocal2,n2,eight
13: PetscMPIInt rank, size
14: PetscInt from(10), to(10)
16: PetscScalar num
17: Vec v1, v2, v3
18: VecScatter scat1, scat2
19: IS fromis, tois
20: n=8
21: nlocal=2
22: PetscCallA(PetscInitialize(ierr))
23: PetscCallMPIA(MPI_COMM_RANK(PETSC_COMM_WORLD,rank,ierr))
24: PetscCallMPIA(MPI_COMM_SIZE(PETSC_COMM_WORLD,size,ierr))
25: if (size.ne.4) then
26: print *, 'Four processor test'
27: stop
28: end if
30: nlocal2 = 2*nlocal
31: n2 = 2*n
32: PetscCallA(VecCreateMPI(PETSC_COMM_WORLD,nlocal2,n2,v1,ierr))
33: PetscCallA(VecCreateMPI(PETSC_COMM_WORLD,nlocal,n,v2,ierr))
34: PetscCallA(VecCreateSeq(PETSC_COMM_SELF,n,v3,ierr))
36: num=2.0
37: row = 1
38: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
39: row = 5
40: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
41: row = 9
42: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
43: row = 13
44: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
45: num=1.0
46: row = 15
47: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
48: row = 3
49: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
50: row = 7
51: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
52: row = 11
53: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
55: PetscCallA(VecAssemblyBegin(v1,ierr))
56: PetscCallA(VecAssemblyEnd(v1,ierr))
58: num=0.0
59: PetscCallA(VecScale(v2,num,ierr))
60: PetscCallA(VecScale(v3,num,ierr))
62: from(1)=1
63: from(2)=5
64: from(3)=9
65: from(4)=13
66: from(5)=3
67: from(6)=7
68: from(7)=11
69: from(8)=15
70: to(1)=0
71: to(2)=0
72: to(3)=0
73: to(4)=0
74: to(5)=7
75: to(6)=7
76: to(7)=7
77: to(8)=7
79: eight = 8
80: PetscCallA(ISCreateGeneral(PETSC_COMM_SELF,eight,from,PETSC_COPY_VALUES,fromis,ierr))
81: PetscCallA(ISCreateGeneral(PETSC_COMM_SELF,eight,to,PETSC_COPY_VALUES,tois,ierr))
82: PetscCallA(VecScatterCreate(v1,fromis,v2,tois,scat1,ierr))
83: PetscCallA(VecScatterCreate(v1,fromis,v3,tois,scat2,ierr))
84: PetscCallA(ISDestroy(fromis,ierr))
85: PetscCallA(ISDestroy(tois,ierr))
87: PetscCallA(VecScatterBegin(scat1,v1,v2,ADD_VALUES,SCATTER_FORWARD,ierr))
88: PetscCallA(VecScatterEnd(scat1,v1,v2,ADD_VALUES,SCATTER_FORWARD,ierr))
90: PetscCallA(VecScatterBegin(scat2,v1,v3,ADD_VALUES,SCATTER_FORWARD,ierr))
91: PetscCallA(VecScatterEnd(scat2,v1,v3,ADD_VALUES,SCATTER_FORWARD,ierr))
93: PetscCallA(PetscObjectSetName(v1, 'V1',ierr))
94: PetscCallA(VecView(v1,PETSC_VIEWER_STDOUT_WORLD,ierr))
96: PetscCallA(PetscObjectSetName(v2, 'V2',ierr))
97: PetscCallA(VecView(v2,PETSC_VIEWER_STDOUT_WORLD,ierr))
99: if (rank.eq.0) then
100: PetscCallA(PetscObjectSetName(v3, 'V3',ierr))
101: PetscCallA(VecView(v3,PETSC_VIEWER_STDOUT_SELF,ierr))
102: end if
104: PetscCallA(VecScatterDestroy(scat1,ierr))
105: PetscCallA(VecScatterDestroy(scat2,ierr))
106: PetscCallA(VecDestroy(v1,ierr))
107: PetscCallA(VecDestroy(v2,ierr))
108: PetscCallA(VecDestroy(v3,ierr))
110: PetscCallA(PetscFinalize(ierr))
112: end
114: !/*TEST
115: !
116: ! test:
117: ! nsize: 4
118: !
119: !TEST*/