jeudi 14 janvier 2016

Random number generator in PGI Fortran not so random

The following code just generates a simple triple of random numbers:

program testrand

integer, parameter :: nz = 160, nf = 160, nlt = 90

real :: tmpidx(3)
integer :: idxarr(3), idx1, idx2, idx3, seed_size, ticks
integer, allocatable :: seed(:)

call random_seed(size=seed_size)
allocate(seed(seed_size))

call system_clock(count=ticks)
seed = ticks+37*(/(i-1, i=1,seed_size)/)
call random_seed(put=seed)

deallocate(seed)

call random_number(tmpidx)
idxarr = tmpidx * (/nz, nf, nlt/)
idx1 = max(1,idxarr(1))
idx2 = max(1,idxarr(2))
idx3 = max(1,idxarr(3))
print *,idx1, idx2, idx3

end program

I compile this with gfortran and run a few times and I get:

> gfortran testrand.f90
> ./a.out
          74          98          86
> ./a.out
         113           3          10
> ./a.out
          44         104          27

Looks pretty random. Now I compile with PGI Fortran and run a few times:

> pgf90 testrand.f90
> ./a.out
            1            1            1
> ./a.out
            1            1            1
> ./a.out
            1            1            1

Of course, there's no way to be completely sure, but I suspect this is not random. :) Anyone know what is going on here? Anyone know the right way to get random numbers with PGI Fortran?




Aucun commentaire:

Enregistrer un commentaire