1 1.1 mrg /* Implementation of the IRAND, RAND, and SRAND intrinsics. 2 1.1.1.4 mrg Copyright (C) 2004-2024 Free Software Foundation, Inc. 3 1.1 mrg Contributed by Steven G. Kargl <kargls (at) comcast.net>. 4 1.1 mrg 5 1.1 mrg This file is part of the GNU Fortran 95 runtime library (libgfortran). 6 1.1 mrg 7 1.1 mrg Libgfortran is free software; you can redistribute it and/or 8 1.1 mrg modify it under the terms of the GNU General Public 9 1.1 mrg License as published by the Free Software Foundation; either 10 1.1 mrg version 3 of the License, or (at your option) any later version. 11 1.1 mrg 12 1.1 mrg Libgfortran is distributed in the hope that it will be useful, 13 1.1 mrg but WITHOUT ANY WARRANTY; without even the implied warranty of 14 1.1 mrg MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 15 1.1 mrg GNU General Public License for more details. 16 1.1 mrg 17 1.1 mrg Under Section 7 of GPL version 3, you are granted additional 18 1.1 mrg permissions described in the GCC Runtime Library Exception, version 19 1.1 mrg 3.1, as published by the Free Software Foundation. 20 1.1 mrg 21 1.1 mrg You should have received a copy of the GNU General Public License and 22 1.1 mrg a copy of the GCC Runtime Library Exception along with this program; 23 1.1 mrg see the files COPYING3 and COPYING.RUNTIME respectively. If not, see 24 1.1 mrg <http://www.gnu.org/licenses/>. */ 25 1.1 mrg 26 1.1 mrg /* Simple multiplicative congruent algorithm. 27 1.1 mrg The period of this generator is approximately 2^31-1, which means that 28 1.1 mrg it should not be used for anything serious. The implementation here 29 1.1 mrg is based of an algorithm from S.K. Park and K.W. Miller, Comm. ACM, 30 1.1 mrg 31, 1192-1201 (1988). It is also provided solely for compatibility 31 1.1 mrg with G77. */ 32 1.1 mrg 33 1.1 mrg #include "libgfortran.h" 34 1.1 mrg #include <gthr.h> 35 1.1 mrg 36 1.1 mrg #define GFC_RAND_A 16807 37 1.1 mrg #define GFC_RAND_M 2147483647 38 1.1 mrg #define GFC_RAND_M1 (GFC_RAND_M - 1) 39 1.1 mrg 40 1.1 mrg static GFC_UINTEGER_8 rand_seed = 1; 41 1.1 mrg #ifdef __GTHREAD_MUTEX_INIT 42 1.1 mrg static __gthread_mutex_t rand_seed_lock = __GTHREAD_MUTEX_INIT; 43 1.1 mrg #else 44 1.1 mrg static __gthread_mutex_t rand_seed_lock; 45 1.1 mrg #endif 46 1.1 mrg 47 1.1 mrg 48 1.1 mrg /* Set the seed of the irand generator. Note 0 is a bad seed. */ 49 1.1 mrg 50 1.1 mrg static void 51 1.1 mrg srand_internal (GFC_INTEGER_8 i) 52 1.1 mrg { 53 1.1 mrg rand_seed = i ? i : 123459876; 54 1.1 mrg } 55 1.1 mrg 56 1.1 mrg extern void PREFIX(srand) (GFC_INTEGER_4 *i); 57 1.1 mrg export_proto_np(PREFIX(srand)); 58 1.1 mrg 59 1.1 mrg void 60 1.1 mrg PREFIX(srand) (GFC_INTEGER_4 *i) 61 1.1 mrg { 62 1.1 mrg __gthread_mutex_lock (&rand_seed_lock); 63 1.1 mrg srand_internal (*i); 64 1.1 mrg __gthread_mutex_unlock (&rand_seed_lock); 65 1.1 mrg } 66 1.1 mrg 67 1.1 mrg /* Return an INTEGER in the range [1,GFC_RAND_M-1]. */ 68 1.1 mrg 69 1.1 mrg extern GFC_INTEGER_4 irand (GFC_INTEGER_4 *); 70 1.1 mrg iexport_proto(irand); 71 1.1 mrg 72 1.1 mrg GFC_INTEGER_4 73 1.1 mrg irand (GFC_INTEGER_4 *i) 74 1.1 mrg { 75 1.1 mrg GFC_INTEGER_4 j; 76 1.1 mrg if (i) 77 1.1 mrg j = *i; 78 1.1 mrg else 79 1.1 mrg j = 0; 80 1.1 mrg 81 1.1 mrg __gthread_mutex_lock (&rand_seed_lock); 82 1.1 mrg 83 1.1 mrg switch (j) 84 1.1 mrg { 85 1.1 mrg /* Return the next RN. */ 86 1.1 mrg case 0: 87 1.1 mrg break; 88 1.1 mrg 89 1.1 mrg /* Reset the RN sequence to system-dependent sequence and return the 90 1.1 mrg first value. */ 91 1.1 mrg case 1: 92 1.1 mrg srand_internal (0); 93 1.1 mrg break; 94 1.1 mrg 95 1.1 mrg /* Seed the RN sequence with j and return the first value. */ 96 1.1 mrg default: 97 1.1 mrg srand_internal (j); 98 1.1 mrg break; 99 1.1 mrg } 100 1.1 mrg 101 1.1 mrg rand_seed = GFC_RAND_A * rand_seed % GFC_RAND_M; 102 1.1 mrg j = (GFC_INTEGER_4) rand_seed; 103 1.1 mrg 104 1.1 mrg __gthread_mutex_unlock (&rand_seed_lock); 105 1.1 mrg 106 1.1 mrg return j; 107 1.1 mrg } 108 1.1 mrg iexport(irand); 109 1.1 mrg 110 1.1 mrg 111 1.1 mrg /* Return a random REAL in the range [0,1). */ 112 1.1 mrg 113 1.1 mrg extern GFC_REAL_4 PREFIX(rand) (GFC_INTEGER_4 *i); 114 1.1 mrg export_proto_np(PREFIX(rand)); 115 1.1 mrg 116 1.1 mrg GFC_REAL_4 117 1.1 mrg PREFIX(rand) (GFC_INTEGER_4 *i) 118 1.1 mrg { 119 1.1 mrg GFC_UINTEGER_4 mask; 120 1.1 mrg #if GFC_REAL_4_RADIX == 2 121 1.1 mrg mask = ~ (GFC_UINTEGER_4) 0u << (32 - GFC_REAL_4_DIGITS + 1); 122 1.1 mrg #elif GFC_REAL_4_RADIX == 16 123 1.1 mrg mask = ~ (GFC_UINTEGER_4) 0u << ((8 - GFC_REAL_4_DIGITS) * 4 + 1); 124 1.1 mrg #else 125 1.1 mrg #error "GFC_REAL_4_RADIX has unknown value" 126 1.1 mrg #endif 127 1.1 mrg return ((GFC_UINTEGER_4) (irand(i) -1) & mask) * (GFC_REAL_4) 0x1.p-31f; 128 1.1 mrg } 129 1.1 mrg 130 1.1 mrg #ifndef __GTHREAD_MUTEX_INIT 131 1.1 mrg static void __attribute__((constructor)) 132 1.1 mrg init (void) 133 1.1 mrg { 134 1.1 mrg __GTHREAD_MUTEX_INIT_FUNCTION (&rand_seed_lock); 135 1.1 mrg } 136 1.1 mrg #endif 137