Skip to content

Commit f71eb55

Browse files
authored
Merge pull request RustPython#5438 from key262yek/update_random_from_v3.12.3
Update random from v3.12.7
2 parents edcc0b7 + 4f5b266 commit f71eb55

File tree

2 files changed

+236
-93
lines changed

2 files changed

+236
-93
lines changed

Lib/random.py

Lines changed: 142 additions & 50 deletions
Original file line numberDiff line numberDiff line change
@@ -32,6 +32,11 @@
3232
circular uniform
3333
von Mises
3434
35+
discrete distributions
36+
----------------------
37+
binomial
38+
39+
3540
General notes on the underlying Mersenne Twister core generator:
3641
3742
* The period is 2**19937-1.
@@ -49,6 +54,7 @@
4954
from math import log as _log, exp as _exp, pi as _pi, e as _e, ceil as _ceil
5055
from math import sqrt as _sqrt, acos as _acos, cos as _cos, sin as _sin
5156
from math import tau as TWOPI, floor as _floor, isfinite as _isfinite
57+
from math import lgamma as _lgamma, fabs as _fabs, log2 as _log2
5258
try:
5359
from os import urandom as _urandom
5460
except ImportError:
@@ -72,7 +78,7 @@ def _urandom(*args, **kwargs):
7278

7379
try:
7480
# hashlib is pretty heavy to load, try lean internal module first
75-
from _sha512 import sha512 as _sha512
81+
from _sha2 import sha512 as _sha512
7682
except ImportError:
7783
# fallback to official implementation
7884
from hashlib import sha512 as _sha512
@@ -81,6 +87,7 @@ def _urandom(*args, **kwargs):
8187
"Random",
8288
"SystemRandom",
8389
"betavariate",
90+
"binomialvariate",
8491
"choice",
8592
"choices",
8693
"expovariate",
@@ -249,7 +256,7 @@ def _randbelow_with_getrandbits(self, n):
249256
"Return a random int in the range [0,n). Defined for n > 0."
250257

251258
getrandbits = self.getrandbits
252-
k = n.bit_length() # don't use (n-1) here because n can be 1
259+
k = n.bit_length()
253260
r = getrandbits(k) # 0 <= r < 2**k
254261
while r >= n:
255262
r = getrandbits(k)
@@ -304,58 +311,25 @@ def randrange(self, start, stop=None, step=_ONE):
304311

305312
# This code is a bit messy to make it fast for the
306313
# common case while still doing adequate error checking.
307-
try:
308-
istart = _index(start)
309-
except TypeError:
310-
istart = int(start)
311-
if istart != start:
312-
_warn('randrange() will raise TypeError in the future',
313-
DeprecationWarning, 2)
314-
raise ValueError("non-integer arg 1 for randrange()")
315-
_warn('non-integer arguments to randrange() have been deprecated '
316-
'since Python 3.10 and will be removed in a subsequent '
317-
'version',
318-
DeprecationWarning, 2)
314+
istart = _index(start)
319315
if stop is None:
320316
# We don't check for "step != 1" because it hasn't been
321317
# type checked and converted to an integer yet.
322318
if step is not _ONE:
323-
raise TypeError('Missing a non-None stop argument')
319+
raise TypeError("Missing a non-None stop argument")
324320
if istart > 0:
325321
return self._randbelow(istart)
326322
raise ValueError("empty range for randrange()")
327323

328-
# stop argument supplied.
329-
try:
330-
istop = _index(stop)
331-
except TypeError:
332-
istop = int(stop)
333-
if istop != stop:
334-
_warn('randrange() will raise TypeError in the future',
335-
DeprecationWarning, 2)
336-
raise ValueError("non-integer stop for randrange()")
337-
_warn('non-integer arguments to randrange() have been deprecated '
338-
'since Python 3.10 and will be removed in a subsequent '
339-
'version',
340-
DeprecationWarning, 2)
324+
# Stop argument supplied.
325+
istop = _index(stop)
341326
width = istop - istart
342-
try:
343-
istep = _index(step)
344-
except TypeError:
345-
istep = int(step)
346-
if istep != step:
347-
_warn('randrange() will raise TypeError in the future',
348-
DeprecationWarning, 2)
349-
raise ValueError("non-integer step for randrange()")
350-
_warn('non-integer arguments to randrange() have been deprecated '
351-
'since Python 3.10 and will be removed in a subsequent '
352-
'version',
353-
DeprecationWarning, 2)
327+
istep = _index(step)
354328
# Fast path.
355329
if istep == 1:
356330
if width > 0:
357331
return istart + self._randbelow(width)
358-
raise ValueError("empty range for randrange() (%d, %d, %d)" % (istart, istop, width))
332+
raise ValueError(f"empty range in randrange({start}, {stop})")
359333

360334
# Non-unit step argument supplied.
361335
if istep > 0:
@@ -365,7 +339,7 @@ def randrange(self, start, stop=None, step=_ONE):
365339
else:
366340
raise ValueError("zero step for randrange()")
367341
if n <= 0:
368-
raise ValueError("empty range for randrange()")
342+
raise ValueError(f"empty range in randrange({start}, {stop}, {step})")
369343
return istart + istep * self._randbelow(n)
370344

371345
def randint(self, a, b):
@@ -531,7 +505,14 @@ def choices(self, population, weights=None, *, cum_weights=None, k=1):
531505
## -------------------- real-valued distributions -------------------
532506

533507
def uniform(self, a, b):
534-
"Get a random number in the range [a, b) or [a, b] depending on rounding."
508+
"""Get a random number in the range [a, b) or [a, b] depending on rounding.
509+
510+
The mean (expected value) and variance of the random variable are:
511+
512+
E[X] = (a + b) / 2
513+
Var[X] = (b - a) ** 2 / 12
514+
515+
"""
535516
return a + (b - a) * self.random()
536517

537518
def triangular(self, low=0.0, high=1.0, mode=None):
@@ -542,6 +523,11 @@ def triangular(self, low=0.0, high=1.0, mode=None):
542523
543524
http://en.wikipedia.org/wiki/Triangular_distribution
544525
526+
The mean (expected value) and variance of the random variable are:
527+
528+
E[X] = (low + high + mode) / 3
529+
Var[X] = (low**2 + high**2 + mode**2 - low*high - low*mode - high*mode) / 18
530+
545531
"""
546532
u = self.random()
547533
try:
@@ -623,7 +609,7 @@ def lognormvariate(self, mu, sigma):
623609
"""
624610
return _exp(self.normalvariate(mu, sigma))
625611

626-
def expovariate(self, lambd):
612+
def expovariate(self, lambd=1.0):
627613
"""Exponential distribution.
628614
629615
lambd is 1.0 divided by the desired mean. It should be
@@ -632,12 +618,15 @@ def expovariate(self, lambd):
632618
positive infinity if lambd is positive, and from negative
633619
infinity to 0 if lambd is negative.
634620
635-
"""
636-
# lambd: rate lambd = 1/mean
637-
# ('lambda' is a Python reserved word)
621+
The mean (expected value) and variance of the random variable are:
638622
623+
E[X] = 1 / lambd
624+
Var[X] = 1 / lambd ** 2
625+
626+
"""
639627
# we use 1-random() instead of random() to preclude the
640628
# possibility of taking the log of zero.
629+
641630
return -_log(1.0 - self.random()) / lambd
642631

643632
def vonmisesvariate(self, mu, kappa):
@@ -693,8 +682,12 @@ def gammavariate(self, alpha, beta):
693682
pdf(x) = --------------------------------------
694683
math.gamma(alpha) * beta ** alpha
695684
685+
The mean (expected value) and variance of the random variable are:
686+
687+
E[X] = alpha * beta
688+
Var[X] = alpha * beta ** 2
689+
696690
"""
697-
# alpha > 0, beta > 0, mean is alpha*beta, variance is alpha*beta**2
698691

699692
# Warning: a few older sources define the gamma distribution in terms
700693
# of alpha > -1.0
@@ -753,6 +746,11 @@ def betavariate(self, alpha, beta):
753746
Conditions on the parameters are alpha > 0 and beta > 0.
754747
Returned values range between 0 and 1.
755748
749+
The mean (expected value) and variance of the random variable are:
750+
751+
E[X] = alpha / (alpha + beta)
752+
Var[X] = alpha * beta / ((alpha + beta)**2 * (alpha + beta + 1))
753+
756754
"""
757755
## See
758756
## http://mail.python.org/pipermail/python-bugs-list/2001-January/003752.html
@@ -793,6 +791,97 @@ def weibullvariate(self, alpha, beta):
793791
return alpha * (-_log(u)) ** (1.0 / beta)
794792

795793

794+
## -------------------- discrete distributions ---------------------
795+
796+
def binomialvariate(self, n=1, p=0.5):
797+
"""Binomial random variable.
798+
799+
Gives the number of successes for *n* independent trials
800+
with the probability of success in each trial being *p*:
801+
802+
sum(random() < p for i in range(n))
803+
804+
Returns an integer in the range: 0 <= X <= n
805+
806+
The mean (expected value) and variance of the random variable are:
807+
808+
E[X] = n * p
809+
Var[x] = n * p * (1 - p)
810+
811+
"""
812+
# Error check inputs and handle edge cases
813+
if n < 0:
814+
raise ValueError("n must be non-negative")
815+
if p <= 0.0 or p >= 1.0:
816+
if p == 0.0:
817+
return 0
818+
if p == 1.0:
819+
return n
820+
raise ValueError("p must be in the range 0.0 <= p <= 1.0")
821+
822+
random = self.random
823+
824+
# Fast path for a common case
825+
if n == 1:
826+
return _index(random() < p)
827+
828+
# Exploit symmetry to establish: p <= 0.5
829+
if p > 0.5:
830+
return n - self.binomialvariate(n, 1.0 - p)
831+
832+
if n * p < 10.0:
833+
# BG: Geometric method by Devroye with running time of O(np).
834+
# https://dl.acm.org/doi/pdf/10.1145/42372.42381
835+
x = y = 0
836+
c = _log2(1.0 - p)
837+
if not c:
838+
return x
839+
while True:
840+
y += _floor(_log2(random()) / c) + 1
841+
if y > n:
842+
return x
843+
x += 1
844+
845+
# BTRS: Transformed rejection with squeeze method by Wolfgang Hörmann
846+
# https://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.47.8407&rep=rep1&type=pdf
847+
assert n*p >= 10.0 and p <= 0.5
848+
setup_complete = False
849+
850+
spq = _sqrt(n * p * (1.0 - p)) # Standard deviation of the distribution
851+
b = 1.15 + 2.53 * spq
852+
a = -0.0873 + 0.0248 * b + 0.01 * p
853+
c = n * p + 0.5
854+
vr = 0.92 - 4.2 / b
855+
856+
while True:
857+
858+
u = random()
859+
u -= 0.5
860+
us = 0.5 - _fabs(u)
861+
k = _floor((2.0 * a / us + b) * u + c)
862+
if k < 0 or k > n:
863+
continue
864+
865+
# The early-out "squeeze" test substantially reduces
866+
# the number of acceptance condition evaluations.
867+
v = random()
868+
if us >= 0.07 and v <= vr:
869+
return k
870+
871+
# Acceptance-rejection test.
872+
# Note, the original paper erroneously omits the call to log(v)
873+
# when comparing to the log of the rescaled binomial distribution.
874+
if not setup_complete:
875+
alpha = (2.83 + 5.1 / b) * spq
876+
lpq = _log(p / (1.0 - p))
877+
m = _floor((n + 1) * p) # Mode of the distribution
878+
h = _lgamma(m + 1) + _lgamma(n - m + 1)
879+
setup_complete = True # Only needs to be done once
880+
v *= alpha / (a / (us * us) + b)
881+
if _log(v) <= h - _lgamma(k + 1) - _lgamma(n - k + 1) + (k - m) * lpq:
882+
return k
883+
884+
796885
## ------------------------------------------------------------------
797886
## --------------- Operating System Random Source ------------------
798887

@@ -859,6 +948,7 @@ def _notimplemented(self, *args, **kwds):
859948
gammavariate = _inst.gammavariate
860949
gauss = _inst.gauss
861950
betavariate = _inst.betavariate
951+
binomialvariate = _inst.binomialvariate
862952
paretovariate = _inst.paretovariate
863953
weibullvariate = _inst.weibullvariate
864954
getstate = _inst.getstate
@@ -883,15 +973,17 @@ def _test_generator(n, func, args):
883973
low = min(data)
884974
high = max(data)
885975

886-
print(f'{t1 - t0:.3f} sec, {n} times {func.__name__}')
976+
print(f'{t1 - t0:.3f} sec, {n} times {func.__name__}{args!r}')
887977
print('avg %g, stddev %g, min %g, max %g\n' % (xbar, sigma, low, high))
888978

889979

890-
def _test(N=2000):
980+
def _test(N=10_000):
891981
_test_generator(N, random, ())
892982
_test_generator(N, normalvariate, (0.0, 1.0))
893983
_test_generator(N, lognormvariate, (0.0, 1.0))
894984
_test_generator(N, vonmisesvariate, (0.0, 1.0))
985+
_test_generator(N, binomialvariate, (15, 0.60))
986+
_test_generator(N, binomialvariate, (100, 0.75))
895987
_test_generator(N, gammavariate, (0.01, 1.0))
896988
_test_generator(N, gammavariate, (0.1, 1.0))
897989
_test_generator(N, gammavariate, (0.1, 2.0))

0 commit comments

Comments
 (0)