Actual source code: ex1f.F90

petsc-3.8.4 2018-03-24
Report Typos and Errors
  1: !
  2: !
  3: !  Formatted test for IS general routines
  4: !
  5:       program main
  6:  #include <petsc/finclude/petscis.h>
  7:       use petscis
  8:       implicit none

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

 17:        call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
 18:        if (ierr .ne. 0) then
 19:          print*,'Unable to initialize PETSc'
 20:          stop
 21:        endif
 22:        call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
 23:        call MPI_Comm_size(PETSC_COMM_WORLD,size,ierr)

 25: !     Test IS of size 0

 27:        n = 0
 28:        call ISCreateGeneral(PETSC_COMM_SELF,n,indices,PETSC_COPY_VALUES,is,ierr);CHKERRA(ierr);
 29:        call ISGetLocalSize(is,n,ierr);CHKERRA(ierr);
 30:        if (n .ne. 0) then
 31:          SETERRA(PETSC_COMM_SELF,1,'Error getting size of zero IS')
 32:        endif
 33:        call ISDestroy(is,ierr)


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

 51: !     Check identity and permutation

 53:       call ISPermutation(is,flag,ierr);CHKERRA(ierr)
 54:       if (flag) then
 55:          SETERRA(PETSC_COMM_SELF,1,'Error checking permutation')
 56:       endif
 57:       call ISIdentity(is,flag,ierr);CHKERRA(ierr)
 58:       if (.not. flag) then
 59:          SETERRA(PETSC_COMM_SELF,1,'Error checking identity')
 60:       endif
 61:       call ISSetPermutation(is,ierr);CHKERRA(ierr)
 62:       call ISSetIdentity(is,ierr);CHKERRA(ierr)
 63:       call ISPermutation(is,flag,ierr);CHKERRA(ierr)
 64:       if (.not. flag) then
 65:          SETERRA(PETSC_COMM_SELF,1,'Error checking permutation second time')
 66:       endif
 67:       call ISIdentity(is,flag,ierr);CHKERRA(ierr)
 68:       if (.not. flag) then
 69:          SETERRA(PETSC_COMM_SELF,1,'Error checking identity second time')
 70:       endif

 72: !     Check equality of index sets

 74:       call ISEqual(is,is,flag,ierr);CHKERRA(ierr)
 75:       if (.not. flag) then
 76:          SETERRA(PETSC_COMM_SELF,1,'Error checking equal')
 77:       endif

 79: !     Sorting

 81:       call ISSort(is,ierr);CHKERRA(ierr)
 82:       call ISSorted(is,flag,ierr);CHKERRA(ierr)
 83:       if (.not. flag) then
 84:          SETERRA(PETSC_COMM_SELF,1,'Error checking sorted')
 85:       endif

 87: !     Thinks it is a different type?

 89:       call PetscObjectTypeCompare(is,ISSTRIDE,flag,ierr);CHKERRA(ierr)
 90:       if (flag) then
 91:          SETERRA(PETSC_COMM_SELF,1,'Error checking stride')
 92:       endif
 93:       call PetscObjectTypeCompare(is,ISBLOCK,flag,ierr);CHKERRA(ierr)
 94:       if (flag) then
 95:          SETERRA(PETSC_COMM_SELF,1,'Error checking block')
 96:       endif

 98:       call ISDestroy(is,ierr);CHKERRA(ierr)

100: !     Inverting permutation

102:       do 30, i=1,n
103:         indices(i) = n - i
104:  30   continue

106:       call ISCreateGeneral(PETSC_COMM_SELF,n,indices,PETSC_COPY_VALUES,is,ierr);CHKERRA(ierr)
107:       call ISSetPermutation(is,ierr);CHKERRA(ierr)
108:       call ISInvertPermutation(is,PETSC_DECIDE,newis,ierr);CHKERRA(ierr)
109:       call ISGetIndices(newis,ii,iis,ierr);CHKERRA(ierr)
110:       do 40, i=1,n
111:         if (ii(iis+i) .ne. n - i) then
112:           SETERRA(PETSC_COMM_SELF,1,'Error getting permutation indices')
113:        endif
114:  40   continue
115:       call ISRestoreIndices(newis,ii,iis,ierr);CHKERRA(ierr)
116:       call ISDestroy(newis,ierr);CHKERRA(ierr)
117:       call ISDestroy(is,ierr);CHKERRA(ierr)
118:       call PetscFinalize(ierr)
119:       end