|
111
|
1 ------------------------------------------------------------------------------
|
|
|
2 -- --
|
|
|
3 -- GNAT RUN-TIME COMPONENTS --
|
|
|
4 -- --
|
|
|
5 -- S Y S T E M . R A N D O M _ N U M B E R S --
|
|
|
6 -- --
|
|
|
7 -- B o d y --
|
|
|
8 -- --
|
|
131
|
9 -- Copyright (C) 2007-2018, Free Software Foundation, Inc. --
|
|
111
|
10 -- --
|
|
|
11 -- GNAT is free software; you can redistribute it and/or modify it under --
|
|
|
12 -- terms of the GNU General Public License as published by the Free Soft- --
|
|
|
13 -- ware Foundation; either version 3, or (at your option) any later ver- --
|
|
|
14 -- sion. GNAT is distributed in the hope that it will be useful, but WITH- --
|
|
|
15 -- OUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY --
|
|
|
16 -- or FITNESS FOR A PARTICULAR PURPOSE. --
|
|
|
17 -- --
|
|
|
18 -- As a special exception under Section 7 of GPL version 3, you are granted --
|
|
|
19 -- additional permissions described in the GCC Runtime Library Exception, --
|
|
|
20 -- version 3.1, as published by the Free Software Foundation. --
|
|
|
21 -- --
|
|
|
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/>. --
|
|
|
26 -- --
|
|
|
27 -- GNAT was originally developed by the GNAT team at New York University. --
|
|
|
28 -- Extensive contributions were provided by Ada Core Technologies Inc. --
|
|
|
29 -- --
|
|
|
30 ------------------------------------------------------------------------------
|
|
|
31
|
|
|
32 ------------------------------------------------------------------------------
|
|
|
33 -- --
|
|
|
34 -- The implementation here is derived from a C-program for MT19937, with --
|
|
|
35 -- initialization improved 2002/1/26. As required, the following notice is --
|
|
|
36 -- copied from the original program. --
|
|
|
37 -- --
|
|
|
38 -- Copyright (C) 1997 - 2002, Makoto Matsumoto and Takuji Nishimura, --
|
|
|
39 -- All rights reserved. --
|
|
|
40 -- --
|
|
|
41 -- Redistribution and use in source and binary forms, with or without --
|
|
|
42 -- modification, are permitted provided that the following conditions --
|
|
|
43 -- are met: --
|
|
|
44 -- --
|
|
|
45 -- 1. Redistributions of source code must retain the above copyright --
|
|
|
46 -- notice, this list of conditions and the following disclaimer. --
|
|
|
47 -- --
|
|
|
48 -- 2. Redistributions in binary form must reproduce the above copyright --
|
|
|
49 -- notice, this list of conditions and the following disclaimer in the --
|
|
|
50 -- documentation and/or other materials provided with the distribution.--
|
|
|
51 -- --
|
|
|
52 -- 3. The names of its contributors may not be used to endorse or promote --
|
|
|
53 -- products derived from this software without specific prior written --
|
|
|
54 -- permission. --
|
|
|
55 -- --
|
|
|
56 -- THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS --
|
|
|
57 -- "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT --
|
|
|
58 -- LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR --
|
|
|
59 -- A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT --
|
|
|
60 -- OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, --
|
|
|
61 -- SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED --
|
|
|
62 -- TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR --
|
|
|
63 -- PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF --
|
|
|
64 -- LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING --
|
|
|
65 -- NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS --
|
|
|
66 -- SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. --
|
|
|
67 -- --
|
|
|
68 ------------------------------------------------------------------------------
|
|
|
69
|
|
|
70 ------------------------------------------------------------------------------
|
|
|
71 -- --
|
|
|
72 -- This is an implementation of the Mersenne Twister, twisted generalized --
|
|
|
73 -- feedback shift register of rational normal form, with state-bit --
|
|
|
74 -- reflection and tempering. This version generates 32-bit integers with a --
|
|
|
75 -- period of 2**19937 - 1 (a Mersenne prime, hence the name). For --
|
|
|
76 -- applications requiring more than 32 bits (up to 64), we concatenate two --
|
|
|
77 -- 32-bit numbers. --
|
|
|
78 -- --
|
|
|
79 -- See http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html for --
|
|
|
80 -- details. --
|
|
|
81 -- --
|
|
|
82 -- In contrast to the original code, we do not generate random numbers in --
|
|
|
83 -- batches of N. Measurement seems to show this has very little if any --
|
|
|
84 -- effect on performance, and it may be marginally better for real-time --
|
|
|
85 -- applications with hard deadlines. --
|
|
|
86 -- --
|
|
|
87 ------------------------------------------------------------------------------
|
|
|
88
|
|
|
89 with Ada.Unchecked_Conversion;
|
|
|
90
|
|
|
91 with System.Random_Seed;
|
|
|
92
|
|
|
93 with Interfaces; use Interfaces;
|
|
|
94
|
|
|
95 use Ada;
|
|
|
96
|
|
|
97 package body System.Random_Numbers with
|
|
|
98 SPARK_Mode => Off
|
|
|
99 is
|
|
|
100 Image_Numeral_Length : constant := Max_Image_Width / N;
|
|
|
101
|
|
|
102 subtype Image_String is String (1 .. Max_Image_Width);
|
|
|
103
|
|
|
104 ----------------------------
|
|
|
105 -- Algorithmic Parameters --
|
|
|
106 ----------------------------
|
|
|
107
|
|
|
108 Lower_Mask : constant := 2**31 - 1;
|
|
|
109 Upper_Mask : constant := 2**31;
|
|
|
110
|
|
|
111 Matrix_A : constant array (State_Val range 0 .. 1) of State_Val
|
|
|
112 := (0, 16#9908b0df#);
|
|
|
113 -- The twist transformation is represented by a matrix of the form
|
|
|
114 --
|
|
|
115 -- [ 0 I(31) ]
|
|
|
116 -- [ _a ]
|
|
|
117 --
|
|
|
118 -- where 0 is a 31x31 block of 0s, I(31) is the 31x31 identity matrix and
|
|
|
119 -- _a is a particular bit row-vector, represented here by a 32-bit integer.
|
|
|
120 -- If integer x represents a row vector of bits (with x(0), the units bit,
|
|
|
121 -- last), then
|
|
|
122 -- x * A = [0 x(31..1)] xor Matrix_A(x(0)).
|
|
|
123
|
|
|
124 U : constant := 11;
|
|
|
125 S : constant := 7;
|
|
|
126 B_Mask : constant := 16#9d2c5680#;
|
|
|
127 T : constant := 15;
|
|
|
128 C_Mask : constant := 16#efc60000#;
|
|
|
129 L : constant := 18;
|
|
|
130 -- The tempering shifts and bit masks, in the order applied
|
|
|
131
|
|
|
132 Seed0 : constant := 5489;
|
|
|
133 -- Default seed, used to initialize the state vector when Reset not called
|
|
|
134
|
|
|
135 Seed1 : constant := 19650218;
|
|
|
136 -- Seed used to initialize the state vector when calling Reset with an
|
|
|
137 -- initialization vector.
|
|
|
138
|
|
|
139 Mult0 : constant := 1812433253;
|
|
|
140 -- Multiplier for a modified linear congruential generator used to
|
|
|
141 -- initialize the state vector when calling Reset with a single integer
|
|
|
142 -- seed.
|
|
|
143
|
|
|
144 Mult1 : constant := 1664525;
|
|
|
145 Mult2 : constant := 1566083941;
|
|
|
146 -- Multipliers for two modified linear congruential generators used to
|
|
|
147 -- initialize the state vector when calling Reset with an initialization
|
|
|
148 -- vector.
|
|
|
149
|
|
|
150 -----------------------
|
|
|
151 -- Local Subprograms --
|
|
|
152 -----------------------
|
|
|
153
|
|
|
154 procedure Init (Gen : Generator; Initiator : Unsigned_32);
|
|
|
155 -- Perform a default initialization of the state of Gen. The resulting
|
|
|
156 -- state is identical for identical values of Initiator.
|
|
|
157
|
|
|
158 procedure Insert_Image
|
|
|
159 (S : in out Image_String;
|
|
|
160 Index : Integer;
|
|
|
161 V : State_Val);
|
|
|
162 -- Insert image of V into S, in the Index'th 11-character substring
|
|
|
163
|
|
|
164 function Extract_Value (S : String; Index : Integer) return State_Val;
|
|
|
165 -- Treat S as a sequence of 11-character decimal numerals and return
|
|
|
166 -- the result of converting numeral #Index (numbering from 0)
|
|
|
167
|
|
|
168 function To_Unsigned is
|
|
|
169 new Unchecked_Conversion (Integer_32, Unsigned_32);
|
|
|
170 function To_Unsigned is
|
|
|
171 new Unchecked_Conversion (Integer_64, Unsigned_64);
|
|
|
172
|
|
|
173 ------------
|
|
|
174 -- Random --
|
|
|
175 ------------
|
|
|
176
|
|
|
177 function Random (Gen : Generator) return Unsigned_32 is
|
|
|
178 G : Generator renames Gen.Writable.Self.all;
|
|
|
179 Y : State_Val;
|
|
|
180 I : Integer; -- should avoid use of identifier I ???
|
|
|
181
|
|
|
182 begin
|
|
|
183 I := G.I;
|
|
|
184
|
|
|
185 if I < N - M then
|
|
|
186 Y := (G.S (I) and Upper_Mask) or (G.S (I + 1) and Lower_Mask);
|
|
|
187 Y := G.S (I + M) xor Shift_Right (Y, 1) xor Matrix_A (Y and 1);
|
|
|
188 I := I + 1;
|
|
|
189
|
|
|
190 elsif I < N - 1 then
|
|
|
191 Y := (G.S (I) and Upper_Mask) or (G.S (I + 1) and Lower_Mask);
|
|
|
192 Y := G.S (I + (M - N))
|
|
|
193 xor Shift_Right (Y, 1)
|
|
|
194 xor Matrix_A (Y and 1);
|
|
|
195 I := I + 1;
|
|
|
196
|
|
|
197 elsif I = N - 1 then
|
|
|
198 Y := (G.S (I) and Upper_Mask) or (G.S (0) and Lower_Mask);
|
|
|
199 Y := G.S (M - 1) xor Shift_Right (Y, 1) xor Matrix_A (Y and 1);
|
|
|
200 I := 0;
|
|
|
201
|
|
|
202 else
|
|
|
203 Init (G, Seed0);
|
|
|
204 return Random (Gen);
|
|
|
205 end if;
|
|
|
206
|
|
|
207 G.S (G.I) := Y;
|
|
|
208 G.I := I;
|
|
|
209
|
|
|
210 Y := Y xor Shift_Right (Y, U);
|
|
|
211 Y := Y xor (Shift_Left (Y, S) and B_Mask);
|
|
|
212 Y := Y xor (Shift_Left (Y, T) and C_Mask);
|
|
|
213 Y := Y xor Shift_Right (Y, L);
|
|
|
214
|
|
|
215 return Y;
|
|
|
216 end Random;
|
|
|
217
|
|
|
218 generic
|
|
|
219 type Unsigned is mod <>;
|
|
|
220 type Real is digits <>;
|
|
|
221 with function Random (G : Generator) return Unsigned is <>;
|
|
|
222 function Random_Float_Template (Gen : Generator) return Real;
|
|
|
223 pragma Inline (Random_Float_Template);
|
|
|
224 -- Template for a random-number generator implementation that delivers
|
|
|
225 -- values of type Real in the range [0 .. 1], using values from Gen,
|
|
|
226 -- assuming that Unsigned is large enough to hold the bits of a mantissa
|
|
|
227 -- for type Real.
|
|
|
228
|
|
|
229 ---------------------------
|
|
|
230 -- Random_Float_Template --
|
|
|
231 ---------------------------
|
|
|
232
|
|
|
233 function Random_Float_Template (Gen : Generator) return Real is
|
|
|
234
|
|
|
235 pragma Compile_Time_Error
|
|
|
236 (Unsigned'Last <= 2**(Real'Machine_Mantissa - 1),
|
|
|
237 "insufficiently large modular type used to hold mantissa");
|
|
|
238
|
|
|
239 begin
|
|
|
240 -- This code generates random floating-point numbers from unsigned
|
|
|
241 -- integers. Assuming that Real'Machine_Radix = 2, it can deliver all
|
|
|
242 -- machine values of type Real (as implied by Real'Machine_Mantissa and
|
|
|
243 -- Real'Machine_Emin), which is not true of the standard method (to
|
|
|
244 -- which we fall back for nonbinary radix): computing Real(<random
|
|
|
245 -- integer>) / (<max random integer>+1). To do so, we first extract an
|
|
|
246 -- (M-1)-bit significand (where M is Real'Machine_Mantissa), and then
|
|
|
247 -- decide on a normalized exponent by repeated coin flips, decrementing
|
|
|
248 -- from 0 as long as we flip heads (1 bits). This process yields the
|
|
|
249 -- proper geometric distribution for the exponent: in a uniformly
|
|
|
250 -- distributed set of floating-point numbers, 1/2 of them will be in
|
|
|
251 -- (0.5, 1], 1/4 will be in (0.25, 0.5], and so forth. It makes a
|
|
|
252 -- further adjustment at binade boundaries (see comments below) to give
|
|
|
253 -- the effect of selecting a uniformly distributed real deviate in
|
|
|
254 -- [0..1] and then rounding to the nearest representable floating-point
|
|
|
255 -- number. The algorithm attempts to be stingy with random integers. In
|
|
|
256 -- the worst case, it can consume roughly -Real'Machine_Emin/32 32-bit
|
|
|
257 -- integers, but this case occurs with probability around
|
|
|
258 -- 2**Machine_Emin, and the expected number of calls to integer-valued
|
|
|
259 -- Random is 1. For another discussion of the issues addressed by this
|
|
|
260 -- process, see Allen Downey's unpublished paper at
|
|
|
261 -- http://allendowney.com/research/rand/downey07randfloat.pdf.
|
|
|
262
|
|
|
263 if Real'Machine_Radix /= 2 then
|
|
|
264 return Real'Machine
|
|
|
265 (Real (Unsigned'(Random (Gen))) * 2.0**(-Unsigned'Size));
|
|
|
266
|
|
|
267 else
|
|
|
268 declare
|
|
|
269 type Bit_Count is range 0 .. 4;
|
|
|
270
|
|
|
271 subtype T is Real'Base;
|
|
|
272
|
|
|
273 Trailing_Ones : constant array (Unsigned_32 range 0 .. 15)
|
|
|
274 of Bit_Count :=
|
|
|
275 (2#00000# => 0, 2#00001# => 1, 2#00010# => 0, 2#00011# => 2,
|
|
|
276 2#00100# => 0, 2#00101# => 1, 2#00110# => 0, 2#00111# => 3,
|
|
|
277 2#01000# => 0, 2#01001# => 1, 2#01010# => 0, 2#01011# => 2,
|
|
|
278 2#01100# => 0, 2#01101# => 1, 2#01110# => 0, 2#01111# => 4);
|
|
|
279
|
|
|
280 Pow_Tab : constant array (Bit_Count range 0 .. 3) of Real
|
|
|
281 := (0 => 2.0**(0 - T'Machine_Mantissa),
|
|
|
282 1 => 2.0**(-1 - T'Machine_Mantissa),
|
|
|
283 2 => 2.0**(-2 - T'Machine_Mantissa),
|
|
|
284 3 => 2.0**(-3 - T'Machine_Mantissa));
|
|
|
285
|
|
|
286 Extra_Bits : constant Natural :=
|
|
|
287 (Unsigned'Size - T'Machine_Mantissa + 1);
|
|
|
288 -- Random bits left over after selecting mantissa
|
|
|
289
|
|
|
290 Mantissa : Unsigned;
|
|
|
291
|
|
|
292 X : Real; -- Scaled mantissa
|
|
|
293 R : Unsigned_32; -- Supply of random bits
|
|
|
294 R_Bits : Natural; -- Number of bits left in R
|
|
|
295 K : Bit_Count; -- Next decrement to exponent
|
|
|
296
|
|
|
297 begin
|
|
|
298 Mantissa := Random (Gen) / 2**Extra_Bits;
|
|
|
299 R := Unsigned_32 (Mantissa mod 2**Extra_Bits);
|
|
|
300 R_Bits := Extra_Bits;
|
|
|
301 X := Real (2**(T'Machine_Mantissa - 1) + Mantissa); -- Exact
|
|
|
302
|
|
|
303 if Extra_Bits < 4 and then R < 2 ** Extra_Bits - 1 then
|
|
|
304
|
|
|
305 -- We got lucky and got a zero in our few extra bits
|
|
|
306
|
|
|
307 K := Trailing_Ones (R);
|
|
|
308
|
|
|
309 else
|
|
|
310 Find_Zero : loop
|
|
|
311
|
|
|
312 -- R has R_Bits unprocessed random bits, a multiple of 4.
|
|
|
313 -- X needs to be halved for each trailing one bit. The
|
|
|
314 -- process stops as soon as a 0 bit is found. If R_Bits
|
|
|
315 -- becomes zero, reload R.
|
|
|
316
|
|
|
317 -- Process 4 bits at a time for speed: the two iterations
|
|
|
318 -- on average with three tests each was still too slow,
|
|
|
319 -- probably because the branches are not predictable.
|
|
|
320 -- This loop now will only execute once 94% of the cases,
|
|
|
321 -- doing more bits at a time will not help.
|
|
|
322
|
|
|
323 while R_Bits >= 4 loop
|
|
|
324 K := Trailing_Ones (R mod 16);
|
|
|
325
|
|
|
326 exit Find_Zero when K < 4; -- Exits 94% of the time
|
|
|
327
|
|
|
328 R_Bits := R_Bits - 4;
|
|
|
329 X := X / 16.0;
|
|
|
330 R := R / 16;
|
|
|
331 end loop;
|
|
|
332
|
|
|
333 -- Do not allow us to loop endlessly even in the (very
|
|
|
334 -- unlikely) case that Random (Gen) keeps yielding all ones.
|
|
|
335
|
|
|
336 exit Find_Zero when X = 0.0;
|
|
|
337 R := Random (Gen);
|
|
|
338 R_Bits := 32;
|
|
|
339 end loop Find_Zero;
|
|
|
340 end if;
|
|
|
341
|
|
|
342 -- K has the count of trailing ones not reflected yet in X. The
|
|
|
343 -- following multiplication takes care of that, as well as the
|
|
|
344 -- correction to move the radix point to the left of the mantissa.
|
|
|
345 -- Doing it at the end avoids repeated rounding errors in the
|
|
|
346 -- exceedingly unlikely case of ever having a subnormal result.
|
|
|
347
|
|
|
348 X := X * Pow_Tab (K);
|
|
|
349
|
|
|
350 -- The smallest value in each binade is rounded to by 0.75 of
|
|
|
351 -- the span of real numbers as its next larger neighbor, and
|
|
|
352 -- 1.0 is rounded to by half of the span of real numbers as its
|
|
|
353 -- next smaller neighbor. To account for this, when we encounter
|
|
|
354 -- the smallest number in a binade, we substitute the smallest
|
|
|
355 -- value in the next larger binade with probability 1/2.
|
|
|
356
|
|
|
357 if Mantissa = 0 and then Unsigned_32'(Random (Gen)) mod 2 = 0 then
|
|
|
358 X := 2.0 * X;
|
|
|
359 end if;
|
|
|
360
|
|
|
361 return X;
|
|
|
362 end;
|
|
|
363 end if;
|
|
|
364 end Random_Float_Template;
|
|
|
365
|
|
|
366 ------------
|
|
|
367 -- Random --
|
|
|
368 ------------
|
|
|
369
|
|
|
370 function Random (Gen : Generator) return Float is
|
|
|
371 function F is new Random_Float_Template (Unsigned_32, Float);
|
|
|
372 begin
|
|
|
373 return F (Gen);
|
|
|
374 end Random;
|
|
|
375
|
|
|
376 function Random (Gen : Generator) return Long_Float is
|
|
|
377 function F is new Random_Float_Template (Unsigned_64, Long_Float);
|
|
|
378 begin
|
|
|
379 return F (Gen);
|
|
|
380 end Random;
|
|
|
381
|
|
|
382 function Random (Gen : Generator) return Unsigned_64 is
|
|
|
383 begin
|
|
|
384 return Shift_Left (Unsigned_64 (Unsigned_32'(Random (Gen))), 32)
|
|
|
385 or Unsigned_64 (Unsigned_32'(Random (Gen)));
|
|
|
386 end Random;
|
|
|
387
|
|
|
388 ---------------------
|
|
|
389 -- Random_Discrete --
|
|
|
390 ---------------------
|
|
|
391
|
|
|
392 function Random_Discrete
|
|
|
393 (Gen : Generator;
|
|
|
394 Min : Result_Subtype := Default_Min;
|
|
|
395 Max : Result_Subtype := Result_Subtype'Last) return Result_Subtype
|
|
|
396 is
|
|
|
397 begin
|
|
|
398 if Max = Min then
|
|
|
399 return Max;
|
|
|
400
|
|
|
401 elsif Max < Min then
|
|
|
402 raise Constraint_Error;
|
|
|
403
|
|
|
404 elsif Result_Subtype'Base'Size > 32 then
|
|
|
405 declare
|
|
|
406 -- In the 64-bit case, we have to be careful, since not all 64-bit
|
|
|
407 -- unsigned values are representable in GNAT's root_integer type.
|
|
|
408 -- Ignore different-size warnings here since GNAT's handling
|
|
|
409 -- is correct.
|
|
|
410
|
|
|
411 pragma Warnings ("Z");
|
|
|
412 function Conv_To_Unsigned is
|
|
|
413 new Unchecked_Conversion (Result_Subtype'Base, Unsigned_64);
|
|
|
414 function Conv_To_Result is
|
|
|
415 new Unchecked_Conversion (Unsigned_64, Result_Subtype'Base);
|
|
|
416 pragma Warnings ("z");
|
|
|
417
|
|
|
418 N : constant Unsigned_64 :=
|
|
|
419 Conv_To_Unsigned (Max) - Conv_To_Unsigned (Min) + 1;
|
|
|
420
|
|
|
421 X, Slop : Unsigned_64;
|
|
|
422
|
|
|
423 begin
|
|
|
424 if N = 0 then
|
|
|
425 return Conv_To_Result (Conv_To_Unsigned (Min) + Random (Gen));
|
|
|
426
|
|
|
427 else
|
|
|
428 Slop := Unsigned_64'Last rem N + 1;
|
|
|
429
|
|
|
430 loop
|
|
|
431 X := Random (Gen);
|
|
|
432 exit when Slop = N or else X <= Unsigned_64'Last - Slop;
|
|
|
433 end loop;
|
|
|
434
|
|
|
435 return Conv_To_Result (Conv_To_Unsigned (Min) + X rem N);
|
|
|
436 end if;
|
|
|
437 end;
|
|
|
438
|
|
|
439 elsif Result_Subtype'Pos (Max) - Result_Subtype'Pos (Min) =
|
|
|
440 2 ** 32 - 1
|
|
|
441 then
|
|
|
442 return Result_Subtype'Val
|
|
|
443 (Result_Subtype'Pos (Min) + Unsigned_32'Pos (Random (Gen)));
|
|
|
444 else
|
|
|
445 declare
|
|
|
446 N : constant Unsigned_32 :=
|
|
|
447 Unsigned_32 (Result_Subtype'Pos (Max) -
|
|
|
448 Result_Subtype'Pos (Min) + 1);
|
|
|
449 Slop : constant Unsigned_32 := Unsigned_32'Last rem N + 1;
|
|
|
450 X : Unsigned_32;
|
|
|
451
|
|
|
452 begin
|
|
|
453 loop
|
|
|
454 X := Random (Gen);
|
|
|
455 exit when Slop = N or else X <= Unsigned_32'Last - Slop;
|
|
|
456 end loop;
|
|
|
457
|
|
|
458 return
|
|
|
459 Result_Subtype'Val
|
|
|
460 (Result_Subtype'Pos (Min) + Unsigned_32'Pos (X rem N));
|
|
|
461 end;
|
|
|
462 end if;
|
|
|
463 end Random_Discrete;
|
|
|
464
|
|
|
465 ------------------
|
|
|
466 -- Random_Float --
|
|
|
467 ------------------
|
|
|
468
|
|
|
469 function Random_Float (Gen : Generator) return Result_Subtype is
|
|
|
470 begin
|
|
|
471 if Result_Subtype'Base'Digits > Float'Digits then
|
|
|
472 return Result_Subtype'Machine (Result_Subtype
|
|
|
473 (Long_Float'(Random (Gen))));
|
|
|
474 else
|
|
|
475 return Result_Subtype'Machine (Result_Subtype
|
|
|
476 (Float'(Random (Gen))));
|
|
|
477 end if;
|
|
|
478 end Random_Float;
|
|
|
479
|
|
|
480 -----------
|
|
|
481 -- Reset --
|
|
|
482 -----------
|
|
|
483
|
|
|
484 procedure Reset (Gen : Generator) is
|
|
|
485 begin
|
|
|
486 Init (Gen, Unsigned_32'Mod (Random_Seed.Get_Seed));
|
|
|
487 end Reset;
|
|
|
488
|
|
|
489 procedure Reset (Gen : Generator; Initiator : Integer_32) is
|
|
|
490 begin
|
|
|
491 Init (Gen, To_Unsigned (Initiator));
|
|
|
492 end Reset;
|
|
|
493
|
|
|
494 procedure Reset (Gen : Generator; Initiator : Unsigned_32) is
|
|
|
495 begin
|
|
|
496 Init (Gen, Initiator);
|
|
|
497 end Reset;
|
|
|
498
|
|
|
499 procedure Reset (Gen : Generator; Initiator : Integer) is
|
|
|
500 begin
|
|
|
501 -- This is probably an unnecessary precaution against future change, but
|
|
|
502 -- since the test is a static expression, no extra code is involved.
|
|
|
503
|
|
|
504 if Integer'Size <= 32 then
|
|
|
505 Init (Gen, To_Unsigned (Integer_32 (Initiator)));
|
|
|
506
|
|
|
507 else
|
|
|
508 declare
|
|
|
509 Initiator1 : constant Unsigned_64 :=
|
|
|
510 To_Unsigned (Integer_64 (Initiator));
|
|
|
511 Init0 : constant Unsigned_32 :=
|
|
|
512 Unsigned_32 (Initiator1 mod 2 ** 32);
|
|
|
513 Init1 : constant Unsigned_32 :=
|
|
|
514 Unsigned_32 (Shift_Right (Initiator1, 32));
|
|
|
515 begin
|
|
|
516 Reset (Gen, Initialization_Vector'(Init0, Init1));
|
|
|
517 end;
|
|
|
518 end if;
|
|
|
519 end Reset;
|
|
|
520
|
|
|
521 procedure Reset (Gen : Generator; Initiator : Initialization_Vector) is
|
|
|
522 G : Generator renames Gen.Writable.Self.all;
|
|
|
523 I, J : Integer;
|
|
|
524
|
|
|
525 begin
|
|
|
526 Init (G, Seed1);
|
|
|
527 I := 1;
|
|
|
528 J := 0;
|
|
|
529
|
|
|
530 if Initiator'Length > 0 then
|
|
|
531 for K in reverse 1 .. Integer'Max (N, Initiator'Length) loop
|
|
|
532 G.S (I) :=
|
|
|
533 (G.S (I) xor ((G.S (I - 1)
|
|
|
534 xor Shift_Right (G.S (I - 1), 30)) * Mult1))
|
|
|
535 + Initiator (J + Initiator'First) + Unsigned_32 (J);
|
|
|
536
|
|
|
537 I := I + 1;
|
|
|
538 J := J + 1;
|
|
|
539
|
|
|
540 if I >= N then
|
|
|
541 G.S (0) := G.S (N - 1);
|
|
|
542 I := 1;
|
|
|
543 end if;
|
|
|
544
|
|
|
545 if J >= Initiator'Length then
|
|
|
546 J := 0;
|
|
|
547 end if;
|
|
|
548 end loop;
|
|
|
549 end if;
|
|
|
550
|
|
|
551 for K in reverse 1 .. N - 1 loop
|
|
|
552 G.S (I) :=
|
|
|
553 (G.S (I) xor ((G.S (I - 1)
|
|
|
554 xor Shift_Right (G.S (I - 1), 30)) * Mult2))
|
|
|
555 - Unsigned_32 (I);
|
|
|
556 I := I + 1;
|
|
|
557
|
|
|
558 if I >= N then
|
|
|
559 G.S (0) := G.S (N - 1);
|
|
|
560 I := 1;
|
|
|
561 end if;
|
|
|
562 end loop;
|
|
|
563
|
|
|
564 G.S (0) := Upper_Mask;
|
|
|
565 end Reset;
|
|
|
566
|
|
|
567 procedure Reset (Gen : Generator; From_State : Generator) is
|
|
|
568 G : Generator renames Gen.Writable.Self.all;
|
|
|
569 begin
|
|
|
570 G.S := From_State.S;
|
|
|
571 G.I := From_State.I;
|
|
|
572 end Reset;
|
|
|
573
|
|
|
574 procedure Reset (Gen : Generator; From_State : State) is
|
|
|
575 G : Generator renames Gen.Writable.Self.all;
|
|
|
576 begin
|
|
|
577 G.I := 0;
|
|
|
578 G.S := From_State;
|
|
|
579 end Reset;
|
|
|
580
|
|
|
581 procedure Reset (Gen : Generator; From_Image : String) is
|
|
|
582 G : Generator renames Gen.Writable.Self.all;
|
|
|
583 begin
|
|
|
584 G.I := 0;
|
|
|
585
|
|
|
586 for J in 0 .. N - 1 loop
|
|
|
587 G.S (J) := Extract_Value (From_Image, J);
|
|
|
588 end loop;
|
|
|
589 end Reset;
|
|
|
590
|
|
|
591 ----------
|
|
|
592 -- Save --
|
|
|
593 ----------
|
|
|
594
|
|
|
595 procedure Save (Gen : Generator; To_State : out State) is
|
|
|
596 Gen2 : Generator;
|
|
|
597
|
|
|
598 begin
|
|
|
599 if Gen.I = N then
|
|
|
600 Init (Gen2, 5489);
|
|
|
601 To_State := Gen2.S;
|
|
|
602
|
|
|
603 else
|
|
|
604 To_State (0 .. N - 1 - Gen.I) := Gen.S (Gen.I .. N - 1);
|
|
|
605 To_State (N - Gen.I .. N - 1) := Gen.S (0 .. Gen.I - 1);
|
|
|
606 end if;
|
|
|
607 end Save;
|
|
|
608
|
|
|
609 -----------
|
|
|
610 -- Image --
|
|
|
611 -----------
|
|
|
612
|
|
|
613 function Image (Of_State : State) return String is
|
|
|
614 Result : Image_String;
|
|
|
615
|
|
|
616 begin
|
|
|
617 Result := (others => ' ');
|
|
|
618
|
|
|
619 for J in Of_State'Range loop
|
|
|
620 Insert_Image (Result, J, Of_State (J));
|
|
|
621 end loop;
|
|
|
622
|
|
|
623 return Result;
|
|
|
624 end Image;
|
|
|
625
|
|
|
626 function Image (Gen : Generator) return String is
|
|
|
627 Result : Image_String;
|
|
|
628
|
|
|
629 begin
|
|
|
630 Result := (others => ' ');
|
|
|
631 for J in 0 .. N - 1 loop
|
|
|
632 Insert_Image (Result, J, Gen.S ((J + Gen.I) mod N));
|
|
|
633 end loop;
|
|
|
634
|
|
|
635 return Result;
|
|
|
636 end Image;
|
|
|
637
|
|
|
638 -----------
|
|
|
639 -- Value --
|
|
|
640 -----------
|
|
|
641
|
|
|
642 function Value (Coded_State : String) return State is
|
|
|
643 Gen : Generator;
|
|
|
644 S : State;
|
|
|
645 begin
|
|
|
646 Reset (Gen, Coded_State);
|
|
|
647 Save (Gen, S);
|
|
|
648 return S;
|
|
|
649 end Value;
|
|
|
650
|
|
|
651 ----------
|
|
|
652 -- Init --
|
|
|
653 ----------
|
|
|
654
|
|
|
655 procedure Init (Gen : Generator; Initiator : Unsigned_32) is
|
|
|
656 G : Generator renames Gen.Writable.Self.all;
|
|
|
657 begin
|
|
|
658 G.S (0) := Initiator;
|
|
|
659
|
|
|
660 for I in 1 .. N - 1 loop
|
|
|
661 G.S (I) :=
|
|
|
662 (G.S (I - 1) xor Shift_Right (G.S (I - 1), 30)) * Mult0
|
|
|
663 + Unsigned_32 (I);
|
|
|
664 end loop;
|
|
|
665
|
|
|
666 G.I := 0;
|
|
|
667 end Init;
|
|
|
668
|
|
|
669 ------------------
|
|
|
670 -- Insert_Image --
|
|
|
671 ------------------
|
|
|
672
|
|
|
673 procedure Insert_Image
|
|
|
674 (S : in out Image_String;
|
|
|
675 Index : Integer;
|
|
|
676 V : State_Val)
|
|
|
677 is
|
|
|
678 Value : constant String := State_Val'Image (V);
|
|
|
679 begin
|
|
|
680 S (Index * 11 + 1 .. Index * 11 + Value'Length) := Value;
|
|
|
681 end Insert_Image;
|
|
|
682
|
|
|
683 -------------------
|
|
|
684 -- Extract_Value --
|
|
|
685 -------------------
|
|
|
686
|
|
|
687 function Extract_Value (S : String; Index : Integer) return State_Val is
|
|
|
688 Start : constant Integer := S'First + Index * Image_Numeral_Length;
|
|
|
689 begin
|
|
|
690 return State_Val'Value (S (Start .. Start + Image_Numeral_Length - 1));
|
|
|
691 end Extract_Value;
|
|
|
692
|
|
|
693 end System.Random_Numbers;
|
