Den 29.03.2011 16:49, skrev Sturla Molden:
This will not work. A boolean array is not compactly stored, but an array of bytes. Only the first bit 0 is 1 with probability p, bits 1 to 7 bits are 1 with probability 0. We thus have to do this 8 times for each byte, shift left by range(8), and combine them with binary or. Also the main use of random bits is crypto, which requires the use of /dev/urandom or CrypGenRandom instead of Mersenne Twister (np.random.rand).
Here's a cleaned-up one for those who might be interested :-) Sturla import numpy as np import os def randombits(n, p, intention='numerical'): """ Returns an array with packed bits drawn from n Bernoulli trials with successrate p. """ assert (intention in ('numerical','crypto')) # number of bytes m = (n >> 3) + 1 if n % 8 else n >> 3 if intention == 'numerical': # Mersenne Twister rflt = np.random.rand(m*8) else: # /dev/urandom on Linux, Apple, et al., # CryptGenRandom on Windows rflt = np.frombuffer(os.urandom(m*64),dtype=np.uint64) rflt = rflt / float(2**64) b = (rflt.reshape((m,8))<p) # pack the bits b = (b<<range(8)).sum(axis=1).astype(np.uint8) # zero the trailing m*8 - n bits b[-1] &= (0xFF >> (m*8 - n)) return b def bitcount(a, bytewise=False): """ Count the number of set bits in an array of np.uint8. """ assert(a.dtype == np.uint8) c = a[:,np.newaxis].repeat(8, axis=1) >> range(8) & 0x01 return (c.sum(axis=1) if bytewise else c.sum()) if __name__ == '__main__': b = randombits(int(1e6), .1, intent='numerical') print bitcount(b) # should be close to 100000 b = randombits(int(1e6), .1, intent='crypto') print bitcount(b) # should be close to 100000