1: !
  2: !    Description:  Creates an index set based on blocks of integers. Views that index set
  3: !    and then destroys it.
  4: !
  5: !/*T
  6: !    Concepts: index sets^manipulating a block index set;
  7: !    Concepts: Fortran90^accessing indices in index set;
  8: !
  9: !T*/
 10: !
 11: !  The following include statements are required for Fortran programs
 12: !  that use PETSc index sets:
 13: !     petscsys.h  - base PETSc routines
 14: !     petscis.h     - index sets (IS objects)
 15: !     petscis.h90   - to allow access to Fortran90 features of index sets
 16: !
 17:       program main
 18:       implicit none

 20: #include <petsc/finclude/petscsys.h>
 21: #include <petsc/finclude/petscis.h>
 22: #include <petsc/finclude/petscis.h90>

 24:       PetscErrorCode ierr
 25:       PetscInt n,bs,issize
 26:       PetscInt inputindices(4)
 27:       PetscInt, pointer :: indices(:)
 28:       IS       set
 29:       PetscBool  isablock;

 31:       n               = 4
 32:       bs              = 3
 33:       inputindices(1) = 0
 34:       inputindices(2) = 1
 35:       inputindices(3) = 3
 36:       inputindices(4) = 4

 38:       call PetscInitialize(PETSC_NULL_CHARACTER,ierr)

 40: !
 41: !    Create a block index set. The index set has 4 blocks each of size 3.
 42: !    The indices are {0,1,2,3,4,5,9,10,11,12,13,14}
 43: !    Note each processor is generating its own index set
 44: !    (in this case they are all identical)
 45: !
 46:       call ISCreateBlock(PETSC_COMM_SELF,bs,n,inputindices,                   &
 47:      &                   PETSC_COPY_VALUES,set,ierr)
 48:       call ISView(set,PETSC_VIEWER_STDOUT_SELF,ierr)

 50: !
 51: !    Extract indices from set.
 52: !
 53:       call ISGetLocalSize(set,issize,ierr)
 54:       call ISGetIndicesF90(set,indices,ierr)
 55:       write(6,100)indices
 56:  100  format(12I3)
 57:       call ISRestoreIndicesF90(set,indices,ierr)

 59: !
 60: !    Extract the block indices. This returns one index per block.
 61: !
 62:       call ISBlockGetIndicesF90(set,indices,ierr)
 63:       write(6,200)indices
 64:  200  format(4I3)
 65:       call ISBlockRestoreIndicesF90(set,indices,ierr)

 67: !
 68: !    Check if this is really a block index set
 69: !
 70:       call PetscObjectTypeCompare(set,ISBLOCK,isablock,ierr)
 71:       if (.not. isablock) then
 72:         write(6,*) 'Index set is not blocked!'
 73:       endif

 75: !
 76: !    Determine the block size of the index set
 77: !
 78:       call ISGetBlockSize(set,bs,ierr)
 79:       if (bs .ne. 3) then
 80:         write(6,*) 'Blocksize != 3'
 81:       endif

 83: !
 84: !    Get the number of blocks
 85: !
 86:       call ISBlockGetLocalSize(set,n,ierr)
 87:       if (n .ne. 4) then
 88:         write(6,*) 'Number of blocks != 4'
 89:       endif

 91:       call ISDestroy(set,ierr)
 92:       call PetscFinalize(ierr)
 93:       end