Correct long standing bugs in the methods for random distributions.
The range of u=random() is [0,1), so log(u) and 1/x can fail. Fix by setting u=1-random() or by reselecting for a usable value. Will backport.
This commit is contained in:
parent
3a57d9de07
commit
73ced7ee99
|
@ -282,7 +282,7 @@ class Random(CoreGenerator):
|
|||
random = self.random
|
||||
while True:
|
||||
u1 = random()
|
||||
u2 = random()
|
||||
u2 = 1.0 - random()
|
||||
z = NV_MAGICCONST*(u1-0.5)/u2
|
||||
zz = z*z/4.0
|
||||
if zz <= -_log(u2):
|
||||
|
@ -422,7 +422,9 @@ class Random(CoreGenerator):
|
|||
|
||||
while True:
|
||||
u1 = random()
|
||||
u2 = random()
|
||||
if not 1e-7 < u1 < .9999999:
|
||||
continue
|
||||
u2 = 1.0 - random()
|
||||
v = _log(u1/(1.0-u1))/ainv
|
||||
x = alpha*_exp(v)
|
||||
z = u1*u1*u2
|
||||
|
@ -554,7 +556,7 @@ class Random(CoreGenerator):
|
|||
"""Pareto distribution. alpha is the shape parameter."""
|
||||
# Jain, pg. 495
|
||||
|
||||
u = self.random()
|
||||
u = 1.0 - self.random()
|
||||
return 1.0 / pow(u, 1.0/alpha)
|
||||
|
||||
## -------------------- Weibull --------------------
|
||||
|
@ -567,7 +569,7 @@ class Random(CoreGenerator):
|
|||
"""
|
||||
# Jain, pg. 499; bug fix courtesy Bill Arms
|
||||
|
||||
u = self.random()
|
||||
u = 1.0 - self.random()
|
||||
return alpha * pow(-_log(u), 1.0/beta)
|
||||
|
||||
## -------------------- Wichmann-Hill -------------------
|
||||
|
|
Loading…
Reference in New Issue