home *** CD-ROM | disk | FTP | other *** search
/ MacHack 2000 / MacHack 2000.toast / pc / The Hacks / MacHacksBug / Python 1.5.2c1 / Lib / random.pyc (.txt) < prev    next >
Encoding:
Python Compiled Bytecode  |  2000-06-23  |  9.5 KB  |  249 lines
  1. # Source Generated with Decompyle++
  2. # File: in.pyc (Python 1.5)
  3.  
  4. import whrandom
  5. randrange
  6. sin
  7.  
  8. def makeseed(a = None):
  9.     '''Turn a hashable value into three seed values for whrandom.seed().
  10.  
  11. \tNone or no argument returns (0, 0, 0), to seed from current time.
  12.  
  13. \t'''
  14.     if a is None:
  15.         return (0, 0, 0)
  16.     
  17.     a = hash(a)
  18.     (a, x) = divmod(a, 256)
  19.     (a, y) = divmod(a, 256)
  20.     (a, z) = divmod(a, 256)
  21.     if not (x + a) % 256:
  22.         pass
  23.     x = 1
  24.     if not (y + a) % 256:
  25.         pass
  26.     y = 1
  27.     if not (z + a) % 256:
  28.         pass
  29.     z = 1
  30.     return (x, y, z)
  31.  
  32.  
  33. def seed(a = None):
  34.     '''Seed the default generator from any hashable value.
  35.  
  36. \tNone or no argument returns (0, 0, 0) to seed from current time.
  37.  
  38. \t'''
  39.     (x, y, z) = makeseed(a)
  40.     whrandom.seed(x, y, z)
  41.  
  42.  
  43. class generator(whrandom.whrandom):
  44.     '''Random generator class.'''
  45.     
  46.     def __init__(self, a = None):
  47.         '''Constructor.  Seed from current time or hashable value.'''
  48.         self.seed(a)
  49.  
  50.     
  51.     def seed(self, a = None):
  52.         '''Seed the generator from current time or hashable value.'''
  53.         (x, y, z) = makeseed(a)
  54.         whrandom.whrandom.seed(self, x, y, z)
  55.  
  56.  
  57.  
  58. def new_generator(a = None):
  59.     '''Return a new random generator instance.'''
  60.     return generator(a)
  61.  
  62.  
  63. def verify(name, expected):
  64.     computed = eval(name)
  65.     if abs(computed - expected) > 1e-07:
  66.         raise ValueError, 'computed value for %s deviates too much (computed %g, expected %g)' % (name, computed, expected)
  67.     
  68.  
  69. NV_MAGICCONST = 4 * exp(-0.5) / sqrt(2.0)
  70. verify('NV_MAGICCONST', 1.71552776992)
  71.  
  72. def normalvariate(mu, sigma):
  73.     while 1:
  74.         u1 = random()
  75.         u2 = random()
  76.         z = NV_MAGICCONST * (u1 - 0.5) / u2
  77.         zz = z * z / 4.0
  78.         if zz <= -log(u2):
  79.             break
  80.         
  81.     return mu + z * sigma
  82.  
  83.  
  84. def lognormvariate(mu, sigma):
  85.     return exp(normalvariate(mu, sigma))
  86.  
  87.  
  88. def cunifvariate(mean, arc):
  89.     return (mean + arc * (random() - 0.5)) % pi
  90.  
  91.  
  92. def expovariate(lambd):
  93.     u = random()
  94.     while u <= 1e-07:
  95.         u = random()
  96.     return -log(u) / lambd
  97.  
  98. TWOPI = 2.0 * pi
  99. verify('TWOPI', 6.28318530718)
  100.  
  101. def vonmisesvariate(mu, kappa):
  102.     if kappa <= 1e-06:
  103.         return TWOPI * random()
  104.     
  105.     a = 1.0 + sqrt(1.0 + 4.0 * kappa * kappa)
  106.     b = (a - sqrt(2.0 * a)) / (2.0 * kappa)
  107.     r = (1.0 + b * b) / (2.0 * b)
  108.     while 1:
  109.         u1 = random()
  110.         z = cos(pi * u1)
  111.         f = (1.0 + r * z) / (r + z)
  112.         c = kappa * (r - f)
  113.         u2 = random()
  114.         if u2 >= c * (2.0 - c):
  115.             pass
  116.         if not (u2 > c * exp(1.0 - c)):
  117.             break
  118.         
  119.     u3 = random()
  120.     if u3 > 0.5:
  121.         theta = mu % TWOPI + acos(f)
  122.     else:
  123.         theta = mu % TWOPI - acos(f)
  124.     return theta
  125.  
  126. LOG4 = log(4.0)
  127. verify('LOG4', 1.38629436112)
  128.  
  129. def gammavariate(alpha, beta):
  130.     ainv = sqrt(2.0 * alpha - 1.0)
  131.     return beta * stdgamma(alpha, ainv, alpha - LOG4, alpha + ainv)
  132.  
  133. SG_MAGICCONST = 1.0 + log(4.5)
  134. verify('SG_MAGICCONST', 2.50407739678)
  135.  
  136. def stdgamma(alpha, ainv, bbb, ccc):
  137.     if alpha <= 0.0:
  138.         raise ValueError, 'stdgamma: alpha must be > 0.0'
  139.     
  140.     if alpha > 1.0:
  141.         while 1:
  142.             u1 = random()
  143.             u2 = random()
  144.             v = log(u1 / (1.0 - u1)) / ainv
  145.             x = alpha * exp(v)
  146.             z = u1 * u1 * u2
  147.             r = bbb + ccc * v - x
  148.             if r + SG_MAGICCONST - 4.5 * z >= 0.0 or r >= log(z):
  149.                 return x
  150.             
  151.     elif alpha == 1.0:
  152.         u = random()
  153.         while u <= 1e-07:
  154.             u = random()
  155.         return -log(u)
  156.     else:
  157.         while 1:
  158.             u = random()
  159.             b = (e + alpha) / e
  160.             p = b * u
  161.             if p <= 1.0:
  162.                 x = pow(p, 1.0 / alpha)
  163.             else:
  164.                 x = -log((b - p) / alpha)
  165.             u1 = random()
  166.             if p <= 1.0 and u1 > exp(-x) and p > 1:
  167.                 pass
  168.             if not (u1 > pow(x, alpha - 1.0)):
  169.                 break
  170.             
  171.         return x
  172.  
  173. gauss_next = None
  174.  
  175. def gauss(mu, sigma):
  176.     global gauss_next, gauss_next
  177.     z = gauss_next
  178.     gauss_next = None
  179.     if z is None:
  180.         x2pi = random() * TWOPI
  181.         g2rad = sqrt(-2.0 * log(1.0 - random()))
  182.         z = cos(x2pi) * g2rad
  183.         gauss_next = sin(x2pi) * g2rad
  184.     
  185.     return mu + z * sigma
  186.  
  187.  
  188. def betavariate(alpha, beta):
  189.     y = expovariate(alpha)
  190.     z = expovariate(1.0 / beta)
  191.     return z / (y + z)
  192.  
  193.  
  194. def paretovariate(alpha):
  195.     u = random()
  196.     return 1.0 / pow(u, 1.0 / alpha)
  197.  
  198.  
  199. def weibullvariate(alpha, beta):
  200.     u = random()
  201.     return alpha * pow(-log(u), 1.0 / beta)
  202.  
  203.  
  204. def test(N = 200):
  205.     print 'TWOPI         =', TWOPI
  206.     print 'LOG4          =', LOG4
  207.     print 'NV_MAGICCONST =', NV_MAGICCONST
  208.     print 'SG_MAGICCONST =', SG_MAGICCONST
  209.     test_generator(N, 'random()')
  210.     test_generator(N, 'normalvariate(0.0, 1.0)')
  211.     test_generator(N, 'lognormvariate(0.0, 1.0)')
  212.     test_generator(N, 'cunifvariate(0.0, 1.0)')
  213.     test_generator(N, 'expovariate(1.0)')
  214.     test_generator(N, 'vonmisesvariate(0.0, 1.0)')
  215.     test_generator(N, 'gammavariate(0.5, 1.0)')
  216.     test_generator(N, 'gammavariate(0.9, 1.0)')
  217.     test_generator(N, 'gammavariate(1.0, 1.0)')
  218.     test_generator(N, 'gammavariate(2.0, 1.0)')
  219.     test_generator(N, 'gammavariate(20.0, 1.0)')
  220.     test_generator(N, 'gammavariate(200.0, 1.0)')
  221.     test_generator(N, 'gauss(0.0, 1.0)')
  222.     test_generator(N, 'betavariate(3.0, 3.0)')
  223.     test_generator(N, 'paretovariate(1.0)')
  224.     test_generator(N, 'weibullvariate(1.0, 1.0)')
  225.  
  226.  
  227. def test_generator(n, funccall):
  228.     import time
  229.     print n, 'times', funccall
  230.     code = compile(funccall, funccall, 'eval')
  231.     sum = 0.0
  232.     sqsum = 0.0
  233.     smallest = 10000000000.0
  234.     largest = -10000000000.0
  235.     t0 = time.time()
  236.     for i in range(n):
  237.         x = eval(code)
  238.         sum = sum + x
  239.         sqsum = sqsum + x * x
  240.         smallest = min(x, smallest)
  241.         largest = max(x, largest)
  242.     
  243.     t1 = time.time()
  244.     print round(t1 - t0, 3), 'sec,',
  245.     avg = sum / n
  246.     stddev = sqrt(sqsum / n - avg * avg)
  247.     print 'avg %g, stddev %g, min %g, max %g' % (avg, stddev, smallest, largest)
  248.  
  249.