1 /* Implementation of the RANDOM intrinsics
2 Copyright (C) 2002-2018 Free Software Foundation, Inc.
3 Contributed by Lars Segerlund <seger@linuxmail.org>,
4 Steve Kargl and Janne Blomqvist.
6 This file is part of the GNU Fortran runtime library (libgfortran).
8 Libgfortran is free software; you can redistribute it and/or
9 modify it under the terms of the GNU General Public
10 License as published by the Free Software Foundation; either
11 version 3 of the License, or (at your option) any later version.
13 Ligbfortran is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
18 Under Section 7 of GPL version 3, you are granted additional
19 permissions described in the GCC Runtime Library Exception, version
20 3.1, as published by the Free Software Foundation.
22 You should have received a copy of the GNU General Public License and
23 a copy of the GCC Runtime Library Exception along with this program;
24 see the files COPYING3 and COPYING.RUNTIME respectively. If not, see
25 <http://www.gnu.org/licenses/>. */
30 #include "libgfortran.h"
44 #include <_mingw.h> /* For __MINGW64_VERSION_MAJOR */
47 extern void random_r4 (GFC_REAL_4
*);
48 iexport_proto(random_r4
);
50 extern void random_r8 (GFC_REAL_8
*);
51 iexport_proto(random_r8
);
53 extern void arandom_r4 (gfc_array_r4
*);
54 export_proto(arandom_r4
);
56 extern void arandom_r8 (gfc_array_r8
*);
57 export_proto(arandom_r8
);
59 #ifdef HAVE_GFC_REAL_10
61 extern void random_r10 (GFC_REAL_10
*);
62 iexport_proto(random_r10
);
64 extern void arandom_r10 (gfc_array_r10
*);
65 export_proto(arandom_r10
);
69 #ifdef HAVE_GFC_REAL_16
71 extern void random_r16 (GFC_REAL_16
*);
72 iexport_proto(random_r16
);
74 extern void arandom_r16 (gfc_array_r16
*);
75 export_proto(arandom_r16
);
79 #ifdef __GTHREAD_MUTEX_INIT
80 static __gthread_mutex_t random_lock
= __GTHREAD_MUTEX_INIT
;
82 static __gthread_mutex_t random_lock
;
85 /* Helper routines to map a GFC_UINTEGER_* to the corresponding
86 GFC_REAL_* types in the range of [0,1). If GFC_REAL_*_RADIX are 2
87 or 16, respectively, we mask off the bits that don't fit into the
88 correct GFC_REAL_*, convert to the real type, then multiply by the
93 rnumber_4 (GFC_REAL_4
*f
, GFC_UINTEGER_4 v
)
96 #if GFC_REAL_4_RADIX == 2
97 mask
= ~ (GFC_UINTEGER_4
) 0u << (32 - GFC_REAL_4_DIGITS
);
98 #elif GFC_REAL_4_RADIX == 16
99 mask
= ~ (GFC_UINTEGER_4
) 0u << ((8 - GFC_REAL_4_DIGITS
) * 4);
101 #error "GFC_REAL_4_RADIX has unknown value"
104 *f
= (GFC_REAL_4
) v
* GFC_REAL_4_LITERAL(0x1.p
-32);
108 rnumber_8 (GFC_REAL_8
*f
, GFC_UINTEGER_8 v
)
111 #if GFC_REAL_8_RADIX == 2
112 mask
= ~ (GFC_UINTEGER_8
) 0u << (64 - GFC_REAL_8_DIGITS
);
113 #elif GFC_REAL_8_RADIX == 16
114 mask
= ~ (GFC_UINTEGER_8
) 0u << (16 - GFC_REAL_8_DIGITS
) * 4);
116 #error "GFC_REAL_8_RADIX has unknown value"
119 *f
= (GFC_REAL_8
) v
* GFC_REAL_8_LITERAL(0x1.p
-64);
122 #ifdef HAVE_GFC_REAL_10
125 rnumber_10 (GFC_REAL_10
*f
, GFC_UINTEGER_8 v
)
128 #if GFC_REAL_10_RADIX == 2
129 mask
= ~ (GFC_UINTEGER_8
) 0u << (64 - GFC_REAL_10_DIGITS
);
130 #elif GFC_REAL_10_RADIX == 16
131 mask
= ~ (GFC_UINTEGER_10
) 0u << ((16 - GFC_REAL_10_DIGITS
) * 4);
133 #error "GFC_REAL_10_RADIX has unknown value"
136 *f
= (GFC_REAL_10
) v
* GFC_REAL_10_LITERAL(0x1.p
-64);
140 #ifdef HAVE_GFC_REAL_16
142 /* For REAL(KIND=16), we only need to mask off the lower bits. */
145 rnumber_16 (GFC_REAL_16
*f
, GFC_UINTEGER_8 v1
, GFC_UINTEGER_8 v2
)
148 #if GFC_REAL_16_RADIX == 2
149 mask
= ~ (GFC_UINTEGER_8
) 0u << (128 - GFC_REAL_16_DIGITS
);
150 #elif GFC_REAL_16_RADIX == 16
151 mask
= ~ (GFC_UINTEGER_8
) 0u << ((32 - GFC_REAL_16_DIGITS
) * 4);
153 #error "GFC_REAL_16_RADIX has unknown value"
156 *f
= (GFC_REAL_16
) v1
* GFC_REAL_16_LITERAL(0x1.p
-64)
157 + (GFC_REAL_16
) v2
* GFC_REAL_16_LITERAL(0x1.p
-128);
164 We use the xorshift1024* generator, a fast high-quality generator
167 - passes TestU1 without any failures
169 - provides a "jump" function making it easy to provide many
170 independent parallel streams.
172 - Long period of 2**1024 - 1
174 A description can be found at
176 http://vigna.di.unimi.it/ftp/papers/xorshift.pdf
180 http://arxiv.org/abs/1402.6246
182 The paper includes public domain source code which is the basis for
183 the implementation below.
192 xorshift1024star_state
;
195 /* master_init, njumps, and master_state are the only variables
196 protected by random_lock. */
197 static bool master_init
;
198 static unsigned njumps
; /* How many times we have jumped. */
199 static uint64_t master_state
[] = {
200 0xad63fa1ed3b55f36ULL
, 0xd94473e78978b497ULL
, 0xbc60592a98172477ULL
,
201 0xa3de7c6e81265301ULL
, 0x586640c5e785af27ULL
, 0x7a2a3f63b67ce5eaULL
,
202 0x9fde969f922d9b82ULL
, 0xe6fe34379b3f3822ULL
, 0x6c277eac3e99b6c2ULL
,
203 0x9197290ab0d3f069ULL
, 0xdb227302f6c25576ULL
, 0xee0209aee527fae9ULL
,
204 0x675666a793cd05b9ULL
, 0xd048c99fbc70c20fULL
, 0x775f8c3dba385ef5ULL
,
205 0x625288bc262faf33ULL
209 static __gthread_key_t rand_state_key
;
211 static xorshift1024star_state
*
212 get_rand_state (void)
214 /* For single threaded apps. */
215 static xorshift1024star_state rand_state
;
217 if (__gthread_active_p ())
219 void* p
= __gthread_getspecific (rand_state_key
);
222 p
= xcalloc (1, sizeof (xorshift1024star_state
));
223 __gthread_setspecific (rand_state_key
, p
);
233 xorshift1024star (xorshift1024star_state
* rs
)
236 const uint64_t s0
= rs
->s
[p
];
237 uint64_t s1
= rs
->s
[p
= (p
+ 1) & 15];
239 rs
->s
[p
] = s1
^ s0
^ (s1
>> 11) ^ (s0
>> 30);
241 return rs
->s
[p
] * UINT64_C(1181783497276652981);
245 /* This is the jump function for the generator. It is equivalent to
246 2^512 calls to xorshift1024star(); it can be used to generate 2^512
247 non-overlapping subsequences for parallel computations. */
250 jump (xorshift1024star_state
* rs
)
252 static const uint64_t JUMP
[] = {
253 0x84242f96eca9c41dULL
, 0xa3c65b8776f96855ULL
, 0x5b34a39f070b5837ULL
,
254 0x4489affce4f31a1eULL
, 0x2ffeeb0a48316f40ULL
, 0xdc2d9891fe68c022ULL
,
255 0x3659132bb12fea70ULL
, 0xaac17d8efa43cab8ULL
, 0xc4cb815590989b13ULL
,
256 0x5ee975283d71c93bULL
, 0x691548c86c1bd540ULL
, 0x7910c41d10a1e6a5ULL
,
257 0x0b5fc64563b3e2a8ULL
, 0x047f7684e9fc949dULL
, 0xb99181f2d8f685caULL
,
258 0x284600e3f30e38c3ULL
261 uint64_t t
[16] = { 0 };
262 for(size_t i
= 0; i
< sizeof JUMP
/ sizeof *JUMP
; i
++)
263 for(int b
= 0; b
< 64; b
++)
265 if (JUMP
[i
] & 1ULL << b
)
266 for(int j
= 0; j
< 16; j
++)
267 t
[j
] ^= rs
->s
[(j
+ rs
->p
) & 15];
268 xorshift1024star (rs
);
270 for(int j
= 0; j
< 16; j
++)
271 rs
->s
[(j
+ rs
->p
) & 15] = t
[j
];
275 /* Super-simple LCG generator used in getosrandom () if /dev/urandom
278 #define M 2147483647 /* 2^31 - 1 (A large prime number) */
279 #define A 16807 /* Prime root of M, passes statistical tests and produces a full cycle */
280 #define Q 127773 /* M / A (To avoid overflow on A * seed) */
281 #define R 2836 /* M % A (To avoid overflow on A * seed) */
283 __attribute__((unused
)) static uint32_t
284 lcg_parkmiller(uint32_t seed
)
286 uint32_t hi
= seed
/ Q
;
287 uint32_t lo
= seed
% Q
;
288 int32_t test
= A
* lo
- R
* hi
;
300 /* Get some random bytes from the operating system in order to seed
304 getosrandom (void *buf
, size_t buflen
)
306 /* rand_s is available in MinGW-w64 but not plain MinGW. */
307 #if defined(__MINGW64_VERSION_MAJOR)
308 unsigned int* b
= buf
;
309 for (size_t i
= 0; i
< buflen
/ sizeof (unsigned int); i
++)
314 TODO: When glibc adds a wrapper for the getrandom() system call
315 on Linux, one could use that.
317 TODO: One could use getentropy() on OpenBSD. */
318 int flags
= O_RDONLY
;
322 int fd
= open("/dev/urandom", flags
);
325 int res
= read(fd
, buf
, buflen
);
329 uint32_t seed
= 1234567890;
332 if (gf_gettime (&secs
, &usecs
) == 0)
338 pid_t pid
= getpid();
342 for (size_t i
= 0; i
< buflen
/ sizeof (uint32_t); i
++)
345 seed
= lcg_parkmiller (seed
);
348 #endif /* __MINGW64_VERSION_MAJOR */
352 /* Initialize the random number generator for the current thread,
353 using the master state and the number of times we must jump. */
356 init_rand_state (xorshift1024star_state
* rs
, const bool locked
)
359 __gthread_mutex_lock (&random_lock
);
362 getosrandom (master_state
, sizeof (master_state
));
366 memcpy (&rs
->s
, master_state
, sizeof (master_state
));
367 unsigned n
= njumps
++;
369 __gthread_mutex_unlock (&random_lock
);
370 for (unsigned i
= 0; i
< n
; i
++)
376 /* This function produces a REAL(4) value from the uniform distribution
380 random_r4 (GFC_REAL_4
*x
)
382 xorshift1024star_state
* rs
= get_rand_state();
384 if (unlikely (!rs
->init
))
385 init_rand_state (rs
, false);
386 uint64_t r
= xorshift1024star (rs
);
387 /* Take the higher bits, ensuring that a stream of real(4), real(8),
388 and real(10) will be identical (except for precision). */
389 uint32_t high
= (uint32_t) (r
>> 32);
394 /* This function produces a REAL(8) value from the uniform distribution
398 random_r8 (GFC_REAL_8
*x
)
401 xorshift1024star_state
* rs
= get_rand_state();
403 if (unlikely (!rs
->init
))
404 init_rand_state (rs
, false);
405 r
= xorshift1024star (rs
);
410 #ifdef HAVE_GFC_REAL_10
412 /* This function produces a REAL(10) value from the uniform distribution
416 random_r10 (GFC_REAL_10
*x
)
419 xorshift1024star_state
* rs
= get_rand_state();
421 if (unlikely (!rs
->init
))
422 init_rand_state (rs
, false);
423 r
= xorshift1024star (rs
);
430 /* This function produces a REAL(16) value from the uniform distribution
433 #ifdef HAVE_GFC_REAL_16
436 random_r16 (GFC_REAL_16
*x
)
438 GFC_UINTEGER_8 r1
, r2
;
439 xorshift1024star_state
* rs
= get_rand_state();
441 if (unlikely (!rs
->init
))
442 init_rand_state (rs
, false);
443 r1
= xorshift1024star (rs
);
444 r2
= xorshift1024star (rs
);
445 rnumber_16 (x
, r1
, r2
);
452 /* This function fills a REAL(4) array with values from the uniform
453 distribution with range [0,1). */
456 arandom_r4 (gfc_array_r4
*x
)
458 index_type count
[GFC_MAX_DIMENSIONS
];
459 index_type extent
[GFC_MAX_DIMENSIONS
];
460 index_type stride
[GFC_MAX_DIMENSIONS
];
464 xorshift1024star_state
* rs
= get_rand_state();
468 dim
= GFC_DESCRIPTOR_RANK (x
);
470 for (index_type n
= 0; n
< dim
; n
++)
473 stride
[n
] = GFC_DESCRIPTOR_STRIDE(x
,n
);
474 extent
[n
] = GFC_DESCRIPTOR_EXTENT(x
,n
);
481 if (unlikely (!rs
->init
))
482 init_rand_state (rs
, false);
486 /* random_r4 (dest); */
487 uint64_t r
= xorshift1024star (rs
);
488 uint32_t high
= (uint32_t) (r
>> 32);
489 rnumber_4 (dest
, high
);
491 /* Advance to the next element. */
494 /* Advance to the next source element. */
496 while (count
[n
] == extent
[n
])
498 /* When we get to the end of a dimension, reset it and increment
499 the next dimension. */
501 /* We could precalculate these products, but this is a less
502 frequently used path so probably not worth it. */
503 dest
-= stride
[n
] * extent
[n
];
519 /* This function fills a REAL(8) array with values from the uniform
520 distribution with range [0,1). */
523 arandom_r8 (gfc_array_r8
*x
)
525 index_type count
[GFC_MAX_DIMENSIONS
];
526 index_type extent
[GFC_MAX_DIMENSIONS
];
527 index_type stride
[GFC_MAX_DIMENSIONS
];
531 xorshift1024star_state
* rs
= get_rand_state();
535 dim
= GFC_DESCRIPTOR_RANK (x
);
537 for (index_type n
= 0; n
< dim
; n
++)
540 stride
[n
] = GFC_DESCRIPTOR_STRIDE(x
,n
);
541 extent
[n
] = GFC_DESCRIPTOR_EXTENT(x
,n
);
548 if (unlikely (!rs
->init
))
549 init_rand_state (rs
, false);
553 /* random_r8 (dest); */
554 uint64_t r
= xorshift1024star (rs
);
557 /* Advance to the next element. */
560 /* Advance to the next source element. */
562 while (count
[n
] == extent
[n
])
564 /* When we get to the end of a dimension, reset it and increment
565 the next dimension. */
567 /* We could precalculate these products, but this is a less
568 frequently used path so probably not worth it. */
569 dest
-= stride
[n
] * extent
[n
];
585 #ifdef HAVE_GFC_REAL_10
587 /* This function fills a REAL(10) array with values from the uniform
588 distribution with range [0,1). */
591 arandom_r10 (gfc_array_r10
*x
)
593 index_type count
[GFC_MAX_DIMENSIONS
];
594 index_type extent
[GFC_MAX_DIMENSIONS
];
595 index_type stride
[GFC_MAX_DIMENSIONS
];
599 xorshift1024star_state
* rs
= get_rand_state();
603 dim
= GFC_DESCRIPTOR_RANK (x
);
605 for (index_type n
= 0; n
< dim
; n
++)
608 stride
[n
] = GFC_DESCRIPTOR_STRIDE(x
,n
);
609 extent
[n
] = GFC_DESCRIPTOR_EXTENT(x
,n
);
616 if (unlikely (!rs
->init
))
617 init_rand_state (rs
, false);
621 /* random_r10 (dest); */
622 uint64_t r
= xorshift1024star (rs
);
623 rnumber_10 (dest
, r
);
625 /* Advance to the next element. */
628 /* Advance to the next source element. */
630 while (count
[n
] == extent
[n
])
632 /* When we get to the end of a dimension, reset it and increment
633 the next dimension. */
635 /* We could precalculate these products, but this is a less
636 frequently used path so probably not worth it. */
637 dest
-= stride
[n
] * extent
[n
];
655 #ifdef HAVE_GFC_REAL_16
657 /* This function fills a REAL(16) array with values from the uniform
658 distribution with range [0,1). */
661 arandom_r16 (gfc_array_r16
*x
)
663 index_type count
[GFC_MAX_DIMENSIONS
];
664 index_type extent
[GFC_MAX_DIMENSIONS
];
665 index_type stride
[GFC_MAX_DIMENSIONS
];
669 xorshift1024star_state
* rs
= get_rand_state();
673 dim
= GFC_DESCRIPTOR_RANK (x
);
675 for (index_type n
= 0; n
< dim
; n
++)
678 stride
[n
] = GFC_DESCRIPTOR_STRIDE(x
,n
);
679 extent
[n
] = GFC_DESCRIPTOR_EXTENT(x
,n
);
686 if (unlikely (!rs
->init
))
687 init_rand_state (rs
, false);
691 /* random_r16 (dest); */
692 uint64_t r1
= xorshift1024star (rs
);
693 uint64_t r2
= xorshift1024star (rs
);
694 rnumber_16 (dest
, r1
, r2
);
696 /* Advance to the next element. */
699 /* Advance to the next source element. */
701 while (count
[n
] == extent
[n
])
703 /* When we get to the end of a dimension, reset it and increment
704 the next dimension. */
706 /* We could precalculate these products, but this is a less
707 frequently used path so probably not worth it. */
708 dest
-= stride
[n
] * extent
[n
];
727 /* Number of elements in master_state array. */
728 #define SZU64 (sizeof (master_state) / sizeof (uint64_t))
731 /* Keys for scrambling the seed in order to avoid poor seeds. */
733 static const uint64_t xor_keys
[] = {
734 0xbd0c5b6e50c2df49ULL
, 0xd46061cd46e1df38ULL
, 0xbb4f4d4ed6103544ULL
,
735 0x114a583d0756ad39ULL
, 0x4b5ad8623d0aaab6ULL
, 0x3f2ed7afbe0c0f21ULL
,
736 0xdec83fd65f113445ULL
, 0x3824f8fbc4f10d24ULL
, 0x5d9025af05878911ULL
,
737 0x500bc46b540340e9ULL
, 0x8bd53298e0d00530ULL
, 0x57886e40a952e06aULL
,
738 0x926e76c88e31cdb6ULL
, 0xbd0724dac0a3a5f9ULL
, 0xc5c8981b858ab796ULL
,
739 0xbb12ab2694c2b32cULL
743 /* Since a XOR cipher is symmetric, we need only one routine, and we
744 can use it both for encryption and decryption. */
747 scramble_seed (uint64_t *dest
, const uint64_t *src
)
749 for (size_t i
= 0; i
< SZU64
; i
++)
750 dest
[i
] = src
[i
] ^ xor_keys
[i
];
754 /* random_seed is used to seed the PRNG with either a default
755 set of seeds or user specified set of seeds. random_seed
756 must be called with no argument or exactly one argument. */
759 random_seed_i4 (GFC_INTEGER_4
*size
, gfc_array_i4
*put
, gfc_array_i4
*get
)
761 uint64_t seed
[SZU64
];
762 #define SZ (sizeof (master_state) / sizeof (GFC_INTEGER_4))
764 /* Check that we only have one argument present. */
765 if ((size
? 1 : 0) + (put
? 1 : 0) + (get
? 1 : 0) > 1)
766 runtime_error ("RANDOM_SEED should have at most one argument present.");
771 xorshift1024star_state
* rs
= get_rand_state();
773 /* Return the seed to GET data. */
776 /* If the rank of the array is not 1, abort. */
777 if (GFC_DESCRIPTOR_RANK (get
) != 1)
778 runtime_error ("Array rank of GET is not 1.");
780 /* If the array is too small, abort. */
781 if (GFC_DESCRIPTOR_EXTENT(get
,0) < (index_type
) SZ
+ 1)
782 runtime_error ("Array size of GET is too small.");
785 init_rand_state (rs
, false);
787 /* Unscramble the seed. */
788 scramble_seed (seed
, rs
->s
);
790 /* Then copy it back to the user variable. */
791 for (size_t i
= 0; i
< SZ
; i
++)
792 memcpy (&(get
->base_addr
[(SZ
- 1 - i
) * GFC_DESCRIPTOR_STRIDE(get
,0)]),
793 (unsigned char*) seed
+ i
* sizeof(GFC_UINTEGER_4
),
794 sizeof(GFC_UINTEGER_4
));
796 /* Finally copy the value of p after the seed. */
797 get
->base_addr
[SZ
* GFC_DESCRIPTOR_STRIDE(get
, 0)] = rs
->p
;
802 __gthread_mutex_lock (&random_lock
);
804 /* From the standard: "If no argument is present, the processor assigns
805 a processor-dependent value to the seed." */
806 if (size
== NULL
&& put
== NULL
&& get
== NULL
)
809 init_rand_state (rs
, true);
814 /* If the rank of the array is not 1, abort. */
815 if (GFC_DESCRIPTOR_RANK (put
) != 1)
816 runtime_error ("Array rank of PUT is not 1.");
818 /* If the array is too small, abort. */
819 if (GFC_DESCRIPTOR_EXTENT(put
,0) < (index_type
) SZ
+ 1)
820 runtime_error ("Array size of PUT is too small.");
822 /* We copy the seed given by the user. */
823 for (size_t i
= 0; i
< SZ
; i
++)
824 memcpy ((unsigned char*) seed
+ i
* sizeof(GFC_UINTEGER_4
),
825 &(put
->base_addr
[(SZ
- 1 - i
) * GFC_DESCRIPTOR_STRIDE(put
,0)]),
826 sizeof(GFC_UINTEGER_4
));
828 /* We put it after scrambling the bytes, to paper around users who
829 provide seeds with quality only in the lower or upper part. */
830 scramble_seed (master_state
, seed
);
833 init_rand_state (rs
, true);
835 rs
->p
= put
->base_addr
[SZ
* GFC_DESCRIPTOR_STRIDE(put
, 0)] & 15;
838 __gthread_mutex_unlock (&random_lock
);
842 iexport(random_seed_i4
);
846 random_seed_i8 (GFC_INTEGER_8
*size
, gfc_array_i8
*put
, gfc_array_i8
*get
)
848 uint64_t seed
[SZU64
];
850 /* Check that we only have one argument present. */
851 if ((size
? 1 : 0) + (put
? 1 : 0) + (get
? 1 : 0) > 1)
852 runtime_error ("RANDOM_SEED should have at most one argument present.");
854 #define SZ (sizeof (master_state) / sizeof (GFC_INTEGER_8))
858 xorshift1024star_state
* rs
= get_rand_state();
860 /* Return the seed to GET data. */
863 /* If the rank of the array is not 1, abort. */
864 if (GFC_DESCRIPTOR_RANK (get
) != 1)
865 runtime_error ("Array rank of GET is not 1.");
867 /* If the array is too small, abort. */
868 if (GFC_DESCRIPTOR_EXTENT(get
,0) < (index_type
) SZ
+ 1)
869 runtime_error ("Array size of GET is too small.");
872 init_rand_state (rs
, false);
874 /* Unscramble the seed. */
875 scramble_seed (seed
, rs
->s
);
877 /* This code now should do correct strides. */
878 for (size_t i
= 0; i
< SZ
; i
++)
879 memcpy (&(get
->base_addr
[i
* GFC_DESCRIPTOR_STRIDE(get
,0)]), &seed
[i
],
880 sizeof (GFC_UINTEGER_8
));
882 get
->base_addr
[SZ
* GFC_DESCRIPTOR_STRIDE(get
, 0)] = rs
->p
;
887 __gthread_mutex_lock (&random_lock
);
889 /* From the standard: "If no argument is present, the processor assigns
890 a processor-dependent value to the seed." */
891 if (size
== NULL
&& put
== NULL
&& get
== NULL
)
894 init_rand_state (rs
, true);
899 /* If the rank of the array is not 1, abort. */
900 if (GFC_DESCRIPTOR_RANK (put
) != 1)
901 runtime_error ("Array rank of PUT is not 1.");
903 /* If the array is too small, abort. */
904 if (GFC_DESCRIPTOR_EXTENT(put
,0) < (index_type
) SZ
+ 1)
905 runtime_error ("Array size of PUT is too small.");
907 /* This code now should do correct strides. */
908 for (size_t i
= 0; i
< SZ
; i
++)
909 memcpy (&seed
[i
], &(put
->base_addr
[i
* GFC_DESCRIPTOR_STRIDE(put
,0)]),
910 sizeof (GFC_UINTEGER_8
));
912 scramble_seed (master_state
, seed
);
915 init_rand_state (rs
, true);
916 rs
->p
= put
->base_addr
[SZ
* GFC_DESCRIPTOR_STRIDE(put
, 0)] & 15;
920 __gthread_mutex_unlock (&random_lock
);
923 iexport(random_seed_i8
);
926 #if !defined __GTHREAD_MUTEX_INIT || defined __GTHREADS
927 static void __attribute__((constructor
))
928 constructor_random (void)
930 #ifndef __GTHREAD_MUTEX_INIT
931 __GTHREAD_MUTEX_INIT_FUNCTION (&random_lock
);
933 if (__gthread_active_p ())
934 __gthread_key_create (&rand_state_key
, &free
);
939 static void __attribute__((destructor
))
940 destructor_random (void)
942 if (__gthread_active_p ())
943 __gthread_key_delete (rand_state_key
);