Actual source code: ex1f.F
1: !
2: ! "$Id: ex1f.F,v 1.19 2001/01/17 22:20:50 bsmith Exp $"
3: !
4: ! Formatted test for IS general routines
5: !
6: program main
7: implicit none
8: #include finclude/petsc.h
9: #include finclude/petscis.h
12: integer i,n,ierr,indices(1000),rank,size,ii(1)
13: PetscOffset iis
14: IS is,newis
15: PetscTruth flag
17: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
18: CHKERRQ(ierr)
19: call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
20: call MPI_Comm_size(PETSC_COMM_WORLD,size,ierr)
22: ! Test IS of size 0
24: call ISCreateGeneral(PETSC_COMM_SELF,0,indices,is,ierr)
25: CHKERRQ(ierr)
26: call ISGetLocalSize(is,n,ierr)
27: CHKERRQ(ierr)
28: if (n .ne. 0) then
29: print*, 'Error getting size of zero IS'
30: stop
31: endif
32: call ISDestroy(is,ierr)
35: ! Create large IS and test ISGetIndices(,ierr)
37: n = 1000
38: do 10, i=1,n
39: indices(i) = rank + i
40: 10 continue
41: call ISCreateGeneral(PETSC_COMM_SELF,n,indices,is,ierr)
42: CHKERRQ(ierr)
43: call ISGetIndices(is,ii,iis,ierr)
44: CHKERRQ(ierr)
45: do 20, i=1,n
46: if (ii(i+iis) .ne. indices(i)) then
47: print*, 'Error getting indices'
48: stop
49: endif
50: 20 continue
51: call ISRestoreIndices(is,ii,iis,ierr)
52: CHKERRQ(ierr)
54: ! Check identity and permutation
55:
56: call ISPermutation(is,flag,ierr)
57: CHKERRQ(ierr)
58: if (flag .eq. PETSC_TRUE) then
59: print*, 'Error checking permutation'
60: stop
61: endif
62: call ISIdentity(is,flag,ierr)
63: CHKERRQ(ierr)
64: if (flag .eq. PETSC_TRUE) then
65: print*, 'Error checking identity'
66: stop
67: endif
68: call ISSetPermutation(is,ierr)
69: CHKERRQ(ierr)
70: call ISSetIdentity(is,ierr)
71: CHKERRQ(ierr)
72: call ISPermutation(is,flag,ierr)
73: CHKERRQ(ierr)
74: if (flag .ne. PETSC_TRUE) then
75: print*, 'Error checking permutation second time'
76: stop
77: endif
78: call ISIdentity(is,flag,ierr)
79: CHKERRQ(ierr)
80: if (flag .ne. PETSC_TRUE) then
81: print*, 'Error checking identity second time'
82: stop
83: endif
85: ! Check equality of index sets
87: call ISEqual(is,is,flag,ierr)
88: CHKERRQ(ierr)
89: if (flag .ne. PETSC_TRUE) then
90: print*, 'Error checking equal'
91: stop
92: endif
94: ! Sorting
96: call ISSort(is,ierr)
97: CHKERRQ(ierr)
98: call ISSorted(is,flag,ierr)
99: CHKERRQ(ierr)
100: if (flag .ne. PETSC_TRUE) then
101: print*, 'Error checking sorted'
102: stop
103: endif
105: ! Thinks it is a different type?
107: call ISStride(is,flag,ierr)
108: CHKERRQ(ierr)
109: if (flag .eq. PETSC_TRUE) then
110: print*, 'Error checking stride'
111: stop
112: endif
113: call ISBlock(is,flag,ierr)
114: CHKERRQ(ierr)
115: if (flag .eq. PETSC_TRUE) then
116: print*, 'Error checking block'
117: stop
118: endif
120: call ISDestroy(is,ierr)
121: CHKERRQ(ierr)
123: ! Inverting permutation
125: do 30, i=1,n
126: indices(i) = n - i
127: 30 continue
129: call ISCreateGeneral(PETSC_COMM_SELF,n,indices,is,ierr)
130: CHKERRQ(ierr)
131: call ISSetPermutation(is,ierr)
132: CHKERRQ(ierr)
133: call ISInvertPermutation(is,PETSC_DECIDE,newis,ierr)
134: CHKERRQ(ierr)
135: call ISGetIndices(newis,ii,iis,ierr)
136: CHKERRQ(ierr)
137: do 40, i=1,n
138: if (ii(iis+i) .ne. n - i) then
139: print*, 'Error getting permutation indices'
140: stop
141: endif
142: 40 continue
143: call ISRestoreIndices(newis,ii,iis,ierr)
144: CHKERRQ(ierr)
145: call ISDestroy(newis,ierr)
146: CHKERRQ(ierr)
147: call ISDestroy(is,ierr)
148: CHKERRQ(ierr)
149: call PetscFinalize(ierr)
150: end
151: