Actual source code: ex1f.F90

  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(1004),ii(1)
 12:        PetscMPIInt size,rank
 13:        PetscOffset iis
 14:        IS          is,newis
 15:        PetscBool   flag,compute,permanent

 17:        PetscCallA(PetscInitialize(ierr))
 18:        PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr))
 19:        PetscCallMPIA(MPI_Comm_size(PETSC_COMM_WORLD,size,ierr))

 21: !     Test IS of size 0

 23:        n = 0
 24:        PetscCallA(ISCreateGeneral(PETSC_COMM_SELF,n,indices,PETSC_COPY_VALUES,is,ierr))
 25:        PetscCallA(ISGetLocalSize(is,n,ierr))
 26:        if (n .ne. 0) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error getting size of zero IS'); endif
 27:        PetscCallA(ISDestroy(is,ierr))

 29: !     Create large IS and test ISGetIndices(,ierr))
 30: !     fortran indices start from 1 - but IS indices start from 0
 31:       n = 1000 + rank
 32:       do 10, i=1,n
 33:         indices(i) = rank + i-1
 34:  10   continue
 35:       PetscCallA(ISCreateGeneral(PETSC_COMM_SELF,n,indices,PETSC_COPY_VALUES,is,ierr))
 36:       PetscCallA(ISGetIndices(is,ii,iis,ierr))
 37:       do 20, i=1,n
 38:         if (ii(i+iis) .ne. indices(i)) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error getting indices'); endif
 39:  20   continue
 40:       PetscCallA(ISRestoreIndices(is,ii,iis,ierr))

 42: !     Check identity and permutation

 44:       compute = PETSC_TRUE
 45:       permanent = PETSC_FALSE
 46:       PetscCallA(ISPermutation(is,flag,ierr))
 47:       if (flag) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error checking permutation'); endif
 48:       PetscCallA(ISGetInfo(is,IS_PERMUTATION,IS_LOCAL,compute,flag,ierr))
 49:       !if ((rank .eq. 0) .and. (.not. flag)) SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,"ISGetInfo(IS_PERMUTATION,IS_LOCAL)")
 50:       !if (rank .eq. 0 .and. flag) SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,"ISGetInfo(IS_PERMUTATION,IS_LOCAL)")
 51:       PetscCallA(ISIdentity(is,flag,ierr))
 52:       if ((rank .eq. 0) .and. (.not. flag)) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error checking identity'); endif
 53:       if ((rank .ne. 0) .and. flag) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error checking identity'); endif
 54:       PetscCallA(ISSetInfo(is,IS_PERMUTATION,IS_LOCAL,permanent,PETSC_TRUE,ierr))
 55:       PetscCallA(ISSetInfo(is,IS_IDENTITY,IS_LOCAL,permanent,PETSC_TRUE,ierr))
 56:       PetscCallA(ISGetInfo(is,IS_PERMUTATION,IS_LOCAL,compute,flag,ierr))
 57:       if (.not. flag) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error checking permutation second time'); endif
 58:       PetscCallA(ISGetInfo(is,IS_IDENTITY,IS_LOCAL,compute,flag,ierr))
 59:       if (.not. flag) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error checking identity second time'); endif
 60:       PetscCallA(ISClearInfoCache(is,PETSC_TRUE,ierr))

 62: !     Check equality of index sets

 64:       PetscCallA(ISEqual(is,is,flag,ierr))
 65:       if (.not. flag) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error checking equal'); endif

 67: !     Sorting

 69:       PetscCallA(ISSort(is,ierr))
 70:       PetscCallA(ISSorted(is,flag,ierr))
 71:       if (.not. flag) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error checking sorted'); endif

 73: !     Thinks it is a different type?

 75:       PetscCallA(PetscObjectTypeCompare(is,ISSTRIDE,flag,ierr))
 76:       if (flag) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error checking stride'); endif
 77:       PetscCallA(PetscObjectTypeCompare(is,ISBLOCK,flag,ierr))
 78:       if (flag) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error checking block'); endif

 80:       PetscCallA(ISDestroy(is,ierr))

 82: !     Inverting permutation

 84:       do 30, i=1,n
 85:         indices(i) = n - i
 86:  30   continue

 88:       PetscCallA(ISCreateGeneral(PETSC_COMM_SELF,n,indices,PETSC_COPY_VALUES,is,ierr))
 89:       PetscCallA(ISSetPermutation(is,ierr))
 90:       PetscCallA(ISInvertPermutation(is,PETSC_DECIDE,newis,ierr))
 91:       PetscCallA(ISGetIndices(newis,ii,iis,ierr))
 92:       do 40, i=1,n
 93:         if (ii(iis+i) .ne. n - i) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error getting permutation indices'); endif
 94:  40   continue
 95:       PetscCallA(ISRestoreIndices(newis,ii,iis,ierr))
 96:       PetscCallA(ISDestroy(newis,ierr))
 97:       PetscCallA(ISDestroy(is,ierr))
 98:       PetscCallA(PetscFinalize(ierr))
 99:       end

101: !/*TEST
102: !
103: !   test:
104: !     nsize: {{1 2 3 4 5}}
105: !     output_file: output/ex1_1.out
106: !
107: !TEST*/