Actual source code: ex3f90.F90

  1: !
  2: !    Description:  Creates an index set based on blocks of integers. Views that index set
  3: !    and then destroys it.
  4: !
  5: !
  6:       program main
  7: #include <petsc/finclude/petscis.h>
  8:       use petscis
  9:       implicit none

 11:       PetscErrorCode ierr
 12:       PetscInt n,bs,issize
 13:       PetscInt inputindices(4)
 14:       PetscInt, pointer :: indices(:)
 15:       IS       set
 16:       PetscBool  isablock;

 18:       n               = 4
 19:       bs              = 3
 20:       inputindices(1) = 0
 21:       inputindices(2) = 1
 22:       inputindices(3) = 3
 23:       inputindices(4) = 4

 25:       call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
 26:       if (ierr .ne. 0) then
 27:          print*,'Unable to initialize PETSc'
 28:          stop
 29:        endif
 30: !
 31: !    Create a block index set. The index set has 4 blocks each of size 3.
 32: !    The indices are {0,1,2,3,4,5,9,10,11,12,13,14}
 33: !    Note each processor is generating its own index set
 34: !    (in this case they are all identical)
 35: !
 36:       call ISCreateBlock(PETSC_COMM_SELF,bs,n,inputindices,PETSC_COPY_VALUES,set,ierr);CHKERRA(ierr)
 37:       call ISView(set,PETSC_VIEWER_STDOUT_SELF,ierr);CHKERRA(ierr)

 39: !
 40: !    Extract indices from set.
 41: !
 42:       call ISGetLocalSize(set,issize,ierr);CHKERRA(ierr)
 43:       call ISGetIndicesF90(set,indices,ierr);CHKERRA(ierr)
 44:       write(6,100)indices
 45:  100  format(12I3)
 46:       call ISRestoreIndicesF90(set,indices,ierr);CHKERRA(ierr)

 48: !
 49: !    Extract the block indices. This returns one index per block.
 50: !
 51:       call ISBlockGetIndicesF90(set,indices,ierr);CHKERRA(ierr)
 52:       write(6,200)indices
 53:  200  format(4I3)
 54:       call ISBlockRestoreIndicesF90(set,indices,ierr);CHKERRA(ierr)

 56: !
 57: !    Check if this is really a block index set
 58: !
 59:       call PetscObjectTypeCompare(set,ISBLOCK,isablock,ierr);CHKERRA(ierr)
 60:       if (.not. isablock) then
 61:         write(6,*) 'Index set is not blocked!'
 62:       endif

 64: !
 65: !    Determine the block size of the index set
 66: !
 67:       call ISGetBlockSize(set,bs,ierr);CHKERRA(ierr)
 68:       if (bs .ne. 3) then
 69:         write(6,*) 'Blocksize != 3'
 70:       endif

 72: !
 73: !    Get the number of blocks
 74: !
 75:       call ISBlockGetLocalSize(set,n,ierr);CHKERRA(ierr)
 76:       if (n .ne. 4) then
 77:         write(6,*) 'Number of blocks != 4'
 78:       endif

 80:       call ISDestroy(set,ierr);CHKERRA(ierr)
 81:       call PetscFinalize(ierr)
 82:       end

 84: !/*TEST
 85: !
 86: !   test:
 87: !      filter: sort -b
 88: !      filter_output: sort -b
 89: !
 90: !TEST*/