1 | /* |
---|
2 | * Copyright (c) 1993 Martin Birgmeier |
---|
3 | * All rights reserved. |
---|
4 | * |
---|
5 | * You may redistribute unmodified or modified versions of this source |
---|
6 | * code provided that the above copyright notice and this and the |
---|
7 | * following conditions are retained. |
---|
8 | * |
---|
9 | * This software is provided ``as is'', and comes with no warranties |
---|
10 | * of any kind. I shall in no event be liable for anything that happens |
---|
11 | * to anyone/anything when using this software. |
---|
12 | */ |
---|
13 | |
---|
14 | /* |
---|
15 | FUNCTION |
---|
16 | <<rand48>>, <<drand48>>, <<erand48>>, <<lrand48>>, <<nrand48>>, <<mrand48>>, <<jrand48>>, <<srand48>>, <<seed48>>, <<lcong48>>---pseudo-random number generators and initialization routines |
---|
17 | |
---|
18 | INDEX |
---|
19 | rand48 |
---|
20 | INDEX |
---|
21 | drand48 |
---|
22 | INDEX |
---|
23 | erand48 |
---|
24 | INDEX |
---|
25 | lrand48 |
---|
26 | INDEX |
---|
27 | nrand48 |
---|
28 | INDEX |
---|
29 | mrand48 |
---|
30 | INDEX |
---|
31 | jrand48 |
---|
32 | INDEX |
---|
33 | srand48 |
---|
34 | INDEX |
---|
35 | seed48 |
---|
36 | INDEX |
---|
37 | lcong48 |
---|
38 | |
---|
39 | SYNOPSIS |
---|
40 | #include <stdlib.h> |
---|
41 | double drand48(void); |
---|
42 | double erand48(unsigned short <[xseed]>[3]); |
---|
43 | long lrand48(void); |
---|
44 | long nrand48(unsigned short <[xseed]>[3]); |
---|
45 | long mrand48(void); |
---|
46 | long jrand48(unsigned short <[xseed]>[3]); |
---|
47 | void srand48(long <[seed]>); |
---|
48 | unsigned short *seed48(unsigned short <[xseed]>[3]); |
---|
49 | void lcong48(unsigned short <[p]>[7]); |
---|
50 | |
---|
51 | DESCRIPTION |
---|
52 | The <<rand48>> family of functions generates pseudo-random numbers |
---|
53 | using a linear congruential algorithm working on integers 48 bits in size. |
---|
54 | The particular formula employed is |
---|
55 | r(n+1) = (a * r(n) + c) mod m |
---|
56 | where the default values are |
---|
57 | for the multiplicand a = 0xfdeece66d = 25214903917 and |
---|
58 | the addend c = 0xb = 11. The modulo is always fixed at m = 2 ** 48. |
---|
59 | r(n) is called the seed of the random number generator. |
---|
60 | |
---|
61 | For all the six generator routines described next, the first |
---|
62 | computational step is to perform a single iteration of the algorithm. |
---|
63 | |
---|
64 | <<drand48>> and <<erand48>> |
---|
65 | return values of type double. The full 48 bits of r(n+1) are |
---|
66 | loaded into the mantissa of the returned value, with the exponent set |
---|
67 | such that the values produced lie in the interval [0.0, 1.0]. |
---|
68 | |
---|
69 | <<lrand48>> and <<nrand48>> |
---|
70 | return values of type long in the range |
---|
71 | [0, 2**31-1]. The high-order (31) bits of |
---|
72 | r(n+1) are loaded into the lower bits of the returned value, with |
---|
73 | the topmost (sign) bit set to zero. |
---|
74 | |
---|
75 | <<mrand48>> and <<jrand48>> |
---|
76 | return values of type long in the range |
---|
77 | [-2**31, 2**31-1]. The high-order (32) bits of |
---|
78 | r(n+1) are loaded into the returned value. |
---|
79 | |
---|
80 | <<drand48>>, <<lrand48>>, and <<mrand48>> |
---|
81 | use an internal buffer to store r(n). For these functions |
---|
82 | the initial value of r(0) = 0x1234abcd330e = 20017429951246. |
---|
83 | |
---|
84 | On the other hand, <<erand48>>, <<nrand48>>, and <<jrand48>> |
---|
85 | use a user-supplied buffer to store the seed r(n), |
---|
86 | which consists of an array of 3 shorts, where the zeroth member |
---|
87 | holds the least significant bits. |
---|
88 | |
---|
89 | All functions share the same multiplicand and addend. |
---|
90 | |
---|
91 | <<srand48>> is used to initialize the internal buffer r(n) of |
---|
92 | <<drand48>>, <<lrand48>>, and <<mrand48>> |
---|
93 | such that the 32 bits of the seed value are copied into the upper 32 bits |
---|
94 | of r(n), with the lower 16 bits of r(n) arbitrarily being set to 0x330e. |
---|
95 | Additionally, the constant multiplicand and addend of the algorithm are |
---|
96 | reset to the default values given above. |
---|
97 | |
---|
98 | <<seed48>> also initializes the internal buffer r(n) of |
---|
99 | <<drand48>>, <<lrand48>>, and <<mrand48>>, |
---|
100 | but here all 48 bits of the seed can be specified in an array of 3 shorts, |
---|
101 | where the zeroth member specifies the lowest bits. Again, |
---|
102 | the constant multiplicand and addend of the algorithm are |
---|
103 | reset to the default values given above. |
---|
104 | <<seed48>> returns a pointer to an array of 3 shorts which contains |
---|
105 | the old seed. |
---|
106 | This array is statically allocated, thus its contents are lost after |
---|
107 | each new call to <<seed48>>. |
---|
108 | |
---|
109 | Finally, <<lcong48>> allows full control over the multiplicand and |
---|
110 | addend used in <<drand48>>, <<erand48>>, <<lrand48>>, <<nrand48>>, |
---|
111 | <<mrand48>>, and <<jrand48>>, |
---|
112 | and the seed used in <<drand48>>, <<lrand48>>, and <<mrand48>>. |
---|
113 | An array of 7 shorts is passed as parameter; the first three shorts are |
---|
114 | used to initialize the seed; the second three are used to initialize the |
---|
115 | multiplicand; and the last short is used to initialize the addend. |
---|
116 | It is thus not possible to use values greater than 0xffff as the addend. |
---|
117 | |
---|
118 | Note that all three methods of seeding the random number generator |
---|
119 | always also set the multiplicand and addend for any of the six |
---|
120 | generator calls. |
---|
121 | |
---|
122 | For a more powerful random number generator, see <<random>>. |
---|
123 | |
---|
124 | PORTABILITY |
---|
125 | SUS requires these functions. |
---|
126 | |
---|
127 | No supporting OS subroutines are required. |
---|
128 | */ |
---|
129 | |
---|
130 | #include "rand48.h" |
---|
131 | |
---|
132 | void |
---|
133 | __dorand48 (struct _reent *r, |
---|
134 | unsigned short xseed[3]) |
---|
135 | { |
---|
136 | unsigned long accu; |
---|
137 | unsigned short temp[2]; |
---|
138 | |
---|
139 | _REENT_CHECK_RAND48(r); |
---|
140 | accu = (unsigned long) __rand48_mult[0] * (unsigned long) xseed[0] + |
---|
141 | (unsigned long) __rand48_add; |
---|
142 | temp[0] = (unsigned short) accu; /* lower 16 bits */ |
---|
143 | accu >>= sizeof(unsigned short) * 8; |
---|
144 | accu += (unsigned long) __rand48_mult[0] * (unsigned long) xseed[1] + |
---|
145 | (unsigned long) __rand48_mult[1] * (unsigned long) xseed[0]; |
---|
146 | temp[1] = (unsigned short) accu; /* middle 16 bits */ |
---|
147 | accu >>= sizeof(unsigned short) * 8; |
---|
148 | accu += __rand48_mult[0] * xseed[2] + __rand48_mult[1] * xseed[1] + __rand48_mult[2] * xseed[0]; |
---|
149 | xseed[0] = temp[0]; |
---|
150 | xseed[1] = temp[1]; |
---|
151 | xseed[2] = (unsigned short) accu; |
---|
152 | } |
---|