Actual source code: ex30f.F

  1: !
  2: !    "$Id: ex30f.F,v 1.6 2001/08/07 03:02:26 balay Exp $";
  3: !
  4: !  Tests parallel to parallel scatter where a to from index are
  5: !  duplicated
  6:       program main
  7:       implicit none

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

 13:       integer ierr, nlocal, n, irank, iprocs
 14:       integer 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:       call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
 23:       call MPI_COMM_RANK(PETSC_COMM_WORLD,irank,ierr)
 24:       call MPI_COMM_SIZE(PETSC_COMM_WORLD,iprocs,ierr)
 25:       if (iprocs.ne.4) then
 26:          print *, 'Four processor test'
 27:          stop
 28:       end if
 29: 
 30:       call VecCreateMPI(PETSC_COMM_WORLD,nlocal*2,n*2,v1,ierr)
 31:       call VecCreateMPI(PETSC_COMM_WORLD,nlocal,n,v2,ierr)
 32:       call VecCreateSeq(PETSC_COMM_SELF,n,v3,ierr)

 34:       num=2.0
 35:       call VecSetValue(v1,1,num,INSERT_VALUES,ierr)
 36:       call VecSetValue(v1,5,num,INSERT_VALUES,ierr)
 37:       call VecSetValue(v1,9,num,INSERT_VALUES,ierr)
 38:       call VecSetValue(v1,13,num,INSERT_VALUES,ierr)
 39:       num=1.0
 40:       call VecSetValue(v1,15,num,INSERT_VALUES,ierr)
 41:       call VecSetValue(v1,3,num,INSERT_VALUES,ierr)
 42:       call VecSetValue(v1,7,num,INSERT_VALUES,ierr)
 43:       call VecSetValue(v1,11,num,INSERT_VALUES,ierr)

 45:       call VecAssemblyBegin(v1,ierr)
 46:       call VecAssemblyEnd(v1,ierr)

 48:       num=0.0
 49:       call VecScale(num,v2,ierr)
 50:       call VecScale(num,v3,ierr)

 52:       from(1)=1
 53:       from(2)=5
 54:       from(3)=9
 55:       from(4)=13
 56:       from(5)=3
 57:       from(6)=7
 58:       from(7)=11
 59:       from(8)=15
 60:       to(1)=0
 61:       to(2)=0
 62:       to(3)=0
 63:       to(4)=0
 64:       to(5)=7
 65:       to(6)=7
 66:       to(7)=7
 67:       to(8)=7

 69:       call ISCreateGeneral(PETSC_COMM_SELF,8,from,fromis,ierr)
 70:       call ISCreateGeneral(PETSC_COMM_SELF,8,to,tois,ierr)
 71:       call VecScatterCreate(v1,fromis,v2,tois,scat1,ierr)
 72:       call VecScatterCreate(v1,fromis,v3,tois,scat2,ierr)
 73:       call ISDestroy(fromis,ierr)
 74:       call ISDestroy(tois,ierr)
 75: 
 76:       call VecScatterBegin(v1,v2,ADD_VALUES,SCATTER_FORWARD,scat1,ierr)
 77:       call VecScatterEnd(v1,v2,ADD_VALUES,SCATTER_FORWARD,scat1,ierr)
 78: 
 79:       call VecScatterBegin(v1,v3,ADD_VALUES,SCATTER_FORWARD,scat2,ierr)
 80:       call VecScatterEnd(v1,v3,ADD_VALUES,SCATTER_FORWARD,scat2,ierr)
 81: 
 82:       if (irank.eq.0) print *, "V1"
 83:       call VecView(v1,PETSC_VIEWER_STDOUT_WORLD,ierr)
 84:       if (irank.eq.0) print *, "V2"
 85:       call VecView(v2,PETSC_VIEWER_STDOUT_WORLD,ierr)
 86:       if (irank.eq.0) then
 87:          print *, "V3"
 88:          call VecView(v3,PETSC_VIEWER_STDOUT_SELF,ierr)
 89:       end if

 91:       call VecScatterDestroy(scat1,ierr)
 92:       call VecScatterDestroy(scat2,ierr)
 93:       call VecDestroy(v1,ierr)
 94:       call VecDestroy(v2,ierr)
 95:       call VecDestroy(v3,ierr)

 97:       call PetscFinalize(ierr)
 98: 
 99:       end
100: