Actual source code: ex1f.F

  1: !
  2: !
  3: !  Formatted test for IS general routines
  4: !
  5:       program main
  6:       implicit none
 7:  #include finclude/petscsys.h
 8:  #include finclude/petscis.h


 11:        PetscErrorCode ierr
 12:        PetscInt i,n,indices(1000),ii(1)
 13:        PetscMPIInt size,rank
 14:        PetscOffset iis
 15:        IS          is,newis
 16:        PetscTruth  flag

 18:        call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
 19:            CHKERRQ(ierr)
 20:        call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
 21:        call MPI_Comm_size(PETSC_COMM_WORLD,size,ierr)

 23: !     Test IS of size 0

 25:        n = 0
 26:        call ISCreateGeneral(PETSC_COMM_SELF,n,indices,is,ierr)
 27:            CHKERRQ(ierr)
 28:        call ISGetLocalSize(is,n,ierr)
 29:            CHKERRQ(ierr)
 30:        if (n .ne. 0) then
 31:          print*, 'Error getting size of zero IS'
 32:          stop
 33:        endif
 34:        call ISDestroy(is,ierr)


 37: !     Create large IS and test ISGetIndices(,ierr)
 38: !     fortran indices start from 1 - but IS indices start from 0
 39:       n = 1000
 40:       do 10, i=1,n
 41:         indices(i) = i-1
 42:  10   continue
 43:       call ISCreateGeneral(PETSC_COMM_SELF,n,indices,is,ierr)
 44:            CHKERRQ(ierr)
 45:       call ISGetIndices(is,ii,iis,ierr)
 46:            CHKERRQ(ierr)
 47:       do 20, i=1,n
 48:         if (ii(i+iis) .ne. indices(i)) then
 49:            print*, 'Error getting indices'
 50:            stop
 51:         endif
 52:  20   continue
 53:       call ISRestoreIndices(is,ii,iis,ierr)
 54:            CHKERRQ(ierr)

 56: !     Check identity and permutation
 57: 
 58:       call ISPermutation(is,flag,ierr)
 59:            CHKERRQ(ierr)
 60:       if (flag) then
 61:          print*, 'Error checking permutation'
 62:          stop
 63:       endif
 64:       call ISIdentity(is,flag,ierr)
 65:            CHKERRQ(ierr)
 66:       if (.not. flag) then
 67:          print*, 'Error checking identity'
 68:          stop
 69:       endif
 70:       call ISSetPermutation(is,ierr)
 71:            CHKERRQ(ierr)
 72:       call ISSetIdentity(is,ierr)
 73:            CHKERRQ(ierr)
 74:       call ISPermutation(is,flag,ierr)
 75:            CHKERRQ(ierr)
 76:       if (.not. flag) then
 77:          print*, 'Error checking permutation second time'
 78:          stop
 79:       endif
 80:       call ISIdentity(is,flag,ierr)
 81:            CHKERRQ(ierr)
 82:       if (.not. flag) then
 83:          print*, 'Error checking identity second time'
 84:          stop
 85:       endif

 87: !     Check equality of index sets

 89:       call ISEqual(is,is,flag,ierr)
 90:            CHKERRQ(ierr)
 91:       if (.not. flag) then
 92:          print*, 'Error checking equal'
 93:          stop
 94:       endif

 96: !     Sorting

 98:       call ISSort(is,ierr)
 99:            CHKERRQ(ierr)
100:       call ISSorted(is,flag,ierr)
101:            CHKERRQ(ierr)
102:       if (.not. flag) then
103:          print*, 'Error checking sorted'
104:          stop
105:       endif

107: !     Thinks it is a different type?

109:       call ISStride(is,flag,ierr)
110:            CHKERRQ(ierr)
111:       if (flag) then
112:          print*, 'Error checking stride'
113:          stop
114:       endif
115:       call ISBlock(is,flag,ierr)
116:            CHKERRQ(ierr)
117:       if (flag) then
118:          print*, 'Error checking block'
119:          stop
120:       endif

122:       call ISDestroy(is,ierr)
123:            CHKERRQ(ierr)

125: !     Inverting permutation

127:       do 30, i=1,n
128:         indices(i) = n - i
129:  30   continue

131:       call ISCreateGeneral(PETSC_COMM_SELF,n,indices,is,ierr)
132:            CHKERRQ(ierr)
133:       call ISSetPermutation(is,ierr)
134:            CHKERRQ(ierr)
135:       call ISInvertPermutation(is,PETSC_DECIDE,newis,ierr)
136:            CHKERRQ(ierr)
137:       call ISGetIndices(newis,ii,iis,ierr)
138:            CHKERRQ(ierr)
139:       do 40, i=1,n
140:         if (ii(iis+i) .ne. n - i) then
141:           print*, 'Error getting permutation indices'
142:           stop
143:        endif
144:  40   continue
145:       call ISRestoreIndices(newis,ii,iis,ierr)
146:            CHKERRQ(ierr)
147:       call ISDestroy(newis,ierr)
148:            CHKERRQ(ierr)
149:       call ISDestroy(is,ierr)
150:            CHKERRQ(ierr)
151:       call PetscFinalize(ierr)
152:       end
153: