Actual source code: ex2f.F90
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, first, step, val
15: IS set
16: PetscInt, pointer :: index(:)
18: PetscCallA(PetscInitialize(ierr))
19: n = 10
20: first = 3
21: step = 2
23: ! Create stride index set, starting at 3 with a stride of 2 Note
24: ! each processor is generating its own index set (in this case they
25: ! are all identical)
27: PetscCallA(ISCreateStride(PETSC_COMM_SELF, n, first, step, set, ierr))
28: PetscCallA(ISView(set, PETSC_VIEWER_STDOUT_SELF, ierr))
30: ! Extract the indices values from the set. Demonstrates how a Fortran
31: ! code can directly access the array storing a PETSc index set with
32: ! ISGetIndices().
34: PetscCallA(ISGetIndices(set, index, ierr))
35: write (6, 20)
36: ! Bug in IRIX64 f90 compiler - write cannot handle
37: ! integer(integer*8) correctly
38: do 10 i = 1, n
39: val = index(i)
40: write (6, 30) val
41: 10 continue
42: 20 format('Printing indices directly')
43: 30 format(i3)
44: PetscCallA(ISRestoreIndices(set, index, ierr))
46: ! Determine information on stride
48: PetscCallA(ISStrideGetInfo(set, first, step, ierr))
49: if (first /= 3 .or. step /= 2) then
50: print *, 'Stride info not correct!'
51: end if
53: PetscCallA(ISDestroy(set, ierr))
54: PetscCallA(PetscFinalize(ierr))
55: end
57: !/*TEST
58: !
59: ! test:
60: !
61: !TEST*/