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:        PetscBool   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,PETSC_COPY_VALUES,    &
 27:      &                      is,ierr)
 28:            CHKERRQ(ierr)
 29:        call ISGetLocalSize(is,n,ierr)
 30:            CHKERRQ(ierr)
 31:        if (n .ne. 0) then
 32:          print*, 'Error getting size of zero IS'
 33:          stop
 34:        endif
 35:        call ISDestroy(is,ierr)


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

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

 89: !     Check equality of index sets

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

 98: !     Sorting

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

109: !     Thinks it is a different type?

111:       call PetscTypeCompare(is,ISSTRIDE,flag,ierr)
112:            CHKERRQ(ierr)
113:       if (flag) then
114:          print*, 'Error checking stride'
115:          stop
116:       endif
117:       call PetscTypeCompare(is,ISBLOCK,flag,ierr)
118:            CHKERRQ(ierr)
119:       if (flag) then
120:          print*, 'Error checking block'
121:          stop
122:       endif

124:       call ISDestroy(is,ierr)
125:            CHKERRQ(ierr)

127: !     Inverting permutation

129:       do 30, i=1,n
130:         indices(i) = n - i
131:  30   continue

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