ยปCore Development>Code coverage>Lib/random.py

Python code coverage for Lib/random.py

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