Actual source code: ex2f.F90
petsc-3.8.4 2018-03-24
1: !
2: ! Formatted Test for IS stride routines
3: !
4: program main
5: #include <petsc/finclude/petscis.h>
6: use petscis
7: implicit none
9: PetscErrorCode ierr
10: PetscInt i,n,ii(1),start
11: PetscInt stride,ssize,first
12: IS is
13: PetscBool flag
14: PetscOffset iis
16: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
17: if (ierr .ne. 0) then
18: print*,'Unable to initialize PETSc'
19: stop
20: endif
21:
22: ! Test IS of size 0
23: ssize = 0
24: stride = 0
25: first = 2
26: call ISCreateStride(PETSC_COMM_SELF,ssize,stride,first,is,ierr)
27: call ISGetLocalSize(is,n,ierr)
28: if (n .ne. 0) then
29: SETERRA(PETSC_COMM_SELF,1,'Wrong result from ISCreateStride ')
30: endif
31: call ISStrideGetInfo(is,start,stride,ierr)
32: if (start .ne. 0) then
33: SETERRA(PETSC_COMM_SELF,1,'Wrong result from ISStrideGetInfo ')
34: endif
35: if (stride .ne. 2) then
36: SETERRA(PETSC_COMM_SELF,1,'Wrong result from ISStrideGetInfo ')
37: endif
38: call PetscObjectTypeCompare(is,ISSTRIDE,flag,ierr)
39: if (.not. flag) then
40: SETERRA(PETSC_COMM_SELF,1,'Wrong result from PetscObjectTypeCompare')
41: endif
42: call ISGetIndices(is,ii,iis,ierr)
43: call ISRestoreIndices(is,ii,iis,ierr)
44: call ISDestroy(is,ierr)
46: ! Test ISGetIndices()
48: ssize = 10000
49: stride = -8
50: first = 3
51: call ISCreateStride(PETSC_COMM_SELF,ssize,stride,first,is,ierr)
52: call ISGetLocalSize(is,n,ierr)
53: call ISGetIndices(is,ii,iis,ierr)
54: do 10, i=1,10000
55: if (ii(i+iis) .ne. -11 + 3*i) then
56: SETERRA(PETSC_COMM_SELF,1,'Wrong result from ISGetIndices')
57: endif
58: 10 continue
59: call ISRestoreIndices(is,ii,iis,ierr)
60: call ISDestroy(is,ierr)
62: call PetscFinalize(ierr)
63: end