Actual source code: ex3f90.F

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

 22:  #include finclude/petsc.h
 23:  #include finclude/petscis.h
 24: #include "finclude/petscis.h90"

 26:       integer n,ierr,bs,issize,inputindices(4)
 27:       integer, pointer :: indices(:)
 28:       IS       set
 29:       PetscTruth isablock;

 31:       n               = 4
 32:       bs              = 3
 33:       inputindices(1) = 0
 34:       inputindices(2) = 3
 35:       inputindices(3) = 9
 36:       inputindices(4) = 12
 37: 
 38:       call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
 39: 
 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,set,ierr)
 47:       call ISView(set,PETSC_VIEWER_STDOUT_SELF,ierr)

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

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

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

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

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

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