blob: a6454f520df0a321495b5da6f8dade6e4c44036f [file] [log] [blame]
Olivier Deprezf4ef2d02021-04-20 13:36:24 +02001"""Random variable generators.
2
3 bytes
4 -----
5 uniform bytes (values between 0 and 255)
6
7 integers
8 --------
9 uniform within range
10
11 sequences
12 ---------
13 pick random element
14 pick random sample
15 pick weighted random sample
16 generate random permutation
17
18 distributions on the real line:
19 ------------------------------
20 uniform
21 triangular
22 normal (Gaussian)
23 lognormal
24 negative exponential
25 gamma
26 beta
27 pareto
28 Weibull
29
30 distributions on the circle (angles 0 to 2pi)
31 ---------------------------------------------
32 circular uniform
33 von Mises
34
35General notes on the underlying Mersenne Twister core generator:
36
37* The period is 2**19937-1.
38* It is one of the most extensively tested generators in existence.
39* The random() method is implemented in C, executes in a single Python step,
40 and is, therefore, threadsafe.
41
42"""
43
44# Translated by Guido van Rossum from C source provided by
45# Adrian Baddeley. Adapted by Raymond Hettinger for use with
46# the Mersenne Twister and os.urandom() core generators.
47
48from warnings import warn as _warn
49from math import log as _log, exp as _exp, pi as _pi, e as _e, ceil as _ceil
50from math import sqrt as _sqrt, acos as _acos, cos as _cos, sin as _sin
51from math import tau as TWOPI, floor as _floor
52from os import urandom as _urandom
53from _collections_abc import Set as _Set, Sequence as _Sequence
54from itertools import accumulate as _accumulate, repeat as _repeat
55from bisect import bisect as _bisect
56import os as _os
57import _random
58
59try:
60 # hashlib is pretty heavy to load, try lean internal module first
61 from _sha512 import sha512 as _sha512
62except ImportError:
63 # fallback to official implementation
64 from hashlib import sha512 as _sha512
65
66__all__ = [
67 "Random",
68 "SystemRandom",
69 "betavariate",
70 "choice",
71 "choices",
72 "expovariate",
73 "gammavariate",
74 "gauss",
75 "getrandbits",
76 "getstate",
77 "lognormvariate",
78 "normalvariate",
79 "paretovariate",
80 "randint",
81 "random",
82 "randrange",
83 "sample",
84 "seed",
85 "setstate",
86 "shuffle",
87 "triangular",
88 "uniform",
89 "vonmisesvariate",
90 "weibullvariate",
91]
92
93NV_MAGICCONST = 4 * _exp(-0.5) / _sqrt(2.0)
94LOG4 = _log(4.0)
95SG_MAGICCONST = 1.0 + _log(4.5)
96BPF = 53 # Number of bits in a float
97RECIP_BPF = 2 ** -BPF
98
99
100class Random(_random.Random):
101 """Random number generator base class used by bound module functions.
102
103 Used to instantiate instances of Random to get generators that don't
104 share state.
105
106 Class Random can also be subclassed if you want to use a different basic
107 generator of your own devising: in that case, override the following
108 methods: random(), seed(), getstate(), and setstate().
109 Optionally, implement a getrandbits() method so that randrange()
110 can cover arbitrarily large ranges.
111
112 """
113
114 VERSION = 3 # used by getstate/setstate
115
116 def __init__(self, x=None):
117 """Initialize an instance.
118
119 Optional argument x controls seeding, as for Random.seed().
120 """
121
122 self.seed(x)
123 self.gauss_next = None
124
125 def seed(self, a=None, version=2):
126 """Initialize internal state from a seed.
127
128 The only supported seed types are None, int, float,
129 str, bytes, and bytearray.
130
131 None or no argument seeds from current time or from an operating
132 system specific randomness source if available.
133
134 If *a* is an int, all bits are used.
135
136 For version 2 (the default), all of the bits are used if *a* is a str,
137 bytes, or bytearray. For version 1 (provided for reproducing random
138 sequences from older versions of Python), the algorithm for str and
139 bytes generates a narrower range of seeds.
140
141 """
142
143 if version == 1 and isinstance(a, (str, bytes)):
144 a = a.decode('latin-1') if isinstance(a, bytes) else a
145 x = ord(a[0]) << 7 if a else 0
146 for c in map(ord, a):
147 x = ((1000003 * x) ^ c) & 0xFFFFFFFFFFFFFFFF
148 x ^= len(a)
149 a = -2 if x == -1 else x
150
151 elif version == 2 and isinstance(a, (str, bytes, bytearray)):
152 if isinstance(a, str):
153 a = a.encode()
154 a += _sha512(a).digest()
155 a = int.from_bytes(a, 'big')
156
157 elif not isinstance(a, (type(None), int, float, str, bytes, bytearray)):
158 _warn('Seeding based on hashing is deprecated\n'
159 'since Python 3.9 and will be removed in a subsequent '
160 'version. The only \n'
161 'supported seed types are: None, '
162 'int, float, str, bytes, and bytearray.',
163 DeprecationWarning, 2)
164
165 super().seed(a)
166 self.gauss_next = None
167
168 def getstate(self):
169 """Return internal state; can be passed to setstate() later."""
170 return self.VERSION, super().getstate(), self.gauss_next
171
172 def setstate(self, state):
173 """Restore internal state from object returned by getstate()."""
174 version = state[0]
175 if version == 3:
176 version, internalstate, self.gauss_next = state
177 super().setstate(internalstate)
178 elif version == 2:
179 version, internalstate, self.gauss_next = state
180 # In version 2, the state was saved as signed ints, which causes
181 # inconsistencies between 32/64-bit systems. The state is
182 # really unsigned 32-bit ints, so we convert negative ints from
183 # version 2 to positive longs for version 3.
184 try:
185 internalstate = tuple(x % (2 ** 32) for x in internalstate)
186 except ValueError as e:
187 raise TypeError from e
188 super().setstate(internalstate)
189 else:
190 raise ValueError("state with version %s passed to "
191 "Random.setstate() of version %s" %
192 (version, self.VERSION))
193
194
195 ## -------------------------------------------------------
196 ## ---- Methods below this point do not need to be overridden or extended
197 ## ---- when subclassing for the purpose of using a different core generator.
198
199
200 ## -------------------- pickle support -------------------
201
202 # Issue 17489: Since __reduce__ was defined to fix #759889 this is no
203 # longer called; we leave it here because it has been here since random was
204 # rewritten back in 2001 and why risk breaking something.
205 def __getstate__(self): # for pickle
206 return self.getstate()
207
208 def __setstate__(self, state): # for pickle
209 self.setstate(state)
210
211 def __reduce__(self):
212 return self.__class__, (), self.getstate()
213
214
215 ## ---- internal support method for evenly distributed integers ----
216
217 def __init_subclass__(cls, /, **kwargs):
218 """Control how subclasses generate random integers.
219
220 The algorithm a subclass can use depends on the random() and/or
221 getrandbits() implementation available to it and determines
222 whether it can generate random integers from arbitrarily large
223 ranges.
224 """
225
226 for c in cls.__mro__:
227 if '_randbelow' in c.__dict__:
228 # just inherit it
229 break
230 if 'getrandbits' in c.__dict__:
231 cls._randbelow = cls._randbelow_with_getrandbits
232 break
233 if 'random' in c.__dict__:
234 cls._randbelow = cls._randbelow_without_getrandbits
235 break
236
237 def _randbelow_with_getrandbits(self, n):
238 "Return a random int in the range [0,n). Returns 0 if n==0."
239
240 if not n:
241 return 0
242 getrandbits = self.getrandbits
243 k = n.bit_length() # don't use (n-1) here because n can be 1
244 r = getrandbits(k) # 0 <= r < 2**k
245 while r >= n:
246 r = getrandbits(k)
247 return r
248
249 def _randbelow_without_getrandbits(self, n, maxsize=1<<BPF):
250 """Return a random int in the range [0,n). Returns 0 if n==0.
251
252 The implementation does not use getrandbits, but only random.
253 """
254
255 random = self.random
256 if n >= maxsize:
257 _warn("Underlying random() generator does not supply \n"
258 "enough bits to choose from a population range this large.\n"
259 "To remove the range limitation, add a getrandbits() method.")
260 return _floor(random() * n)
261 if n == 0:
262 return 0
263 rem = maxsize % n
264 limit = (maxsize - rem) / maxsize # int(limit * maxsize) % n == 0
265 r = random()
266 while r >= limit:
267 r = random()
268 return _floor(r * maxsize) % n
269
270 _randbelow = _randbelow_with_getrandbits
271
272
273 ## --------------------------------------------------------
274 ## ---- Methods below this point generate custom distributions
275 ## ---- based on the methods defined above. They do not
276 ## ---- directly touch the underlying generator and only
277 ## ---- access randomness through the methods: random(),
278 ## ---- getrandbits(), or _randbelow().
279
280
281 ## -------------------- bytes methods ---------------------
282
283 def randbytes(self, n):
284 """Generate n random bytes."""
285 return self.getrandbits(n * 8).to_bytes(n, 'little')
286
287
288 ## -------------------- integer methods -------------------
289
290 def randrange(self, start, stop=None, step=1):
291 """Choose a random item from range(start, stop[, step]).
292
293 This fixes the problem with randint() which includes the
294 endpoint; in Python this is usually not what you want.
295
296 """
297
298 # This code is a bit messy to make it fast for the
299 # common case while still doing adequate error checking.
300 istart = int(start)
301 if istart != start:
302 raise ValueError("non-integer arg 1 for randrange()")
303 if stop is None:
304 if istart > 0:
305 return self._randbelow(istart)
306 raise ValueError("empty range for randrange()")
307
308 # stop argument supplied.
309 istop = int(stop)
310 if istop != stop:
311 raise ValueError("non-integer stop for randrange()")
312 width = istop - istart
313 if step == 1 and width > 0:
314 return istart + self._randbelow(width)
315 if step == 1:
316 raise ValueError("empty range for randrange() (%d, %d, %d)" % (istart, istop, width))
317
318 # Non-unit step argument supplied.
319 istep = int(step)
320 if istep != step:
321 raise ValueError("non-integer step for randrange()")
322 if istep > 0:
323 n = (width + istep - 1) // istep
324 elif istep < 0:
325 n = (width + istep + 1) // istep
326 else:
327 raise ValueError("zero step for randrange()")
328
329 if n <= 0:
330 raise ValueError("empty range for randrange()")
331
332 return istart + istep * self._randbelow(n)
333
334 def randint(self, a, b):
335 """Return random integer in range [a, b], including both end points.
336 """
337
338 return self.randrange(a, b+1)
339
340
341 ## -------------------- sequence methods -------------------
342
343 def choice(self, seq):
344 """Choose a random element from a non-empty sequence."""
345 # raises IndexError if seq is empty
346 return seq[self._randbelow(len(seq))]
347
348 def shuffle(self, x, random=None):
349 """Shuffle list x in place, and return None.
350
351 Optional argument random is a 0-argument function returning a
352 random float in [0.0, 1.0); if it is the default None, the
353 standard random.random will be used.
354
355 """
356
357 if random is None:
358 randbelow = self._randbelow
359 for i in reversed(range(1, len(x))):
360 # pick an element in x[:i+1] with which to exchange x[i]
361 j = randbelow(i + 1)
362 x[i], x[j] = x[j], x[i]
363 else:
364 _warn('The *random* parameter to shuffle() has been deprecated\n'
365 'since Python 3.9 and will be removed in a subsequent '
366 'version.',
367 DeprecationWarning, 2)
368 floor = _floor
369 for i in reversed(range(1, len(x))):
370 # pick an element in x[:i+1] with which to exchange x[i]
371 j = floor(random() * (i + 1))
372 x[i], x[j] = x[j], x[i]
373
374 def sample(self, population, k, *, counts=None):
375 """Chooses k unique random elements from a population sequence or set.
376
377 Returns a new list containing elements from the population while
378 leaving the original population unchanged. The resulting list is
379 in selection order so that all sub-slices will also be valid random
380 samples. This allows raffle winners (the sample) to be partitioned
381 into grand prize and second place winners (the subslices).
382
383 Members of the population need not be hashable or unique. If the
384 population contains repeats, then each occurrence is a possible
385 selection in the sample.
386
387 Repeated elements can be specified one at a time or with the optional
388 counts parameter. For example:
389
390 sample(['red', 'blue'], counts=[4, 2], k=5)
391
392 is equivalent to:
393
394 sample(['red', 'red', 'red', 'red', 'blue', 'blue'], k=5)
395
396 To choose a sample from a range of integers, use range() for the
397 population argument. This is especially fast and space efficient
398 for sampling from a large population:
399
400 sample(range(10000000), 60)
401
402 """
403
404 # Sampling without replacement entails tracking either potential
405 # selections (the pool) in a list or previous selections in a set.
406
407 # When the number of selections is small compared to the
408 # population, then tracking selections is efficient, requiring
409 # only a small set and an occasional reselection. For
410 # a larger number of selections, the pool tracking method is
411 # preferred since the list takes less space than the
412 # set and it doesn't suffer from frequent reselections.
413
414 # The number of calls to _randbelow() is kept at or near k, the
415 # theoretical minimum. This is important because running time
416 # is dominated by _randbelow() and because it extracts the
417 # least entropy from the underlying random number generators.
418
419 # Memory requirements are kept to the smaller of a k-length
420 # set or an n-length list.
421
422 # There are other sampling algorithms that do not require
423 # auxiliary memory, but they were rejected because they made
424 # too many calls to _randbelow(), making them slower and
425 # causing them to eat more entropy than necessary.
426
427 if isinstance(population, _Set):
428 _warn('Sampling from a set deprecated\n'
429 'since Python 3.9 and will be removed in a subsequent version.',
430 DeprecationWarning, 2)
431 population = tuple(population)
432 if not isinstance(population, _Sequence):
433 raise TypeError("Population must be a sequence. For dicts or sets, use sorted(d).")
434 n = len(population)
435 if counts is not None:
436 cum_counts = list(_accumulate(counts))
437 if len(cum_counts) != n:
438 raise ValueError('The number of counts does not match the population')
439 total = cum_counts.pop()
440 if not isinstance(total, int):
441 raise TypeError('Counts must be integers')
442 if total <= 0:
443 raise ValueError('Total of counts must be greater than zero')
444 selections = sample(range(total), k=k)
445 bisect = _bisect
446 return [population[bisect(cum_counts, s)] for s in selections]
447 randbelow = self._randbelow
448 if not 0 <= k <= n:
449 raise ValueError("Sample larger than population or is negative")
450 result = [None] * k
451 setsize = 21 # size of a small set minus size of an empty list
452 if k > 5:
453 setsize += 4 ** _ceil(_log(k * 3, 4)) # table size for big sets
454 if n <= setsize:
455 # An n-length list is smaller than a k-length set.
456 # Invariant: non-selected at pool[0 : n-i]
457 pool = list(population)
458 for i in range(k):
459 j = randbelow(n - i)
460 result[i] = pool[j]
461 pool[j] = pool[n - i - 1] # move non-selected item into vacancy
462 else:
463 selected = set()
464 selected_add = selected.add
465 for i in range(k):
466 j = randbelow(n)
467 while j in selected:
468 j = randbelow(n)
469 selected_add(j)
470 result[i] = population[j]
471 return result
472
473 def choices(self, population, weights=None, *, cum_weights=None, k=1):
474 """Return a k sized list of population elements chosen with replacement.
475
476 If the relative weights or cumulative weights are not specified,
477 the selections are made with equal probability.
478
479 """
480 random = self.random
481 n = len(population)
482 if cum_weights is None:
483 if weights is None:
484 floor = _floor
485 n += 0.0 # convert to float for a small speed improvement
486 return [population[floor(random() * n)] for i in _repeat(None, k)]
487 cum_weights = list(_accumulate(weights))
488 elif weights is not None:
489 raise TypeError('Cannot specify both weights and cumulative weights')
490 if len(cum_weights) != n:
491 raise ValueError('The number of weights does not match the population')
492 total = cum_weights[-1] + 0.0 # convert to float
493 if total <= 0.0:
494 raise ValueError('Total of weights must be greater than zero')
495 bisect = _bisect
496 hi = n - 1
497 return [population[bisect(cum_weights, random() * total, 0, hi)]
498 for i in _repeat(None, k)]
499
500
501 ## -------------------- real-valued distributions -------------------
502
503 def uniform(self, a, b):
504 "Get a random number in the range [a, b) or [a, b] depending on rounding."
505 return a + (b - a) * self.random()
506
507 def triangular(self, low=0.0, high=1.0, mode=None):
508 """Triangular distribution.
509
510 Continuous distribution bounded by given lower and upper limits,
511 and having a given mode value in-between.
512
513 http://en.wikipedia.org/wiki/Triangular_distribution
514
515 """
516 u = self.random()
517 try:
518 c = 0.5 if mode is None else (mode - low) / (high - low)
519 except ZeroDivisionError:
520 return low
521 if u > c:
522 u = 1.0 - u
523 c = 1.0 - c
524 low, high = high, low
525 return low + (high - low) * _sqrt(u * c)
526
527 def normalvariate(self, mu, sigma):
528 """Normal distribution.
529
530 mu is the mean, and sigma is the standard deviation.
531
532 """
533 # Uses Kinderman and Monahan method. Reference: Kinderman,
534 # A.J. and Monahan, J.F., "Computer generation of random
535 # variables using the ratio of uniform deviates", ACM Trans
536 # Math Software, 3, (1977), pp257-260.
537
538 random = self.random
539 while True:
540 u1 = random()
541 u2 = 1.0 - random()
542 z = NV_MAGICCONST * (u1 - 0.5) / u2
543 zz = z * z / 4.0
544 if zz <= -_log(u2):
545 break
546 return mu + z * sigma
547
548 def gauss(self, mu, sigma):
549 """Gaussian distribution.
550
551 mu is the mean, and sigma is the standard deviation. This is
552 slightly faster than the normalvariate() function.
553
554 Not thread-safe without a lock around calls.
555
556 """
557 # When x and y are two variables from [0, 1), uniformly
558 # distributed, then
559 #
560 # cos(2*pi*x)*sqrt(-2*log(1-y))
561 # sin(2*pi*x)*sqrt(-2*log(1-y))
562 #
563 # are two *independent* variables with normal distribution
564 # (mu = 0, sigma = 1).
565 # (Lambert Meertens)
566 # (corrected version; bug discovered by Mike Miller, fixed by LM)
567
568 # Multithreading note: When two threads call this function
569 # simultaneously, it is possible that they will receive the
570 # same return value. The window is very small though. To
571 # avoid this, you have to use a lock around all calls. (I
572 # didn't want to slow this down in the serial case by using a
573 # lock here.)
574
575 random = self.random
576 z = self.gauss_next
577 self.gauss_next = None
578 if z is None:
579 x2pi = random() * TWOPI
580 g2rad = _sqrt(-2.0 * _log(1.0 - random()))
581 z = _cos(x2pi) * g2rad
582 self.gauss_next = _sin(x2pi) * g2rad
583
584 return mu + z * sigma
585
586 def lognormvariate(self, mu, sigma):
587 """Log normal distribution.
588
589 If you take the natural logarithm of this distribution, you'll get a
590 normal distribution with mean mu and standard deviation sigma.
591 mu can have any value, and sigma must be greater than zero.
592
593 """
594 return _exp(self.normalvariate(mu, sigma))
595
596 def expovariate(self, lambd):
597 """Exponential distribution.
598
599 lambd is 1.0 divided by the desired mean. It should be
600 nonzero. (The parameter would be called "lambda", but that is
601 a reserved word in Python.) Returned values range from 0 to
602 positive infinity if lambd is positive, and from negative
603 infinity to 0 if lambd is negative.
604
605 """
606 # lambd: rate lambd = 1/mean
607 # ('lambda' is a Python reserved word)
608
609 # we use 1-random() instead of random() to preclude the
610 # possibility of taking the log of zero.
611 return -_log(1.0 - self.random()) / lambd
612
613 def vonmisesvariate(self, mu, kappa):
614 """Circular data distribution.
615
616 mu is the mean angle, expressed in radians between 0 and 2*pi, and
617 kappa is the concentration parameter, which must be greater than or
618 equal to zero. If kappa is equal to zero, this distribution reduces
619 to a uniform random angle over the range 0 to 2*pi.
620
621 """
622 # Based upon an algorithm published in: Fisher, N.I.,
623 # "Statistical Analysis of Circular Data", Cambridge
624 # University Press, 1993.
625
626 # Thanks to Magnus Kessler for a correction to the
627 # implementation of step 4.
628
629 random = self.random
630 if kappa <= 1e-6:
631 return TWOPI * random()
632
633 s = 0.5 / kappa
634 r = s + _sqrt(1.0 + s * s)
635
636 while True:
637 u1 = random()
638 z = _cos(_pi * u1)
639
640 d = z / (r + z)
641 u2 = random()
642 if u2 < 1.0 - d * d or u2 <= (1.0 - d) * _exp(d):
643 break
644
645 q = 1.0 / r
646 f = (q + z) / (1.0 + q * z)
647 u3 = random()
648 if u3 > 0.5:
649 theta = (mu + _acos(f)) % TWOPI
650 else:
651 theta = (mu - _acos(f)) % TWOPI
652
653 return theta
654
655 def gammavariate(self, alpha, beta):
656 """Gamma distribution. Not the gamma function!
657
658 Conditions on the parameters are alpha > 0 and beta > 0.
659
660 The probability distribution function is:
661
662 x ** (alpha - 1) * math.exp(-x / beta)
663 pdf(x) = --------------------------------------
664 math.gamma(alpha) * beta ** alpha
665
666 """
667 # alpha > 0, beta > 0, mean is alpha*beta, variance is alpha*beta**2
668
669 # Warning: a few older sources define the gamma distribution in terms
670 # of alpha > -1.0
671 if alpha <= 0.0 or beta <= 0.0:
672 raise ValueError('gammavariate: alpha and beta must be > 0.0')
673
674 random = self.random
675 if alpha > 1.0:
676
677 # Uses R.C.H. Cheng, "The generation of Gamma
678 # variables with non-integral shape parameters",
679 # Applied Statistics, (1977), 26, No. 1, p71-74
680
681 ainv = _sqrt(2.0 * alpha - 1.0)
682 bbb = alpha - LOG4
683 ccc = alpha + ainv
684
685 while 1:
686 u1 = random()
687 if not 1e-7 < u1 < 0.9999999:
688 continue
689 u2 = 1.0 - random()
690 v = _log(u1 / (1.0 - u1)) / ainv
691 x = alpha * _exp(v)
692 z = u1 * u1 * u2
693 r = bbb + ccc * v - x
694 if r + SG_MAGICCONST - 4.5 * z >= 0.0 or r >= _log(z):
695 return x * beta
696
697 elif alpha == 1.0:
698 # expovariate(1/beta)
699 return -_log(1.0 - random()) * beta
700
701 else:
702 # alpha is between 0 and 1 (exclusive)
703 # Uses ALGORITHM GS of Statistical Computing - Kennedy & Gentle
704 while True:
705 u = random()
706 b = (_e + alpha) / _e
707 p = b * u
708 if p <= 1.0:
709 x = p ** (1.0 / alpha)
710 else:
711 x = -_log((b - p) / alpha)
712 u1 = random()
713 if p > 1.0:
714 if u1 <= x ** (alpha - 1.0):
715 break
716 elif u1 <= _exp(-x):
717 break
718 return x * beta
719
720 def betavariate(self, alpha, beta):
721 """Beta distribution.
722
723 Conditions on the parameters are alpha > 0 and beta > 0.
724 Returned values range between 0 and 1.
725
726 """
727 ## See
728 ## http://mail.python.org/pipermail/python-bugs-list/2001-January/003752.html
729 ## for Ivan Frohne's insightful analysis of why the original implementation:
730 ##
731 ## def betavariate(self, alpha, beta):
732 ## # Discrete Event Simulation in C, pp 87-88.
733 ##
734 ## y = self.expovariate(alpha)
735 ## z = self.expovariate(1.0/beta)
736 ## return z/(y+z)
737 ##
738 ## was dead wrong, and how it probably got that way.
739
740 # This version due to Janne Sinkkonen, and matches all the std
741 # texts (e.g., Knuth Vol 2 Ed 3 pg 134 "the beta distribution").
742 y = self.gammavariate(alpha, 1.0)
743 if y:
744 return y / (y + self.gammavariate(beta, 1.0))
745 return 0.0
746
747 def paretovariate(self, alpha):
748 """Pareto distribution. alpha is the shape parameter."""
749 # Jain, pg. 495
750
751 u = 1.0 - self.random()
752 return 1.0 / u ** (1.0 / alpha)
753
754 def weibullvariate(self, alpha, beta):
755 """Weibull distribution.
756
757 alpha is the scale parameter and beta is the shape parameter.
758
759 """
760 # Jain, pg. 499; bug fix courtesy Bill Arms
761
762 u = 1.0 - self.random()
763 return alpha * (-_log(u)) ** (1.0 / beta)
764
765
766## ------------------------------------------------------------------
767## --------------- Operating System Random Source ------------------
768
769
770class SystemRandom(Random):
771 """Alternate random number generator using sources provided
772 by the operating system (such as /dev/urandom on Unix or
773 CryptGenRandom on Windows).
774
775 Not available on all systems (see os.urandom() for details).
776
777 """
778
779 def random(self):
780 """Get the next random number in the range [0.0, 1.0)."""
781 return (int.from_bytes(_urandom(7), 'big') >> 3) * RECIP_BPF
782
783 def getrandbits(self, k):
784 """getrandbits(k) -> x. Generates an int with k random bits."""
785 if k < 0:
786 raise ValueError('number of bits must be non-negative')
787 numbytes = (k + 7) // 8 # bits / 8 and rounded up
788 x = int.from_bytes(_urandom(numbytes), 'big')
789 return x >> (numbytes * 8 - k) # trim excess bits
790
791 def randbytes(self, n):
792 """Generate n random bytes."""
793 # os.urandom(n) fails with ValueError for n < 0
794 # and returns an empty bytes string for n == 0.
795 return _urandom(n)
796
797 def seed(self, *args, **kwds):
798 "Stub method. Not used for a system random number generator."
799 return None
800
801 def _notimplemented(self, *args, **kwds):
802 "Method should not be called for a system random number generator."
803 raise NotImplementedError('System entropy source does not have state.')
804 getstate = setstate = _notimplemented
805
806
807# ----------------------------------------------------------------------
808# Create one instance, seeded from current time, and export its methods
809# as module-level functions. The functions share state across all uses
810# (both in the user's code and in the Python libraries), but that's fine
811# for most programs and is easier for the casual user than making them
812# instantiate their own Random() instance.
813
814_inst = Random()
815seed = _inst.seed
816random = _inst.random
817uniform = _inst.uniform
818triangular = _inst.triangular
819randint = _inst.randint
820choice = _inst.choice
821randrange = _inst.randrange
822sample = _inst.sample
823shuffle = _inst.shuffle
824choices = _inst.choices
825normalvariate = _inst.normalvariate
826lognormvariate = _inst.lognormvariate
827expovariate = _inst.expovariate
828vonmisesvariate = _inst.vonmisesvariate
829gammavariate = _inst.gammavariate
830gauss = _inst.gauss
831betavariate = _inst.betavariate
832paretovariate = _inst.paretovariate
833weibullvariate = _inst.weibullvariate
834getstate = _inst.getstate
835setstate = _inst.setstate
836getrandbits = _inst.getrandbits
837randbytes = _inst.randbytes
838
839
840## ------------------------------------------------------
841## ----------------- test program -----------------------
842
843def _test_generator(n, func, args):
844 from statistics import stdev, fmean as mean
845 from time import perf_counter
846
847 t0 = perf_counter()
848 data = [func(*args) for i in range(n)]
849 t1 = perf_counter()
850
851 xbar = mean(data)
852 sigma = stdev(data, xbar)
853 low = min(data)
854 high = max(data)
855
856 print(f'{t1 - t0:.3f} sec, {n} times {func.__name__}')
857 print('avg %g, stddev %g, min %g, max %g\n' % (xbar, sigma, low, high))
858
859
860def _test(N=2000):
861 _test_generator(N, random, ())
862 _test_generator(N, normalvariate, (0.0, 1.0))
863 _test_generator(N, lognormvariate, (0.0, 1.0))
864 _test_generator(N, vonmisesvariate, (0.0, 1.0))
865 _test_generator(N, gammavariate, (0.01, 1.0))
866 _test_generator(N, gammavariate, (0.1, 1.0))
867 _test_generator(N, gammavariate, (0.1, 2.0))
868 _test_generator(N, gammavariate, (0.5, 1.0))
869 _test_generator(N, gammavariate, (0.9, 1.0))
870 _test_generator(N, gammavariate, (1.0, 1.0))
871 _test_generator(N, gammavariate, (2.0, 1.0))
872 _test_generator(N, gammavariate, (20.0, 1.0))
873 _test_generator(N, gammavariate, (200.0, 1.0))
874 _test_generator(N, gauss, (0.0, 1.0))
875 _test_generator(N, betavariate, (3.0, 3.0))
876 _test_generator(N, triangular, (0.0, 1.0, 1.0 / 3.0))
877
878
879## ------------------------------------------------------
880## ------------------ fork support ---------------------
881
882if hasattr(_os, "fork"):
883 _os.register_at_fork(after_in_child=_inst.seed)
884
885
886if __name__ == '__main__':
887 _test()