From b152f5a2b33b251ab1874a43d97ce73d11eec0a4 Mon Sep 17 00:00:00 2001 From: Janne Blomqvist Date: Thu, 11 Aug 2016 11:58:55 +0300 Subject: Replace KISS PRNG with xorshift1024* using per-thread state. frontend: 2016-08-11 Janne Blomqvist * check.c (gfc_check_random_seed): Use new seed size in check. * intrinsic.texi (RANDOM_NUMBER): Updated documentation. (RANDOM_SEED): Likewise. testsuite: 2016-08-11 Janne Blomqvist * gfortran.dg/random_7.f90: Take into account that the last seed value is the special p value. * gfortran.dg/random_seed_1.f90: Seed size is now constant. libgfortran: 2016-08-11 Janne Blomqvist * intrinsics/random.c: Replace KISS with xorshift1024* using per-thread state. * runtime/main.c (init): Don't call random_seed_i4. From-SVN: r239356 --- gcc/fortran/ChangeLog | 6 ++ gcc/fortran/check.c | 20 +++--- gcc/fortran/intrinsic.texi | 99 +++++++++-------------------- gcc/testsuite/ChangeLog | 6 ++ gcc/testsuite/gfortran.dg/random_7.f90 | 4 +- gcc/testsuite/gfortran.dg/random_seed_1.f90 | 18 +----- 6 files changed, 57 insertions(+), 96 deletions(-) (limited to 'gcc') diff --git a/gcc/fortran/ChangeLog b/gcc/fortran/ChangeLog index ef8fc17..4b9b739 100644 --- a/gcc/fortran/ChangeLog +++ b/gcc/fortran/ChangeLog @@ -1,3 +1,9 @@ +2016-08-11 Janne Blomqvist + + * check.c (gfc_check_random_seed): Use new seed size in check. + * intrinsic.texi (RANDOM_NUMBER): Updated documentation. + (RANDOM_SEED): Likewise. + 2016-08-08 Jakub Jelinek PR fortran/72716 diff --git a/gcc/fortran/check.c b/gcc/fortran/check.c index 288957a..ff5e80b 100644 --- a/gcc/fortran/check.c +++ b/gcc/fortran/check.c @@ -5527,16 +5527,14 @@ gfc_check_random_number (gfc_expr *harvest) bool gfc_check_random_seed (gfc_expr *size, gfc_expr *put, gfc_expr *get) { - unsigned int nargs = 0, kiss_size; + unsigned int nargs = 0, seed_size; locus *where = NULL; mpz_t put_size, get_size; - bool have_gfc_real_16; /* Try and mimic HAVE_GFC_REAL_16 in libgfortran. */ - have_gfc_real_16 = gfc_validate_kind (BT_REAL, 16, true) != -1; - - /* Keep the number of bytes in sync with kiss_size in - libgfortran/intrinsics/random.c. */ - kiss_size = (have_gfc_real_16 ? 48 : 32) / gfc_default_integer_kind; + /* Keep the number of bytes in sync with master_state in + libgfortran/intrinsics/random.c. +1 due to the integer p which is + part of the state too. */ + seed_size = 128 / gfc_default_integer_kind + 1; if (size != NULL) { @@ -5579,11 +5577,11 @@ gfc_check_random_seed (gfc_expr *size, gfc_expr *put, gfc_expr *get) return false; if (gfc_array_size (put, &put_size) - && mpz_get_ui (put_size) < kiss_size) + && mpz_get_ui (put_size) < seed_size) gfc_error ("Size of %qs argument of %qs intrinsic at %L " "too small (%i/%i)", gfc_current_intrinsic_arg[1]->name, gfc_current_intrinsic, - where, (int) mpz_get_ui (put_size), kiss_size); + where, (int) mpz_get_ui (put_size), seed_size); } if (get != NULL) @@ -5611,11 +5609,11 @@ gfc_check_random_seed (gfc_expr *size, gfc_expr *put, gfc_expr *get) return false; if (gfc_array_size (get, &get_size) - && mpz_get_ui (get_size) < kiss_size) + && mpz_get_ui (get_size) < seed_size) gfc_error ("Size of %qs argument of %qs intrinsic at %L " "too small (%i/%i)", gfc_current_intrinsic_arg[2]->name, gfc_current_intrinsic, - where, (int) mpz_get_ui (get_size), kiss_size); + where, (int) mpz_get_ui (get_size), seed_size); } /* RANDOM_SEED may not have more than one non-optional argument. */ diff --git a/gcc/fortran/intrinsic.texi b/gcc/fortran/intrinsic.texi index ae5d814..67b0231 100644 --- a/gcc/fortran/intrinsic.texi +++ b/gcc/fortran/intrinsic.texi @@ -11126,23 +11126,16 @@ end program test_rand Returns a single pseudorandom number or an array of pseudorandom numbers from the uniform distribution over the range @math{ 0 \leq x < 1}. -The runtime-library implements George Marsaglia's KISS (Keep It Simple -Stupid) random number generator (RNG). This RNG combines: -@enumerate -@item The congruential generator @math{x(n) = 69069 \cdot x(n-1) + 1327217885} -with a period of @math{2^{32}}, -@item A 3-shift shift-register generator with a period of @math{2^{32} - 1}, -@item Two 16-bit multiply-with-carry generators with a period of -@math{597273182964842497 > 2^{59}}. -@end enumerate -The overall period exceeds @math{2^{123}}. +The runtime-library implements the xorshift1024* random number +generator (RNG). This generator has a period of @math{2^{1024} - 1}, +and when using multiple threads up to @math{2^{512}} threads can each +generate @math{2^{512}} random numbers before any aliasing occurs. + +Note that in a multi-threaded program (e.g. using OpenMP directives), +each thread will have its own random number state. For details of the +seeding procedure, see the documentation for the @code{RANDOM_SEED} +intrinsic. -Please note, this RNG is thread safe if used within OpenMP directives, -i.e., its state will be consistent while called from multiple threads. -However, the KISS generator does not create random numbers in parallel -from multiple sources, but in sequence from a single source. If an -OpenMP-enabled application heavily relies on random numbers, one should -consider employing a dedicated parallel random number generator instead. @item @emph{Standard}: Fortran 95 and later @@ -11184,12 +11177,23 @@ end program Restarts or queries the state of the pseudorandom number generator used by @code{RANDOM_NUMBER}. -If @code{RANDOM_SEED} is called without arguments, it is initialized -to a default state. The example below shows how to initialize the -random seed with a varying seed in order to ensure a different random -number sequence for each invocation of the program. Note that setting -any of the seed values to zero should be avoided as it can result in -poor quality random numbers being generated. +If @code{RANDOM_SEED} is called without arguments, it is seeded with +random data retrieved from the operating system. + +As an extension to the Fortran standard, the GFortran +@code{RANDOM_NUMBER} supports multiple threads. Each thread in a +multi-threaded program has its own seed. When @code{RANDOM_SEED} is +called either without arguments or with the @var{PUT} argument, the +given seed is copied into a master seed as well as the seed of the +current thread. When a new thread uses @code{RANDOM_NUMBER} for the +first time, the seed is copied from the master seed, and forwarded +@math{N * 2^{512}} steps to guarantee that the random stream does not +alias any other stream in the system, where @var{N} is the number of +threads that have used @code{RANDOM_NUMBER} so far during the program +execution. + +Note that setting any of the seed values to zero should be avoided as +it can result in poor quality random numbers being generated. @item @emph{Standard}: Fortran 95 and later @@ -11217,57 +11221,16 @@ the @var{SIZE} argument. @item @emph{Example}: @smallexample -subroutine init_random_seed() - use iso_fortran_env, only: int64 +program test_random_seed implicit none integer, allocatable :: seed(:) - integer :: i, n, un, istat, dt(8), pid - integer(int64) :: t + integer :: n call random_seed(size = n) allocate(seed(n)) - ! First try if the OS provides a random number generator - open(newunit=un, file="/dev/urandom", access="stream", & - form="unformatted", action="read", status="old", iostat=istat) - if (istat == 0) then - read(un) seed - close(un) - else - ! Fallback to XOR:ing the current time and pid. The PID is - ! useful in case one launches multiple instances of the same - ! program in parallel. - 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) -contains - ! This simple PRNG might not be good enough for real work, but is - ! sufficient for seeding a better PRNG. - 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 + call random_seed(get=seed) + write (*, *) seed +end program test_random_seed @end smallexample @item @emph{See also}: diff --git a/gcc/testsuite/ChangeLog b/gcc/testsuite/ChangeLog index 6798ff8..7bca99a 100644 --- a/gcc/testsuite/ChangeLog +++ b/gcc/testsuite/ChangeLog @@ -1,3 +1,9 @@ +2016-08-11 Janne Blomqvist + + * gfortran.dg/random_7.f90: Take into account that the last seed + value is the special p value. + * gfortran.dg/random_seed_1.f90: Seed size is now constant. + 2016-08-11 Richard Biener * gcc.dg/tree-ssa/ssa-dom-thread-7.c: Adjust. diff --git a/gcc/testsuite/gfortran.dg/random_7.f90 b/gcc/testsuite/gfortran.dg/random_7.f90 index 6435a34..aafe346d 100644 --- a/gcc/testsuite/gfortran.dg/random_7.f90 +++ b/gcc/testsuite/gfortran.dg/random_7.f90 @@ -10,8 +10,8 @@ program trs seed(:) = huge(seed) / 17 call test_random_seed(put=seed) call test_random_seed(get=check) - print *, seed - print *, check + ! In the current implementation seed(17) is special + seed(17) = check(17) if (any (seed /= check)) call abort contains subroutine test_random_seed(size, put, get) diff --git a/gcc/testsuite/gfortran.dg/random_seed_1.f90 b/gcc/testsuite/gfortran.dg/random_seed_1.f90 index ccbcf00..39c81ce 100644 --- a/gcc/testsuite/gfortran.dg/random_seed_1.f90 +++ b/gcc/testsuite/gfortran.dg/random_seed_1.f90 @@ -3,10 +3,6 @@ ! Emit a diagnostic for too small PUT array at compile time ! See PR fortran/37159 -! Possible improvement: -! Provide a separate testcase for systems that support REAL(16), -! to test the minimum size of 12 (instead of 8). -! ! Updated to check for arrays of unexpected size, ! this also works for -fdefault-integer-8. ! @@ -14,19 +10,11 @@ PROGRAM random_seed_1 IMPLICIT NONE - ! Find out what the's largest kind size - INTEGER, PARAMETER :: k1 = kind (0.d0) - INTEGER, PARAMETER :: & - k2 = max (k1, selected_real_kind (precision (0._k1) + 1)) - INTEGER, PARAMETER :: & - k3 = max (k2, selected_real_kind (precision (0._k2) + 1)) - INTEGER, PARAMETER :: & - k4 = max (k3, selected_real_kind (precision (0._k3) + 1)) - - INTEGER, PARAMETER :: nbytes = MERGE(48, 32, k4 == 16) + INTEGER, PARAMETER :: nbytes = 128 + ! +1 due to the special 'p' value in xorshift1024* ! '+1' to avoid out-of-bounds warnings - INTEGER, PARAMETER :: n = nbytes / KIND(n) + 1 + INTEGER, PARAMETER :: n = nbytes / KIND(n) + 2 INTEGER, DIMENSION(n) :: seed ! Get seed, array too small -- cgit v1.1