Actual source code: ex30f.F

  1: !
  2: !
  3: !  Tests parallel to parallel scatter where a to from index are
  4: !  duplicated
  5:       program main
  6:       implicit none

 8:  #include include/finclude/petsc.h
 9:  #include include/finclude/petscis.h
 10:  #include include/finclude/petscvec.h

 12:       PetscErrorCode ierr
 13:       PetscInt  nlocal, n, row, nlocal2,n2,eight
 14:       PetscMPIInt rank, size
 15:       PetscInt from(10), to(10)

 17:       PetscScalar num
 18:       Vec v1, v2, v3
 19:       VecScatter scat1, scat2
 20:       IS fromis, tois
 21:       n=8
 22:       nlocal=2
 23:       call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
 24:       call MPI_COMM_RANK(PETSC_COMM_WORLD,rank,ierr)
 25:       call MPI_COMM_SIZE(PETSC_COMM_WORLD,size,ierr)
 26:       if (size.ne.4) then
 27:          print *, 'Four processor test'
 28:          stop
 29:       end if
 30: 
 31:       nlocal2 = 2*nlocal
 32:       n2      = 2*n
 33:       call VecCreateMPI(PETSC_COMM_WORLD,nlocal2,n2,v1,ierr)
 34:       call VecCreateMPI(PETSC_COMM_WORLD,nlocal,n,v2,ierr)
 35:       call VecCreateSeq(PETSC_COMM_SELF,n,v3,ierr)

 37:       num=2.0
 38:       row = 1
 39:       call VecSetValue(v1,row,num,INSERT_VALUES,ierr)
 40:       row = 5
 41:       call VecSetValue(v1,row,num,INSERT_VALUES,ierr)
 42:       row = 9
 43:       call VecSetValue(v1,row,num,INSERT_VALUES,ierr)
 44:       row = 13
 45:       call VecSetValue(v1,row,num,INSERT_VALUES,ierr)
 46:       num=1.0
 47:       row = 15
 48:       call VecSetValue(v1,row,num,INSERT_VALUES,ierr)
 49:       row = 3
 50:       call VecSetValue(v1,row,num,INSERT_VALUES,ierr)
 51:       row = 7
 52:       call VecSetValue(v1,row,num,INSERT_VALUES,ierr)
 53:       row = 11
 54:       call VecSetValue(v1,row,num,INSERT_VALUES,ierr)

 56:       call VecAssemblyBegin(v1,ierr)
 57:       call VecAssemblyEnd(v1,ierr)

 59:       num=0.0
 60:       call VecScale(v2,num,ierr)
 61:       call 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:       call ISCreateGeneral(PETSC_COMM_SELF,eight,from,fromis,ierr)
 82:       call ISCreateGeneral(PETSC_COMM_SELF,eight,to,tois,ierr)
 83:       call VecScatterCreate(v1,fromis,v2,tois,scat1,ierr)
 84:       call VecScatterCreate(v1,fromis,v3,tois,scat2,ierr)
 85:       call ISDestroy(fromis,ierr)
 86:       call ISDestroy(tois,ierr)
 87: 
 88:       call VecScatterBegin(v1,v2,ADD_VALUES,SCATTER_FORWARD,scat1,ierr)
 89:       call VecScatterEnd(v1,v2,ADD_VALUES,SCATTER_FORWARD,scat1,ierr)
 90: 
 91:       call VecScatterBegin(v1,v3,ADD_VALUES,SCATTER_FORWARD,scat2,ierr)
 92:       call VecScatterEnd(v1,v3,ADD_VALUES,SCATTER_FORWARD,scat2,ierr)
 93: 
 94:       if (rank.eq.0) print *, "V1"
 95:       call VecView(v1,PETSC_VIEWER_STDOUT_WORLD,ierr)
 96:       if (rank.eq.0) print *, "V2"
 97:       call VecView(v2,PETSC_VIEWER_STDOUT_WORLD,ierr)
 98:       if (rank.eq.0) then
 99:          print *, "V3"
100:          call VecView(v3,PETSC_VIEWER_STDOUT_SELF,ierr)
101:       end if

103:       call VecScatterDestroy(scat1,ierr)
104:       call VecScatterDestroy(scat2,ierr)
105:       call VecDestroy(v1,ierr)
106:       call VecDestroy(v2,ierr)
107:       call VecDestroy(v3,ierr)

109:       call PetscFinalize(ierr)
110: 
111:       end
112: