-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathperformance_test.f90
95 lines (79 loc) · 3.42 KB
/
performance_test.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
program test_m_random
use m_random
implicit none
integer, parameter :: dp = kind(0.0d0)
integer, parameter :: i8 = selected_int_kind(18)
integer, parameter :: n_samples = 10*1000*1000
integer :: nn, rng_seed
real(dp) :: mean, variance
real(dp), parameter :: poisson_lambda = 15.0_dp
real(dp) :: time_start, time_end
real(dp), allocatable :: rand_results(:)
type(RNG_t) :: rng
allocate(rand_results(n_samples))
print *, "Testing implementation of m_random.f90"
print *, "This is just checking whether everything works, and by no means"
print *, "a test of the 'randomness' of the pseudo random number generator."
print *, "For these tests, ", n_samples, " values are used"
call system_clock(count=rng_seed)
call cpu_time(time_start)
call random_number(rand_results)
call cpu_time(time_end)
mean = sum(rand_results) / n_samples
variance = sum((rand_results - mean)**2) / n_samples
print *, ""
print *, "For uniform random numbers (built-in), the result is:"
print *, "nanoseconds per number (upper bound)", 1.0e9_dp * (time_end - time_start) / n_samples
print *, "mean/<mean>", mean/0.5_dp
print *, "std dev/<std dev>", sqrt(variance)*sqrt(12.0_dp)
call rng%set_seed([89732_i8, 1892342989_i8])
call cpu_time(time_start)
do nn = 1, n_samples
rand_results(nn) = rng%unif_01()
end do
call cpu_time(time_end)
mean = sum(rand_results) / n_samples
variance = sum((rand_results - mean)**2) / n_samples
print *, ""
print *, "For uniform random numbers (unif_01), the result is:"
print *, "nanoseconds per number (upper bound)", 1.0e9_dp * (time_end - time_start) / n_samples
print *, "mean/<mean>", mean/0.5_dp
print *, "std dev/<std dev>", sqrt(variance)*sqrt(12.0_dp)
call cpu_time(time_start)
do nn = 1, n_samples
rand_results(nn) = rng%normal()
end do
call cpu_time(time_end)
rand_results = rand_results + 1
mean = sum(rand_results) / n_samples
variance = sum((rand_results - mean)**2) / n_samples
print *, ""
print *, "For normal/Gaussian random numbers, the result is:"
print *, "nanoseconds per number (upper bound)", 1.0e9_dp * (time_end - time_start) / n_samples
print *, "mean/<mean>", mean/1.0_dp ! Above we add one to RNG_normal()
print *, "std dev/<std dev>", sqrt(variance)
call cpu_time(time_start)
do nn = 1, n_samples
rand_results(nn) = rng%poisson(poisson_lambda)
end do
call cpu_time(time_end)
mean = sum(rand_results) / n_samples
variance = sum((rand_results - mean)**2) / n_samples
print *, ""
print *, "For Poisson random numbers, the result is:"
print *, "nanoseconds per number (upper bound)", 1.0e9_dp * (time_end - time_start) / n_samples
print *, "mean/<mean>", mean/poisson_lambda ! Above we add one to RNG_normal()
print *, "std dev/<std dev>", sqrt(variance/poisson_lambda)
call cpu_time(time_start)
do nn = 1, n_samples
rand_results(nn) = rng%exponential(1.0_dp)
end do
call cpu_time(time_end)
mean = sum(rand_results) / n_samples
variance = sum((rand_results - mean)**2) / n_samples
print *, ""
print *, "For exponential random numbers, the result is:"
print *, "nanoseconds per number (upper bound)", 1.0e9_dp * (time_end - time_start) / n_samples
print *, "mean/<mean>", mean
print *, "std dev/<std dev>", sqrt(variance)
end program test_m_random