I've been playing around with various random seed generators. Here is a simple one:
subroutine init_random_seed()
integer :: i, n, clock
integer, dimension(:), allocatable :: seed
call random_seed(size = n)
allocate(seed(n))
call system_clock(count=clock)
seed = clock + 37 * (/ (i - 1, i = 1, n) /)
call random_seed(put = seed)
deallocate(seed)
end
...and a more robust one:
SUBROUTINE init_random_seed()
USE ISO_Fortran_env, ONLY: INT64
IMPLICIT NONE
INTEGER, ALLOCATABLE :: seed(:)
INTEGER :: i, n, un, istat, dt(8), pid
INTEGER(INT64) :: t
CALL RANDOM_SEED(size = n)
ALLOCATE(seed(n))
OPEN(newunit=un, file='/dev/urandom', access='stream', status='old', action='read', form='unformatted', iostat=istat)
IF (istat == 0) THEN
READ(un) seed
CLOSE(un)
ELSE
CALL SYSTEM_CLOCK(t)
IF (t == 0) THEN
CALL DATE_AND_TIME(values = dt)
t = (dt(1) - 1970) * 365_INT64 * 24 * 60 * 60 * 1000 + dt(2) * 31_INT64 * 24 * 60 * 60 * 1000 + dt(3) * 24_INT64 * 60 * 60 &
* 1000 + dt(5) * 60 * 60 * 1000 + dt(6) * 60 * 1000 + dt(7) * 1000 + dt(8)
END IF
pid = GETPID()
t = IEOR(t, INT(pid, KIND(t)))
DO i = 1, n
seed(i) = lcg(t)
END DO
END IF
CALL RANDOM_SEED(put = seed)
DEALLOCATE(seed)
CONTAINS
FUNCTION lcg(s)
INTEGER :: lcg
INTEGER(INT64) :: s
IF (s == 0) THEN
s = 104729
ELSE
s = MOD(s, 4294967296_INT64)
END IF
s = MOD(s * 279470273_INT64, 4294967291_INT64)
lcg = INT(MOD(s, INT(HUGE(0), INT64)), KIND(0))
END FUNCTION lcg
END SUBROUTINE init_random_seed
The second one generates higher-quality random numbers, but is comparatively slow. Does anyone see why?
Aucun commentaire:
Enregistrer un commentaire