Actual source code: ex2f.F
1: !
2: ! Description: Creates an index set based on a stride. Views that index set
3: ! and then destroys it.
4: !
5: !/*T
6: ! Concepts: index sets^manipulating a stride index set;
7: ! Concepts: index sets^accessing indices from Fortran
8: !T*/
9: !
10: ! Include petscis.h so we can use PETSc IS objects.
11: !
12: program main
13: implicit none
14: #include finclude/petsc.h
15: #include finclude/petscis.h
17: PetscErrorCode ierr
18: PetscInt i,n,index(1),first,step,val
19: IS set
20: PetscOffset iss
22: #define indices(ib) index(iss + (ib))
24: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
25: n = 10
26: first = 3
27: step = 2
29: ! Create stride index set, starting at 3 with a stride of 2
30: ! Note each processor is generating its own index set
31: ! (in this case they 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 code can
37: ! directly access the array storing a PETSc index set with ISGetIndices().
38: ! The user declares an array (index(1)) and index variable (iss), which are
39: ! then used together to allow the Fortran to directly manipulate the PETSc array
41: call ISGetIndices(set,index,iss,ierr)
42: write(6,20)
43: do 10 i=1,n
44: ! Bug in IRIX64 f90 compiler - write cannot handle integer(integer*8) correctly
45: val = indices(i)
46: write(6,30) val
47: 10 continue
48: 20 format('Printing indices directly')
49: 30 format(i3)
50: call ISRestoreIndices(set,index,iss,ierr)
52: ! Determine information on stride
54: call ISStrideGetInfo(set,first,step,ierr)
55: if (first .ne. 3 .or. step .ne. 2) then
56: print*,'Stride info not correct!'
57: endif
59: call ISDestroy(set,ierr)
60: call PetscFinalize(ierr)
61: end