[Python-checkins] python/nondist/sandbox/twister random.py,1.1,1.2
rhettinger@users.sourceforge.net
rhettinger@users.sourceforge.net
Wed, 18 Dec 2002 10:10:28 -0800
Update of /cvsroot/python/python/nondist/sandbox/twister
In directory sc8-pr-cvs1:/tmp/cvs-serv29268
Modified Files:
random.py
Log Message:
Actually, this checkin in the real Lib/random.py 1.40 used as a basis
for incorporating the MersenneTwister.
Index: random.py
===================================================================
RCS file: /cvsroot/python/python/nondist/sandbox/twister/random.py,v
retrieving revision 1.1
retrieving revision 1.2
diff -C2 -d -r1.1 -r1.2
*** random.py 18 Dec 2002 18:05:01 -0000 1.1
--- random.py 18 Dec 2002 18:10:26 -0000 1.2
***************
*** 108,113 ****
# Adrian Baddeley.
- import MersenneTwister
-
class Random:
"""Random number generator base class used by bound module functions.
--- 108,111 ----
***************
*** 125,149 ****
"""
! VERSION = 1
! def __init__(self, *seeds):
"""Initialize an instance.
! Optional argument "seeds" controls seeding, as for Random.seed().
"""
! if self.__class__ == Random:
! geninstance = MersenneTwister.Random()
! self.random = geninstance.random
! self.seed = geninstance.seed
! self._getstate = geninstance.getstate
! self._setstate = geninstance.setstate
! self._randbelow = geninstance._randbelow
! self.jumpahead = geninstance.jumpahead
self.gauss_next = None
! self.seed(*seeds)
!
def getstate(self):
"""Return internal state; can be passed to setstate() later."""
! return self.VERSION, self._getstate(), self.gauss_next
def setstate(self, state):
--- 123,204 ----
"""
! VERSION = 1 # used by getstate/setstate
! def __init__(self, x=None):
"""Initialize an instance.
! Optional argument x controls seeding, as for Random.seed().
"""
!
! self.seed(x)
!
! ## -------------------- core generator -------------------
!
! # Specific to Wichmann-Hill generator. Subclasses wishing to use a
! # different core generator should override the seed(), random(),
! # getstate(), setstate() and jumpahead() methods.
!
! def seed(self, a=None):
! """Initialize internal state from hashable object.
!
! None or no argument seeds from current time.
!
! If a is not None or an int or long, hash(a) is used instead.
!
! If a is an int or long, a is used directly. Distinct values between
! 0 and 27814431486575L inclusive are guaranteed to yield distinct
! internal states (this guarantee is specific to the default
! Wichmann-Hill generator).
! """
!
! if a is None:
! # Initialize from current time
! import time
! a = long(time.time() * 256)
!
! if type(a) not in (type(3), type(3L)):
! a = hash(a)
!
! a, x = divmod(a, 30268)
! a, y = divmod(a, 30306)
! a, z = divmod(a, 30322)
! self._seed = int(x)+1, int(y)+1, int(z)+1
!
self.gauss_next = None
!
! def random(self):
! """Get the next random number in the range [0.0, 1.0)."""
!
! # Wichman-Hill random number generator.
! #
! # Wichmann, B. A. & Hill, I. D. (1982)
! # Algorithm AS 183:
! # An efficient and portable pseudo-random number generator
! # Applied Statistics 31 (1982) 188-190
! #
! # see also:
! # Correction to Algorithm AS 183
! # Applied Statistics 33 (1984) 123
! #
! # McLeod, A. I. (1985)
! # A remark on Algorithm AS 183
! # Applied Statistics 34 (1985),198-200
!
! # This part is thread-unsafe:
! # BEGIN CRITICAL SECTION
! x, y, z = self._seed
! x = (171 * x) % 30269
! y = (172 * y) % 30307
! z = (170 * z) % 30323
! self._seed = x, y, z
! # END CRITICAL SECTION
!
! # Note: on a platform using IEEE-754 double arithmetic, this can
! # never return 0.0 (asserted by Tim; proof too long for a comment).
! return (x/30269.0 + y/30307.0 + z/30323.0) % 1.0
!
def getstate(self):
"""Return internal state; can be passed to setstate() later."""
! return self.VERSION, self._seed, self.gauss_next
def setstate(self, state):
***************
*** 151,156 ****
version = state[0]
if version == 1:
! version, internalstate, self.gauss_next = state
! self._setstate(internalstate)
else:
raise ValueError("state with version %s passed to "
--- 206,210 ----
version = state[0]
if version == 1:
! version, self._seed, self.gauss_next = state
else:
raise ValueError("state with version %s passed to "
***************
*** 158,161 ****
--- 212,283 ----
(version, self.VERSION))
+ def jumpahead(self, n):
+ """Act as if n calls to random() were made, but quickly.
+
+ n is an int, greater than or equal to 0.
+
+ Example use: If you have 2 threads and know that each will
+ consume no more than a million random numbers, create two Random
+ objects r1 and r2, then do
+ r2.setstate(r1.getstate())
+ r2.jumpahead(1000000)
+ Then r1 and r2 will use guaranteed-disjoint segments of the full
+ period.
+ """
+
+ if not n >= 0:
+ raise ValueError("n must be >= 0")
+ x, y, z = self._seed
+ x = int(x * pow(171, n, 30269)) % 30269
+ y = int(y * pow(172, n, 30307)) % 30307
+ z = int(z * pow(170, n, 30323)) % 30323
+ self._seed = x, y, z
+
+ def __whseed(self, x=0, y=0, z=0):
+ """Set the Wichmann-Hill seed from (x, y, z).
+
+ These must be integers in the range [0, 256).
+ """
+
+ if not type(x) == type(y) == type(z) == int:
+ raise TypeError('seeds must be integers')
+ if not (0 <= x < 256 and 0 <= y < 256 and 0 <= z < 256):
+ raise ValueError('seeds must be in range(0, 256)')
+ if 0 == x == y == z:
+ # Initialize from current time
+ import time
+ t = long(time.time() * 256)
+ t = int((t&0xffffff) ^ (t>>24))
+ t, x = divmod(t, 256)
+ t, y = divmod(t, 256)
+ t, z = divmod(t, 256)
+ # Zero is a poor seed, so substitute 1
+ self._seed = (x or 1, y or 1, z or 1)
+
+ self.gauss_next = None
+
+ def whseed(self, a=None):
+ """Seed from hashable object's hash code.
+
+ None or no argument seeds from current time. It is not guaranteed
+ that objects with distinct hash codes lead to distinct internal
+ states.
+
+ This is obsolete, provided for compatibility with the seed routine
+ used prior to Python 2.1. Use the .seed() method instead.
+ """
+
+ if a is None:
+ self.__whseed()
+ return
+ a = hash(a)
+ a, x = divmod(a, 256)
+ a, y = divmod(a, 256)
+ a, z = divmod(a, 256)
+ x = (x + a) % 256 or 1
+ y = (y + a) % 256 or 1
+ z = (z + a) % 256 or 1
+ self.__whseed(x, y, z)
+
## ---- Methods below this point do not need to be overridden when
## ---- subclassing for the purpose of using a different core generator.
***************
*** 171,178 ****
## -------------------- integer methods -------------------
! def _randbelow(self, n, int=int):
! return int(self.random() * n)
!
! def randrange(self, start, stop=None, step=1, _randbelow=_randbelow, default=None):
"""Choose a random item from range(start, stop[, step]).
--- 293,297 ----
## -------------------- integer methods -------------------
! def randrange(self, start, stop=None, step=1, int=int, default=None):
"""Choose a random item from range(start, stop[, step]).
***************
*** 189,193 ****
if stop is default:
if istart > 0:
! return _randbelow(istart)
raise ValueError, "empty range for randrange()"
--- 308,312 ----
if stop is default:
if istart > 0:
! return int(self.random() * istart)
raise ValueError, "empty range for randrange()"
***************
*** 198,202 ****
if step == 1 and istart < istop:
try:
! return istart + _randbelow(istop - istart)
except OverflowError:
# This can happen if istop-istart > sys.maxint + 1, and
--- 317,321 ----
if step == 1 and istart < istop:
try:
! return istart + int(self.random()*(istop - istart))
except OverflowError:
# This can happen if istop-istart > sys.maxint + 1, and
***************
*** 223,227 ****
if n <= 0:
raise ValueError, "empty range for randrange()"
! return istart + istep*self._randbelow(n)
def randint(self, a, b):
--- 342,346 ----
if n <= 0:
raise ValueError, "empty range for randrange()"
! return istart + istep*int(self.random() * n)
def randint(self, a, b):
***************
*** 237,241 ****
return seq[int(self.random() * len(seq))]
! def shuffle(self, x, _randbelow=None):
"""x, random=random.random -> shuffle list x in place; return None.
--- 356,360 ----
return seq[int(self.random() * len(seq))]
! def shuffle(self, x, random=None, int=int):
"""x, random=random.random -> shuffle list x in place; return None.
***************
*** 249,260 ****
"""
! if _randbelow is None:
! _randbelow = self._randbelow
for i in xrange(len(x)-1, 0, -1):
# pick an element in x[:i+1] with which to exchange x[i]
! j = _randbelow(i+1)
x[i], x[j] = x[j], x[i]
! def sample(self, population, k, _randbelow=None):
"""Chooses k unique random elements from a population sequence.
--- 368,379 ----
"""
! if random is None:
! random = self.random
for i in xrange(len(x)-1, 0, -1):
# pick an element in x[:i+1] with which to exchange x[i]
! j = int(random() * (i+1))
x[i], x[j] = x[j], x[i]
! def sample(self, population, k, random=None, int=int):
"""Chooses k unique random elements from a population sequence.
***************
*** 294,304 ****
if not 0 <= k <= n:
raise ValueError, "sample larger than population"
! if _randbelow is None:
! _randbelow = self._randbelow
result = [None] * k
if n < 6 * k: # if n len list takes less space than a k len dict
pool = list(population)
for i in xrange(k): # invariant: non-selected at [0,n-i)
! j = _randbelow(n-i)
result[i] = pool[j]
pool[j] = pool[n-i-1]
--- 413,423 ----
if not 0 <= k <= n:
raise ValueError, "sample larger than population"
! if random is None:
! random = self.random
result = [None] * k
if n < 6 * k: # if n len list takes less space than a k len dict
pool = list(population)
for i in xrange(k): # invariant: non-selected at [0,n-i)
! j = int(random() * (n-i))
result[i] = pool[j]
pool[j] = pool[n-i-1]
***************
*** 306,312 ****
selected = {}
for i in xrange(k):
! j = _randbelow(n)
while j in selected:
! j = _randbelow(n)
result[i] = selected[j] = population[j]
return result
--- 425,431 ----
selected = {}
for i in xrange(k):
! j = int(random() * n)
while j in selected:
! j = int(random() * n)
result[i] = selected[j] = population[j]
return result
***************
*** 626,780 ****
return alpha * pow(-_log(u), 1.0/beta)
- ## -------------------- Wichmann-Hill -------------------
-
- class WichmannHill(Random):
- # Specific to Wichmann-Hill generator. Subclasses wishing to use a
- # different core generator should override the seed(), random(),
- # getstate(), setstate() and jumpahead() methods.
-
- def seed(self, *seeds):
- """Initialize internal state from hashable object.
-
- None or no argument seeds from current time.
-
- If a is not None or an int or long, hash(a) is used instead.
-
- If a is an int or long, a is used directly. Distinct values between
- 0 and 27814431486575L inclusive are guaranteed to yield distinct
- internal states (this guarantee is specific to the default
- Wichmann-Hill generator).
- """
-
- if len(seeds):
- a = 0
- for s in seeds:
- if not isinstance(s, (int, long)):
- s = hash(s)
- a ^= s
- else:
- # Initialize from current time
- import time
- a = long(time.time() * 256)
-
- a, x = divmod(a, 30268)
- a, y = divmod(a, 30306)
- a, z = divmod(a, 30322)
- self._seed = int(x)+1, int(y)+1, int(z)+1
-
- self.gauss_next = None
-
- def random(self):
- """Get the next random number in the range [0.0, 1.0)."""
-
- # Wichman-Hill random number generator.
- #
- # Wichmann, B. A. & Hill, I. D. (1982)
- # Algorithm AS 183:
- # An efficient and portable pseudo-random number generator
- # Applied Statistics 31 (1982) 188-190
- #
- # see also:
- # Correction to Algorithm AS 183
- # Applied Statistics 33 (1984) 123
- #
- # McLeod, A. I. (1985)
- # A remark on Algorithm AS 183
- # Applied Statistics 34 (1985),198-200
-
- # This part is thread-unsafe:
- # BEGIN CRITICAL SECTION
- x, y, z = self._seed
- x = (171 * x) % 30269
- y = (172 * y) % 30307
- z = (170 * z) % 30323
- self._seed = x, y, z
- # END CRITICAL SECTION
-
- # Note: on a platform using IEEE-754 double arithmetic, this can
- # never return 0.0 (asserted by Tim; proof too long for a comment).
- return (x/30269.0 + y/30307.0 + z/30323.0) % 1.0
-
- def getstate(self):
- """Return internal state; can be passed to setstate() later."""
- return self.VERSION, self._seed, self.gauss_next
-
- def setstate(self, state):
- """Restore internal state from object returned by getstate()."""
- version = state[0]
- if version == 1:
- version, self._seed, self.gauss_next = state
- else:
- raise ValueError("state with version %s passed to "
- "Random.setstate() of version %s" %
- (version, self.VERSION))
-
- def jumpahead(self, n):
- """Act as if n calls to random() were made, but quickly.
-
- n is an int, greater than or equal to 0.
-
- Example use: If you have 2 threads and know that each will
- consume no more than a million random numbers, create two Random
- objects r1 and r2, then do
- r2.setstate(r1.getstate())
- r2.jumpahead(1000000)
- Then r1 and r2 will use guaranteed-disjoint segments of the full
- period.
- """
-
- if not n >= 0:
- raise ValueError("n must be >= 0")
- x, y, z = self._seed
- x = int(x * pow(171, n, 30269)) % 30269
- y = int(y * pow(172, n, 30307)) % 30307
- z = int(z * pow(170, n, 30323)) % 30323
- self._seed = x, y, z
-
- def __whseed(self, x=0, y=0, z=0):
- """Set the Wichmann-Hill seed from (x, y, z).
-
- These must be integers in the range [0, 256).
- """
-
- if not type(x) == type(y) == type(z) == int:
- raise TypeError('seeds must be integers')
- if not (0 <= x < 256 and 0 <= y < 256 and 0 <= z < 256):
- raise ValueError('seeds must be in range(0, 256)')
- if 0 == x == y == z:
- # Initialize from current time
- import time
- t = long(time.time() * 256)
- t = int((t&0xffffff) ^ (t>>24))
- t, x = divmod(t, 256)
- t, y = divmod(t, 256)
- t, z = divmod(t, 256)
- # Zero is a poor seed, so substitute 1
- self._seed = (x or 1, y or 1, z or 1)
-
- self.gauss_next = None
-
- def whseed(self, a=None):
- """Seed from hashable object's hash code.
-
- None or no argument seeds from current time. It is not guaranteed
- that objects with distinct hash codes lead to distinct internal
- states.
-
- This is obsolete, provided for compatibility with the seed routine
- used prior to Python 2.1. Use the .seed() method instead.
- """
-
- if a is None:
- self.__whseed()
- return
- a = hash(a)
- a, x = divmod(a, 256)
- a, y = divmod(a, 256)
- a, z = divmod(a, 256)
- x = (x + a) % 256 or 1
- y = (y + a) % 256 or 1
- z = (z + a) % 256 or 1
- self.__whseed(x, y, z)
-
## -------------------- test program --------------------
--- 745,748 ----
***************
*** 843,859 ****
_test_sample(500)
! if isinstance(_inst, WichmannHill):
! # Test jumpahead.
! s = getstate()
! jumpahead(N)
! r1 = random()
! # now do it the slow way
! setstate(s)
! for i in range(N):
! random()
! r2 = random()
! if r1 != r2:
! raise ValueError("jumpahead test failed " + `(N, r1, r2)`)
! print _inst.__class__
# Create one instance, seeded from current time, and export its methods
--- 811,825 ----
_test_sample(500)
! # Test jumpahead.
! s = getstate()
! jumpahead(N)
! r1 = random()
! # now do it the slow way
! setstate(s)
! for i in range(N):
! random()
! r2 = random()
! if r1 != r2:
! raise ValueError("jumpahead test failed " + `(N, r1, r2)`)
# Create one instance, seeded from current time, and export its methods
***************
*** 862,872 ****
# libraries), but that's fine for most programs and is easier for the
# casual user than making them instantiate their own Random() instance.
-
_inst = Random()
- #_inst = WichmannHill()
seed = _inst.seed
random = _inst.random
uniform = _inst.uniform
- _randbelow = _inst._randbelow
randint = _inst.randint
choice = _inst.choice
--- 828,835 ----
***************
*** 888,895 ****
setstate = _inst.setstate
jumpahead = _inst.jumpahead
! try:
! whseed = _inst.whseed
! except AttributeError:
! pass
if __name__ == '__main__':
--- 851,855 ----
setstate = _inst.setstate
jumpahead = _inst.jumpahead
! whseed = _inst.whseed
if __name__ == '__main__':