Skip to content
Newer
Older
100644 413 lines (360 sloc) 8.97 KB
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
1 module types
2 implicit none
3 private
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
4 public dp, i64
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
5 integer, parameter :: dp=kind(0.d0) ! double precision
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
6 integer, parameter :: i64 = selected_int_kind(18) ! At least 64-bit integer
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
7 end module
8
9
10 module utils
11 ! Various utilities
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
12 use types, only: dp, i64
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
13 implicit none
14 private
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
15 public trace, mean, std, init_random_seed, randn, assert, stop_error, &
16 sysclock2ms
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
17
18 contains
19
20 subroutine stop_error(msg)
21 ! Aborts the program with nonzero exit code
22 !
23 ! The statement "stop msg" will return 0 exit code when compiled using
24 ! gfortran. stop_error() uses the statement "stop 1" which returns an exit code
25 ! 1 and a print statement to print the message.
26 !
27 ! Example
28 ! -------
29 !
30 ! call stop_error("Invalid argument")
31
32 character(len=*) :: msg ! Message to print on stdout
33 print *, msg
34 stop 1
35 end subroutine
36
37 subroutine assert(condition)
38 ! If condition == .false., it aborts the program.
39 !
40 ! Arguments
41 ! ---------
42 !
43 logical, intent(in) :: condition
44 !
45 ! Example
46 ! -------
47 !
48 ! call assert(a == 5)
49
50 if (.not. condition) call stop_error("Assert failed.")
51 end subroutine
52
53 real(dp) function trace(A) result(t)
54 real(dp), intent(in) :: A(:, :)
55 integer :: i
56 t = 0
57 do i = 1, size(A, 1)
58 t = t + A(i, i)
59 end do
60 end function
61
62 real(dp) function mean(x) result(t)
63 real(dp), intent(in) :: x(:)
64 t = sum(x) / size(x)
65 end function
66
67 real(dp) function std(x) result(t)
68 real(dp), intent(in) :: x(:)
69 t = sqrt(mean(abs(x - mean(x))**2))
70 end function
71
72 subroutine init_random_seed()
73 integer :: i, n, clock
74 integer, allocatable :: seed(:)
75 call random_seed(size=n)
76 allocate(seed(n))
77 call system_clock(count=clock)
78 seed = clock + 37 * [ (i - 1, i = 1, n) ]
79 call random_seed(put=seed)
80 end subroutine
81
82 FUNCTION rnorm() RESULT( fn_val )
83
84 ! This subroutine was taken from: http://jblevins.org/mirror/amiller/rnorm.f90
85
86 ! Generate a random normal deviate using the polar method.
87 ! Reference: Marsaglia,G. & Bray,T.A. 'A convenient method for generating
88 ! normal variables', Siam Rev., vol.6, 260-264, 1964.
89
90 IMPLICIT NONE
91 REAL(dp) :: fn_val
92
93 ! Local variables
94
95 REAL(dp) :: u, sum
96 REAL(dp), SAVE :: v, sln
97 LOGICAL, SAVE :: second = .FALSE.
98 REAL(dp), PARAMETER :: one = 1, vsmall = TINY( one )
99
100 IF (second) THEN
101 ! If second, use the second random number generated on last call
102
103 second = .false.
104 fn_val = v*sln
105
106 ELSE
107 ! First call; generate a pair of random normals
108
109 second = .true.
110 DO
111 CALL RANDOM_NUMBER( u )
112 CALL RANDOM_NUMBER( v )
113 u = SCALE( u, 1 ) - one
114 v = SCALE( v, 1 ) - one
115 sum = u*u + v*v + vsmall ! vsmall added to prevent LOG(zero) / zero
116 IF(sum < one) EXIT
117 END DO
118 sln = SQRT(- SCALE( LOG(sum), 1 ) / sum)
119 fn_val = u*sln
120 END IF
121
122 RETURN
123 END FUNCTION rnorm
124
125 subroutine randn(A)
126 real(dp), intent(out) :: A(:, :)
127 integer :: i, j
128 do j = 1, size(A, 2)
129 do i = 1, size(A, 1)
130 A(i, j) = rnorm()
131 end do
132 end do
133 end subroutine
134
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
135 ! Convert a number of clock ticks, as returned by system_clock called
136 ! with integer(i64) arguments, to milliseconds
137 function sysclock2ms(t)
138 integer(i64), intent(in) :: t
139 integer(i64) :: rate
140 real(dp) :: sysclock2ms, r
141 call system_clock(count_rate=rate)
142 r = 1000._dp / rate
143 sysclock2ms = t * r
144 end function sysclock2ms
145
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
146 end module
147
148
149
150 module bench
88946be @certik Add the rest of the benchmarks
certik authored Jun 6, 2012
151 use utils, only: trace, randn, std, mean, stop_error
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
152 use types, only: dp
153 implicit none
154 private
88946be @certik Add the rest of the benchmarks
certik authored Jun 7, 2012
155 public fib, parse_int, quicksort, mandelperf, pisum, randmatstat, randmatmul
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
156
157 contains
158
159 integer recursive function fib(n) result(r)
160 integer, intent(in) :: n
161 if (n < 2) then
162 r = n
163 else
164 r = fib(n-1) + fib(n-2)
165 end if
166 end function
167
88946be @certik Add the rest of the benchmarks
certik authored Jun 7, 2012
168 integer function parse_int(s, base) result(n)
169 character(len=*), intent(in) :: s
170 integer, intent(in) :: base
171 integer :: i, d
172 character :: c
173 n = 0
174 do i = 1, len(s)
175 c = s(i:i)
176 d = 0
177 if (ichar(c) >= ichar('0') .and. ichar(c) <= ichar('9')) then
178 d = ichar(c) - ichar('0')
179 else if (ichar(c) >= ichar('A') .and. ichar(c) <= ichar('Z')) then
180 d = ichar(c) - ichar('A') + 10
181 else if (ichar(c) >= ichar('a') .and. ichar(c) <= ichar('z')) then
182 d = ichar(c) - ichar('a') + 10
183 else
184 call stop_error("parse_int 1")
185 end if
186
187 if (base <= d) call stop_error("parse_int 2")
188 n = n*base + d
189 end do
190
191 end function
192
193
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
194 integer function mandel(z0) result(r)
195 complex(dp), intent(in) :: z0
196 complex(dp) :: c, z
f7832ed @magistere Fix pisum and mandelbrot microbenchmarks
magistere authored Oct 12, 2013
197 integer :: n, maxiter
198 maxiter = 80
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
199 z = z0
200 c = z0
f7832ed @magistere Fix pisum and mandelbrot microbenchmarks
magistere authored Oct 12, 2013
201 do n = 1, maxiter
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
202 if (abs(z) > 2) then
f7832ed @magistere Fix pisum and mandelbrot microbenchmarks
magistere authored Oct 12, 2013
203 r = n-1
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
204 return
205 end if
206 z = z**2 + c
207 end do
f7832ed @magistere Fix pisum and mandelbrot microbenchmarks
magistere authored Oct 12, 2013
208 r = maxiter
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
209 end function
210
211 integer function mandelperf() result(mandel_sum)
4f1fb71 @magistere Fix mandelbrot microbenchmark in C, Fortran and Go
magistere authored Oct 16, 2013
212 integer :: re, im
e4c53c0 @jiahao Add volatile keyword to mandelperf in perf.f90
jiahao authored Nov 14, 2014
213 volatile :: mandel_sum
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
214 mandel_sum = 0
4f1fb71 @magistere Fix mandelbrot microbenchmark in C, Fortran and Go
magistere authored Oct 16, 2013
215 re = -20
216 do while (re <= 5)
217 im = -10
218 do while (im <= 10)
219 mandel_sum = mandel_sum + mandel(cmplx(re/10._dp, im/10._dp, dp))
220 im = im + 1
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
221 end do
4f1fb71 @magistere Fix mandelbrot microbenchmark in C, Fortran and Go
magistere authored Oct 16, 2013
222 re = re + 1
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
223 end do
224 end function
225
88946be @certik Add the rest of the benchmarks
certik authored Jun 7, 2012
226 recursive subroutine quicksort(a, lo0, hi)
227 real(dp), intent(inout) :: a(:)
228 integer, intent(in) :: lo0, hi
229 integer :: i, j, lo
230 real(dp) :: pivot, t
231 lo = lo0
232 i = lo
233 j = hi
234 do while (i < hi)
235 pivot = a((lo+hi)/2)
236 do while (i <= j)
237 do while (a(i) < pivot)
238 i = i + 1
239 end do
240 do while (a(j) > pivot)
241 j = j - 1
242 end do
243 if (i <= j) then
244 t = a(i)
245 a(i) = a(j)
246 a(j) = t
247 i = i + 1
248 j = j - 1
249 end if
250 end do
251 if (lo < j) call quicksort(a, lo, j)
252 lo = i
253 j = hi
254 end do
255 end subroutine
256
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
257 real(dp) function pisum() result(s)
258 integer :: j, k
f7832ed @magistere Fix pisum and mandelbrot microbenchmarks
magistere authored Oct 12, 2013
259 do j = 1, 500
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
260 s = 0
f7832ed @magistere Fix pisum and mandelbrot microbenchmarks
magistere authored Oct 12, 2013
261 do k = 1, 10000
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
262 s = s + 1._dp / k**2
263 end do
264 end do
265 end function
266
267 subroutine randmatstat(t, s1, s2)
268 integer, intent(in) :: t
269 real(dp), intent(out) :: s1, s2
270 real(dp), allocatable, dimension(:, :) :: a, b, c, d, P, Q, X
271 real(dp), allocatable :: v(:), w(:)
272 integer :: n, i
273 n = 5
274 allocate(a(n, n), b(n, n), c(n, n), d(n, n))
275 allocate(P(4*n, n), Q(2*n, 2*n), X(2*n, 2*n))
276 allocate(v(t), w(t))
277 do i = 1, t
278 call randn(a)
279 call randn(b)
280 call randn(c)
281 call randn(d)
282 P(:n, :)=a; P(n+1:2*n, :)=b; P(2*n+1:3*n, :)=c; P(3*n+1:, :)=d
283 Q(:n, :n) = a; Q(n+1:, :n) = b
284 Q(:n, n+1: ) = c; Q(n+1:, n+1: ) = d
285 X = matmul(transpose(P), P)
286 X = matmul(X, X)
287 X = matmul(X, X)
288 v(i) = trace(X)
289 X = matmul(transpose(Q), Q)
290 X = matmul(X, X)
291 X = matmul(X, X)
292 w(i) = trace(X)
293 end do
294 s1 = std(v) / mean(v)
295 s2 = std(w) / mean(w)
296 end subroutine
297
298 subroutine randmatmul(n, C)
299 integer, intent(in) :: n
300 real(dp), intent(out), allocatable :: C(:, :)
301 real(dp), allocatable :: A(:, :), B(:, :)
302 allocate(A(n, n), B(n, n), C(n, n))
303 call random_number(A)
304 call random_number(B)
305 C = matmul(A, B)
306 end subroutine
307
308 end module
309
310 program perf
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
311 use types, only: dp, i64
312 use utils, only: assert, init_random_seed, sysclock2ms
88946be @certik Add the rest of the benchmarks
certik authored Jun 7, 2012
313 use bench, only: fib, parse_int, quicksort, mandelperf, pisum, randmatstat, &
314 randmatmul
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
315 implicit none
316
3b0b665 @certik Increase the number of runs to get nonzero timings
certik authored Jun 7, 2012
317 integer, parameter :: NRUNS = 1000
318 integer :: i, f, n, m, k, k2
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
319 integer(i64) :: t1, t2, tmin
320 real(dp) :: pi, s1, s2
88946be @certik Add the rest of the benchmarks
certik authored Jun 7, 2012
321 real(dp), allocatable :: C(:, :), d(:)
322 character(len=11) :: s
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
323
324 call init_random_seed()
325
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
326 tmin = huge(0_i64)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
327 do i = 1, 5
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
328 call system_clock(t1)
3b0b665 @certik Increase the number of runs to get nonzero timings
certik authored Jun 7, 2012
329 do k = 1, NRUNS
330 f = fib(20)
331 end do
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
332 call system_clock(t2)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
333 if (t2-t1 < tmin) tmin = t2-t1
334 end do
335 call assert(f == 6765)
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
336 print "('fortran,fib,',f0.6)", sysclock2ms(tmin) / NRUNS
88946be @certik Add the rest of the benchmarks
certik authored Jun 7, 2012
337
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
338 tmin = huge(0_i64)
88946be @certik Add the rest of the benchmarks
certik authored Jun 7, 2012
339 do i = 1, 5
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
340 call system_clock(t1)
3b0b665 @certik Increase the number of runs to get nonzero timings
certik authored Jun 7, 2012
341 do k2 = 1, NRUNS
342 do k = 1, 1000
343 call random_number(s1)
344 n = int(s1*huge(n))
345 write(s, '(z0)') n
346 m = parse_int(s(:len_trim(s)), 16)
347 call assert(m == n)
348 end do
88946be @certik Add the rest of the benchmarks
certik authored Jun 7, 2012
349 end do
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
350 call system_clock(t2)
88946be @certik Add the rest of the benchmarks
certik authored Jun 7, 2012
351 if (t2-t1 < tmin) tmin = t2-t1
352 end do
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
353 print "('fortran,parse_int,',f0.6)", sysclock2ms(tmin) / NRUNS
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
354
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
355 tmin = huge(0_i64)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
356 do i = 1, 5
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
357 call system_clock(t1)
3b0b665 @certik Increase the number of runs to get nonzero timings
certik authored Jun 7, 2012
358 do k = 1, NRUNS
359 f = mandelperf()
360 end do
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
361 call system_clock(t2)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
362 if (t2-t1 < tmin) tmin = t2-t1
363 end do
4f1fb71 @magistere Fix mandelbrot microbenchmark in C, Fortran and Go
magistere authored Oct 16, 2013
364 call assert(f == 14791)
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
365 print "('fortran,mandel,',f0.6)", sysclock2ms(tmin) / NRUNS
88946be @certik Add the rest of the benchmarks
certik authored Jun 7, 2012
366
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
367 tmin = huge(0_i64)
88946be @certik Add the rest of the benchmarks
certik authored Jun 7, 2012
368 do i = 1, 5
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
369 call system_clock(t1)
3b0b665 @certik Increase the number of runs to get nonzero timings
certik authored Jun 7, 2012
370 do k = 1, NRUNS
371 allocate(d(5000))
372 call random_number(d)
373 call quicksort(d, 1, size(d))
374 deallocate(d)
375 end do
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
376 call system_clock(t2)
88946be @certik Add the rest of the benchmarks
certik authored Jun 7, 2012
377 if (t2-t1 < tmin) tmin = t2-t1
378 end do
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
379 print "('fortran,quicksort,',f0.6)", sysclock2ms(tmin) / NRUNS
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
380
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
381 tmin = huge(0_i64)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
382 do i = 1, 5
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
383 call system_clock(t1)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
384 pi = pisum()
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
385 call system_clock(t2)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
386 if (t2-t1 < tmin) tmin = t2-t1
387 end do
388 call assert(abs(pi - 1.644834071848065_dp) < 1e-6_dp)
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
389 print "('fortran,pi_sum,',f0.6)", sysclock2ms(tmin)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
390
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
391 tmin = huge(0_i64)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
392 do i = 1, 5
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
393 call system_clock(t1)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
394 call randmatstat(1000, s1, s2)
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
395 call system_clock(t2)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
396 if (t2-t1 < tmin) tmin = t2-t1
397 end do
d07ed3b @StefanKarpinski perf: spruce up Fortran benchmark stuff.
StefanKarpinski authored Jul 23, 2012
398 ! call assert(s1 > 0.5_dp .and. s1 < 1)
399 ! call assert(s2 > 0.5_dp .and. s2 < 1)
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
400 print "('fortran,rand_mat_stat,',f0.6)", sysclock2ms(tmin)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
401
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
402 tmin = huge(0_i64)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
403 do i = 1, 5
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
404 call system_clock(t1)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
405 call randmatmul(1000, C)
406 call assert(C(1, 1) >= 0)
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
407 call system_clock(t2)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
408 if (t2-t1 < tmin) tmin = t2-t1
409 end do
d3e963f @jabl Use system_clock instead of cpu_time
jabl authored Feb 12, 2015
410 print "('fortran,rand_mat_mul,',f0.6)", sysclock2ms(tmin)
ce91f65 @certik Add Fortran version of the benchmarks
certik authored Jun 6, 2012
411
412 end program
Something went wrong with that request. Please try again.