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