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, i1
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: i1 = 1
33: PetscCallA(VecCreateFromOptions(PETSC_COMM_WORLD,PETSC_NULL_CHARACTER,i1,nlocal2,n2,v1,ierr))
34: PetscCallA(VecCreateFromOptions(PETSC_COMM_WORLD,PETSC_NULL_CHARACTER,i1,nlocal,n,v2,ierr))
35: PetscCallA(VecCreateSeq(PETSC_COMM_SELF,n,v3,ierr))
37: num=2.0
38: row = 1
39: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
40: row = 5
41: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
42: row = 9
43: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
44: row = 13
45: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
46: num=1.0
47: row = 15
48: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
49: row = 3
50: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
51: row = 7
52: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
53: row = 11
54: PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
56: PetscCallA(VecAssemblyBegin(v1,ierr))
57: PetscCallA(VecAssemblyEnd(v1,ierr))
59: num=0.0
60: PetscCallA(VecScale(v2,num,ierr))
61: PetscCallA(VecScale(v3,num,ierr))
63: from(1)=1
64: from(2)=5
65: from(3)=9
66: from(4)=13
67: from(5)=3
68: from(6)=7
69: from(7)=11
70: from(8)=15
71: to(1)=0
72: to(2)=0
73: to(3)=0
74: to(4)=0
75: to(5)=7
76: to(6)=7
77: to(7)=7
78: to(8)=7
80: eight = 8
81: PetscCallA(ISCreateGeneral(PETSC_COMM_SELF,eight,from,PETSC_COPY_VALUES,fromis,ierr))
82: PetscCallA(ISCreateGeneral(PETSC_COMM_SELF,eight,to,PETSC_COPY_VALUES,tois,ierr))
83: PetscCallA(VecScatterCreate(v1,fromis,v2,tois,scat1,ierr))
84: PetscCallA(VecScatterCreate(v1,fromis,v3,tois,scat2,ierr))
85: PetscCallA(ISDestroy(fromis,ierr))
86: PetscCallA(ISDestroy(tois,ierr))
88: PetscCallA(VecScatterBegin(scat1,v1,v2,ADD_VALUES,SCATTER_FORWARD,ierr))
89: PetscCallA(VecScatterEnd(scat1,v1,v2,ADD_VALUES,SCATTER_FORWARD,ierr))
91: PetscCallA(VecScatterBegin(scat2,v1,v3,ADD_VALUES,SCATTER_FORWARD,ierr))
92: PetscCallA(VecScatterEnd(scat2,v1,v3,ADD_VALUES,SCATTER_FORWARD,ierr))
94: PetscCallA(PetscObjectSetName(v1, 'V1',ierr))
95: PetscCallA(VecView(v1,PETSC_VIEWER_STDOUT_WORLD,ierr))
97: PetscCallA(PetscObjectSetName(v2, 'V2',ierr))
98: PetscCallA(VecView(v2,PETSC_VIEWER_STDOUT_WORLD,ierr))
100: if (rank.eq.0) then
101: PetscCallA(PetscObjectSetName(v3, 'V3',ierr))
102: PetscCallA(VecView(v3,PETSC_VIEWER_STDOUT_SELF,ierr))
103: end if
105: PetscCallA(VecScatterDestroy(scat1,ierr))
106: PetscCallA(VecScatterDestroy(scat2,ierr))
107: PetscCallA(VecDestroy(v1,ierr))
108: PetscCallA(VecDestroy(v2,ierr))
109: PetscCallA(VecDestroy(v3,ierr))
111: PetscCallA(PetscFinalize(ierr))
113: end
115: !/*TEST
116: !
117: ! test:
118: ! nsize: 4
119: !
120: !TEST*/