Actual source code: ex2f.F
1: !
2: ! Description: Creates an index set based on a stride. Views that
3: ! index set and then destroys it.
4: !
5: !
6: ! Include petscis.h so we can use PETSc IS objects.
7: !
8: program main
9: #include <petsc/finclude/petscis.h>
10: use petscis
11: implicit none
13: PetscErrorCode ierr
14: PetscInt i,n,index(1),first,step,val
15: IS set
16: PetscOffset iss
18: #define indices(ib) index(iss + (ib))
20: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
21: if (ierr .ne. 0) then
22: print*,'Unable to initialize PETSc'
23: stop
24: endif
25: n = 10
26: first = 3
27: step = 2
29: ! Create stride index set, starting at 3 with a stride of 2 Note
30: ! each processor is generating its own index set (in this case they
31: ! are all identical)
33: call ISCreateStride(PETSC_COMM_SELF,n,first,step,set,ierr)
34: call ISView(set,PETSC_VIEWER_STDOUT_SELF,ierr)
36: ! Extract the indice values from the set. Demonstrates how a Fortran
37: ! code can directly access the array storing a PETSc index set with
38: ! ISGetIndices(). The user declares an array (index(1)) and index
39: ! variable (iss), which are then used together to allow the Fortran
40: ! to directly manipulate the PETSc array
42: call ISGetIndices(set,index,iss,ierr)
43: write(6,20)
44: ! Bug in IRIX64 f90 compiler - write cannot handle
45: ! integer(integer*8) correctly
46: do 10 i=1,n
47: val = indices(i)
48: write(6,30) val
49: 10 continue
50: 20 format('Printing indices directly')
51: 30 format(i3)
52: call ISRestoreIndices(set,index,iss,ierr)
54: ! Determine information on stride
56: call ISStrideGetInfo(set,first,step,ierr)
57: if (first .ne. 3 .or. step .ne. 2) then
58: print*,'Stride info not correct!'
59: endif
61: call ISDestroy(set,ierr)
62: call PetscFinalize(ierr)
63: end
65: !/*TEST
66: !
67: ! test:
68: !
69: !TEST*/