forked from AMReX-Codes/amrex
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathsort_i.f90
248 lines (206 loc) · 5.82 KB
/
sort_i.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
module sort_i_module
! Adapted from Meissner, Adams, et.al., and NR.
implicit none
interface sort
module procedure heapsort_i
module procedure heapsort_indirect_i
end interface sort
private less_i
private swap_i
contains
pure function less_i(a,b) result(r)
logical :: r
integer, intent(in) :: a, b
r = a < b
end function less_i
subroutine swap_i(a,b)
integer, intent(inout) :: a, b
integer :: t
t = a; a = b; b = t
end subroutine swap_i
subroutine heapsort_i( array, cmp)
integer, dimension(:), intent(in out) :: array
interface
function cmp(a,b) result(r)
implicit none
logical :: r
integer, intent(in) :: a, b
end function cmp
end interface
optional cmp
if ( present(cmp) ) then
call he_sort_cmp(array,cmp)
else
call he_sort_no_cmp(array)
end if
contains
subroutine he_sort_cmp(array,cmp)
integer, dimension(:), intent(in out) :: array
interface
function cmp(a,b) result(r)
implicit none
logical :: r
integer, intent(in) :: a, b
end function cmp
end interface
integer :: n, i
n = size (array)
do i = n / 2, 1, -1
call peck_cmp(array,i,n,cmp)
end do
do i = n, 2, -1
call swap_i(array(1),array(i))
call peck_cmp(array,1, i-1,cmp)
end do
end subroutine he_sort_cmp
subroutine he_sort_no_cmp(array)
integer, dimension(:), intent(in out) :: array
integer :: n, i
n = size (array)
do i = n / 2, 1, -1
call peck_no_cmp(array,i,n)
end do
do i = n, 2, -1
call swap_i(array(1),array(i))
call peck_no_cmp(array,1, i-1)
end do
end subroutine he_sort_no_cmp
subroutine peck_cmp(array,l,r,cmp)
integer, dimension(:), intent(inout) :: array
integer, intent(in) :: r
integer, intent(in) :: l
interface
function cmp(a,b) result(r)
implicit none
logical :: r
integer, intent(in) :: a, b
end function cmp
end interface
integer :: i
integer :: next
next = array(l)
i = 2*l
do while ( i <= r )
if ( i < r ) then
if ( cmp(array(i), array(i+1)) ) i = i+1
end if
if ( .not. cmp( next, array(i)) ) exit
array(i/2) = array(i)
i = 2*i
end do
array(i/2) = next
end subroutine peck_cmp
subroutine peck_no_cmp(array,l,r)
integer, dimension(:), intent(inout) :: array
integer, intent(in) :: r
integer, intent(in) :: l
integer :: i
integer :: next
next = array(l)
i = 2*l
do while ( i <= r )
if ( i < r ) then
if ( less_i(array(i), array(i+1)) ) i = i+1
end if
if ( .not. less_i( next, array(i)) ) exit
array(i/2) = array(i)
i = 2*i
end do
array(i/2) = next
end subroutine peck_no_cmp
end subroutine heapsort_i
subroutine heapsort_indirect_i( array, iarray, cmp)
integer, dimension(:), intent(in) :: array
integer, dimension(:), intent(out) :: iarray
interface
function cmp(a,b) result(r)
implicit none
logical :: r
integer, intent(in) :: a, b
end function cmp
end interface
optional cmp
if ( present(cmp) ) then
call he_sort_cmp(iarray, cmp)
else
call he_sort_no_cmp(iarray)
end if
contains
subroutine he_sort_no_cmp(iarray)
integer, dimension(:), intent(out) :: iarray
integer :: n, i
n = size (array)
iarray = (/(i,i=1,n)/)
do i = n / 2, 1, -1
call peck_no_cmp(iarray, i, n)
end do
do i = n, 2, -1
call swap_i(iarray(1), iarray(i))
call peck_no_cmp(iarray, 1, i-1)
end do
end subroutine he_sort_no_cmp
subroutine peck_no_cmp(iarray, l, r)
integer, dimension(:), intent(inout) :: iarray
integer, intent(in) :: r
integer, intent(in) :: l
integer :: i
integer :: next
next = iarray(l)
i = 2*l
do while ( i <= r )
if ( i < r ) then
if ( less_i(array(iarray(i)), array(iarray(i+1))) ) i = i+1
end if
if ( .not. less_i( array(next), array(iarray(i))) ) exit
iarray(i/2) = iarray(i)
i = 2*i
end do
iarray(i/2) = next
end subroutine peck_no_cmp
subroutine he_sort_cmp(iarray, cmp)
integer, dimension(:), intent(out) :: iarray
interface
function cmp(a,b) result(r)
implicit none
logical :: r
integer, intent(in) :: a, b
end function cmp
end interface
integer :: n, i
n = size (array)
iarray = (/(i,i=1,n)/)
do i = n / 2, 1, -1
call peck_cmp(iarray, i, n, cmp)
end do
do i = n, 2, -1
call swap_i(iarray(1), iarray(i))
call peck_cmp(iarray, 1, i-1, cmp)
end do
end subroutine he_sort_cmp
subroutine peck_cmp(iarray, l, r, cmp)
integer, dimension(:), intent(inout) :: iarray
integer, intent(in) :: r
integer, intent(in) :: l
interface
function cmp(a,b) result(r)
implicit none
logical :: r
integer, intent(in) :: a, b
end function cmp
end interface
integer :: i
integer :: next
next = iarray(l)
i = 2*l
do while ( i <= r )
if ( i < r ) then
if ( cmp(array(iarray(i)), array(iarray(i+1))) ) i = i+1
end if
if ( .not. cmp( array(next), array(iarray(i))) ) exit
iarray(i/2) = iarray(i)
i = 2*i
end do
iarray(i/2) = next
end subroutine peck_cmp
end subroutine heapsort_indirect_i
end module sort_i_module