test_direct.py 19 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580
  1. import os
  2. from os.path import join
  3. import sys
  4. import numpy as np
  5. from numpy.testing import (assert_equal, assert_allclose, assert_array_equal,
  6. assert_raises)
  7. import pytest
  8. from numpy.random import (
  9. Generator, MT19937, PCG64, PCG64DXSM, Philox, RandomState, SeedSequence,
  10. SFC64, default_rng
  11. )
  12. from numpy.random._common import interface
  13. try:
  14. import cffi # noqa: F401
  15. MISSING_CFFI = False
  16. except ImportError:
  17. MISSING_CFFI = True
  18. try:
  19. import ctypes # noqa: F401
  20. MISSING_CTYPES = False
  21. except ImportError:
  22. MISSING_CTYPES = False
  23. if sys.flags.optimize > 1:
  24. # no docstrings present to inspect when PYTHONOPTIMIZE/Py_OptimizeFlag > 1
  25. # cffi cannot succeed
  26. MISSING_CFFI = True
  27. pwd = os.path.dirname(os.path.abspath(__file__))
  28. def assert_state_equal(actual, target):
  29. for key in actual:
  30. if isinstance(actual[key], dict):
  31. assert_state_equal(actual[key], target[key])
  32. elif isinstance(actual[key], np.ndarray):
  33. assert_array_equal(actual[key], target[key])
  34. else:
  35. assert actual[key] == target[key]
  36. def uint32_to_float32(u):
  37. return ((u >> np.uint32(8)) * (1.0 / 2**24)).astype(np.float32)
  38. def uniform32_from_uint64(x):
  39. x = np.uint64(x)
  40. upper = np.array(x >> np.uint64(32), dtype=np.uint32)
  41. lower = np.uint64(0xffffffff)
  42. lower = np.array(x & lower, dtype=np.uint32)
  43. joined = np.column_stack([lower, upper]).ravel()
  44. return uint32_to_float32(joined)
  45. def uniform32_from_uint53(x):
  46. x = np.uint64(x) >> np.uint64(16)
  47. x = np.uint32(x & np.uint64(0xffffffff))
  48. return uint32_to_float32(x)
  49. def uniform32_from_uint32(x):
  50. return uint32_to_float32(x)
  51. def uniform32_from_uint(x, bits):
  52. if bits == 64:
  53. return uniform32_from_uint64(x)
  54. elif bits == 53:
  55. return uniform32_from_uint53(x)
  56. elif bits == 32:
  57. return uniform32_from_uint32(x)
  58. else:
  59. raise NotImplementedError
  60. def uniform_from_uint(x, bits):
  61. if bits in (64, 63, 53):
  62. return uniform_from_uint64(x)
  63. elif bits == 32:
  64. return uniform_from_uint32(x)
  65. def uniform_from_uint64(x):
  66. return (x >> np.uint64(11)) * (1.0 / 9007199254740992.0)
  67. def uniform_from_uint32(x):
  68. out = np.empty(len(x) // 2)
  69. for i in range(0, len(x), 2):
  70. a = x[i] >> 5
  71. b = x[i + 1] >> 6
  72. out[i // 2] = (a * 67108864.0 + b) / 9007199254740992.0
  73. return out
  74. def uniform_from_dsfmt(x):
  75. return x.view(np.double) - 1.0
  76. def gauss_from_uint(x, n, bits):
  77. if bits in (64, 63):
  78. doubles = uniform_from_uint64(x)
  79. elif bits == 32:
  80. doubles = uniform_from_uint32(x)
  81. else: # bits == 'dsfmt'
  82. doubles = uniform_from_dsfmt(x)
  83. gauss = []
  84. loc = 0
  85. x1 = x2 = 0.0
  86. while len(gauss) < n:
  87. r2 = 2
  88. while r2 >= 1.0 or r2 == 0.0:
  89. x1 = 2.0 * doubles[loc] - 1.0
  90. x2 = 2.0 * doubles[loc + 1] - 1.0
  91. r2 = x1 * x1 + x2 * x2
  92. loc += 2
  93. f = np.sqrt(-2.0 * np.log(r2) / r2)
  94. gauss.append(f * x2)
  95. gauss.append(f * x1)
  96. return gauss[:n]
  97. def test_seedsequence():
  98. from numpy.random.bit_generator import (ISeedSequence,
  99. ISpawnableSeedSequence,
  100. SeedlessSeedSequence)
  101. s1 = SeedSequence(range(10), spawn_key=(1, 2), pool_size=6)
  102. s1.spawn(10)
  103. s2 = SeedSequence(**s1.state)
  104. assert_equal(s1.state, s2.state)
  105. assert_equal(s1.n_children_spawned, s2.n_children_spawned)
  106. # The interfaces cannot be instantiated themselves.
  107. assert_raises(TypeError, ISeedSequence)
  108. assert_raises(TypeError, ISpawnableSeedSequence)
  109. dummy = SeedlessSeedSequence()
  110. assert_raises(NotImplementedError, dummy.generate_state, 10)
  111. assert len(dummy.spawn(10)) == 10
  112. def test_generator_spawning():
  113. """ Test spawning new generators and bit_generators directly.
  114. """
  115. rng = np.random.default_rng()
  116. seq = rng.bit_generator.seed_seq
  117. new_ss = seq.spawn(5)
  118. expected_keys = [seq.spawn_key + (i,) for i in range(5)]
  119. assert [c.spawn_key for c in new_ss] == expected_keys
  120. new_bgs = rng.bit_generator.spawn(5)
  121. expected_keys = [seq.spawn_key + (i,) for i in range(5, 10)]
  122. assert [bg.seed_seq.spawn_key for bg in new_bgs] == expected_keys
  123. new_rngs = rng.spawn(5)
  124. expected_keys = [seq.spawn_key + (i,) for i in range(10, 15)]
  125. found_keys = [rng.bit_generator.seed_seq.spawn_key for rng in new_rngs]
  126. assert found_keys == expected_keys
  127. # Sanity check that streams are actually different:
  128. assert new_rngs[0].uniform() != new_rngs[1].uniform()
  129. def test_non_spawnable():
  130. from numpy.random.bit_generator import ISeedSequence
  131. class FakeSeedSequence:
  132. def generate_state(self, n_words, dtype=np.uint32):
  133. return np.zeros(n_words, dtype=dtype)
  134. ISeedSequence.register(FakeSeedSequence)
  135. rng = np.random.default_rng(FakeSeedSequence())
  136. with pytest.raises(TypeError, match="The underlying SeedSequence"):
  137. rng.spawn(5)
  138. with pytest.raises(TypeError, match="The underlying SeedSequence"):
  139. rng.bit_generator.spawn(5)
  140. class Base:
  141. dtype = np.uint64
  142. data2 = data1 = {}
  143. @classmethod
  144. def setup_class(cls):
  145. cls.bit_generator = PCG64
  146. cls.bits = 64
  147. cls.dtype = np.uint64
  148. cls.seed_error_type = TypeError
  149. cls.invalid_init_types = []
  150. cls.invalid_init_values = []
  151. @classmethod
  152. def _read_csv(cls, filename):
  153. with open(filename) as csv:
  154. seed = csv.readline()
  155. seed = seed.split(',')
  156. seed = [int(s.strip(), 0) for s in seed[1:]]
  157. data = []
  158. for line in csv:
  159. data.append(int(line.split(',')[-1].strip(), 0))
  160. return {'seed': seed, 'data': np.array(data, dtype=cls.dtype)}
  161. def test_raw(self):
  162. bit_generator = self.bit_generator(*self.data1['seed'])
  163. uints = bit_generator.random_raw(1000)
  164. assert_equal(uints, self.data1['data'])
  165. bit_generator = self.bit_generator(*self.data1['seed'])
  166. uints = bit_generator.random_raw()
  167. assert_equal(uints, self.data1['data'][0])
  168. bit_generator = self.bit_generator(*self.data2['seed'])
  169. uints = bit_generator.random_raw(1000)
  170. assert_equal(uints, self.data2['data'])
  171. def test_random_raw(self):
  172. bit_generator = self.bit_generator(*self.data1['seed'])
  173. uints = bit_generator.random_raw(output=False)
  174. assert uints is None
  175. uints = bit_generator.random_raw(1000, output=False)
  176. assert uints is None
  177. def test_gauss_inv(self):
  178. n = 25
  179. rs = RandomState(self.bit_generator(*self.data1['seed']))
  180. gauss = rs.standard_normal(n)
  181. assert_allclose(gauss,
  182. gauss_from_uint(self.data1['data'], n, self.bits))
  183. rs = RandomState(self.bit_generator(*self.data2['seed']))
  184. gauss = rs.standard_normal(25)
  185. assert_allclose(gauss,
  186. gauss_from_uint(self.data2['data'], n, self.bits))
  187. def test_uniform_double(self):
  188. rs = Generator(self.bit_generator(*self.data1['seed']))
  189. vals = uniform_from_uint(self.data1['data'], self.bits)
  190. uniforms = rs.random(len(vals))
  191. assert_allclose(uniforms, vals)
  192. assert_equal(uniforms.dtype, np.float64)
  193. rs = Generator(self.bit_generator(*self.data2['seed']))
  194. vals = uniform_from_uint(self.data2['data'], self.bits)
  195. uniforms = rs.random(len(vals))
  196. assert_allclose(uniforms, vals)
  197. assert_equal(uniforms.dtype, np.float64)
  198. def test_uniform_float(self):
  199. rs = Generator(self.bit_generator(*self.data1['seed']))
  200. vals = uniform32_from_uint(self.data1['data'], self.bits)
  201. uniforms = rs.random(len(vals), dtype=np.float32)
  202. assert_allclose(uniforms, vals)
  203. assert_equal(uniforms.dtype, np.float32)
  204. rs = Generator(self.bit_generator(*self.data2['seed']))
  205. vals = uniform32_from_uint(self.data2['data'], self.bits)
  206. uniforms = rs.random(len(vals), dtype=np.float32)
  207. assert_allclose(uniforms, vals)
  208. assert_equal(uniforms.dtype, np.float32)
  209. def test_repr(self):
  210. rs = Generator(self.bit_generator(*self.data1['seed']))
  211. assert 'Generator' in repr(rs)
  212. assert f'{id(rs):#x}'.upper().replace('X', 'x') in repr(rs)
  213. def test_str(self):
  214. rs = Generator(self.bit_generator(*self.data1['seed']))
  215. assert 'Generator' in str(rs)
  216. assert str(self.bit_generator.__name__) in str(rs)
  217. assert f'{id(rs):#x}'.upper().replace('X', 'x') not in str(rs)
  218. def test_pickle(self):
  219. import pickle
  220. bit_generator = self.bit_generator(*self.data1['seed'])
  221. state = bit_generator.state
  222. bitgen_pkl = pickle.dumps(bit_generator)
  223. reloaded = pickle.loads(bitgen_pkl)
  224. reloaded_state = reloaded.state
  225. assert_array_equal(Generator(bit_generator).standard_normal(1000),
  226. Generator(reloaded).standard_normal(1000))
  227. assert bit_generator is not reloaded
  228. assert_state_equal(reloaded_state, state)
  229. ss = SeedSequence(100)
  230. aa = pickle.loads(pickle.dumps(ss))
  231. assert_equal(ss.state, aa.state)
  232. def test_pickle_preserves_seed_sequence(self):
  233. # GH 26234
  234. # Add explicit test that bit generators preserve seed sequences
  235. import pickle
  236. bit_generator = self.bit_generator(*self.data1['seed'])
  237. ss = bit_generator.seed_seq
  238. bg_plk = pickle.loads(pickle.dumps(bit_generator))
  239. ss_plk = bg_plk.seed_seq
  240. assert_equal(ss.state, ss_plk.state)
  241. assert_equal(ss.pool, ss_plk.pool)
  242. bit_generator.seed_seq.spawn(10)
  243. bg_plk = pickle.loads(pickle.dumps(bit_generator))
  244. ss_plk = bg_plk.seed_seq
  245. assert_equal(ss.state, ss_plk.state)
  246. assert_equal(ss.n_children_spawned, ss_plk.n_children_spawned)
  247. def test_invalid_state_type(self):
  248. bit_generator = self.bit_generator(*self.data1['seed'])
  249. with pytest.raises(TypeError):
  250. bit_generator.state = {'1'}
  251. def test_invalid_state_value(self):
  252. bit_generator = self.bit_generator(*self.data1['seed'])
  253. state = bit_generator.state
  254. state['bit_generator'] = 'otherBitGenerator'
  255. with pytest.raises(ValueError):
  256. bit_generator.state = state
  257. def test_invalid_init_type(self):
  258. bit_generator = self.bit_generator
  259. for st in self.invalid_init_types:
  260. with pytest.raises(TypeError):
  261. bit_generator(*st)
  262. def test_invalid_init_values(self):
  263. bit_generator = self.bit_generator
  264. for st in self.invalid_init_values:
  265. with pytest.raises((ValueError, OverflowError)):
  266. bit_generator(*st)
  267. def test_benchmark(self):
  268. bit_generator = self.bit_generator(*self.data1['seed'])
  269. bit_generator._benchmark(1)
  270. bit_generator._benchmark(1, 'double')
  271. with pytest.raises(ValueError):
  272. bit_generator._benchmark(1, 'int32')
  273. @pytest.mark.skipif(MISSING_CFFI, reason='cffi not available')
  274. def test_cffi(self):
  275. bit_generator = self.bit_generator(*self.data1['seed'])
  276. cffi_interface = bit_generator.cffi
  277. assert isinstance(cffi_interface, interface)
  278. other_cffi_interface = bit_generator.cffi
  279. assert other_cffi_interface is cffi_interface
  280. @pytest.mark.skipif(MISSING_CTYPES, reason='ctypes not available')
  281. def test_ctypes(self):
  282. bit_generator = self.bit_generator(*self.data1['seed'])
  283. ctypes_interface = bit_generator.ctypes
  284. assert isinstance(ctypes_interface, interface)
  285. other_ctypes_interface = bit_generator.ctypes
  286. assert other_ctypes_interface is ctypes_interface
  287. def test_getstate(self):
  288. bit_generator = self.bit_generator(*self.data1['seed'])
  289. state = bit_generator.state
  290. alt_state = bit_generator.__getstate__()
  291. assert isinstance(alt_state, tuple)
  292. assert_state_equal(state, alt_state[0])
  293. assert isinstance(alt_state[1], SeedSequence)
  294. class TestPhilox(Base):
  295. @classmethod
  296. def setup_class(cls):
  297. cls.bit_generator = Philox
  298. cls.bits = 64
  299. cls.dtype = np.uint64
  300. cls.data1 = cls._read_csv(
  301. join(pwd, './data/philox-testset-1.csv'))
  302. cls.data2 = cls._read_csv(
  303. join(pwd, './data/philox-testset-2.csv'))
  304. cls.seed_error_type = TypeError
  305. cls.invalid_init_types = []
  306. cls.invalid_init_values = [(1, None, 1), (-1,), (None, None, 2 ** 257 + 1)]
  307. def test_set_key(self):
  308. bit_generator = self.bit_generator(*self.data1['seed'])
  309. state = bit_generator.state
  310. keyed = self.bit_generator(counter=state['state']['counter'],
  311. key=state['state']['key'])
  312. assert_state_equal(bit_generator.state, keyed.state)
  313. class TestPCG64(Base):
  314. @classmethod
  315. def setup_class(cls):
  316. cls.bit_generator = PCG64
  317. cls.bits = 64
  318. cls.dtype = np.uint64
  319. cls.data1 = cls._read_csv(join(pwd, './data/pcg64-testset-1.csv'))
  320. cls.data2 = cls._read_csv(join(pwd, './data/pcg64-testset-2.csv'))
  321. cls.seed_error_type = (ValueError, TypeError)
  322. cls.invalid_init_types = [(3.2,), ([None],), (1, None)]
  323. cls.invalid_init_values = [(-1,)]
  324. def test_advance_symmetry(self):
  325. rs = Generator(self.bit_generator(*self.data1['seed']))
  326. state = rs.bit_generator.state
  327. step = -0x9e3779b97f4a7c150000000000000000
  328. rs.bit_generator.advance(step)
  329. val_neg = rs.integers(10)
  330. rs.bit_generator.state = state
  331. rs.bit_generator.advance(2**128 + step)
  332. val_pos = rs.integers(10)
  333. rs.bit_generator.state = state
  334. rs.bit_generator.advance(10 * 2**128 + step)
  335. val_big = rs.integers(10)
  336. assert val_neg == val_pos
  337. assert val_big == val_pos
  338. def test_advange_large(self):
  339. rs = Generator(self.bit_generator(38219308213743))
  340. pcg = rs.bit_generator
  341. state = pcg.state["state"]
  342. initial_state = 287608843259529770491897792873167516365
  343. assert state["state"] == initial_state
  344. pcg.advance(sum(2**i for i in (96, 64, 32, 16, 8, 4, 2, 1)))
  345. state = pcg.state["state"]
  346. advanced_state = 135275564607035429730177404003164635391
  347. assert state["state"] == advanced_state
  348. class TestPCG64DXSM(Base):
  349. @classmethod
  350. def setup_class(cls):
  351. cls.bit_generator = PCG64DXSM
  352. cls.bits = 64
  353. cls.dtype = np.uint64
  354. cls.data1 = cls._read_csv(join(pwd, './data/pcg64dxsm-testset-1.csv'))
  355. cls.data2 = cls._read_csv(join(pwd, './data/pcg64dxsm-testset-2.csv'))
  356. cls.seed_error_type = (ValueError, TypeError)
  357. cls.invalid_init_types = [(3.2,), ([None],), (1, None)]
  358. cls.invalid_init_values = [(-1,)]
  359. def test_advance_symmetry(self):
  360. rs = Generator(self.bit_generator(*self.data1['seed']))
  361. state = rs.bit_generator.state
  362. step = -0x9e3779b97f4a7c150000000000000000
  363. rs.bit_generator.advance(step)
  364. val_neg = rs.integers(10)
  365. rs.bit_generator.state = state
  366. rs.bit_generator.advance(2**128 + step)
  367. val_pos = rs.integers(10)
  368. rs.bit_generator.state = state
  369. rs.bit_generator.advance(10 * 2**128 + step)
  370. val_big = rs.integers(10)
  371. assert val_neg == val_pos
  372. assert val_big == val_pos
  373. def test_advange_large(self):
  374. rs = Generator(self.bit_generator(38219308213743))
  375. pcg = rs.bit_generator
  376. state = pcg.state
  377. initial_state = 287608843259529770491897792873167516365
  378. assert state["state"]["state"] == initial_state
  379. pcg.advance(sum(2**i for i in (96, 64, 32, 16, 8, 4, 2, 1)))
  380. state = pcg.state["state"]
  381. advanced_state = 277778083536782149546677086420637664879
  382. assert state["state"] == advanced_state
  383. class TestMT19937(Base):
  384. @classmethod
  385. def setup_class(cls):
  386. cls.bit_generator = MT19937
  387. cls.bits = 32
  388. cls.dtype = np.uint32
  389. cls.data1 = cls._read_csv(join(pwd, './data/mt19937-testset-1.csv'))
  390. cls.data2 = cls._read_csv(join(pwd, './data/mt19937-testset-2.csv'))
  391. cls.seed_error_type = ValueError
  392. cls.invalid_init_types = []
  393. cls.invalid_init_values = [(-1,)]
  394. def test_seed_float_array(self):
  395. assert_raises(TypeError, self.bit_generator, np.array([np.pi]))
  396. assert_raises(TypeError, self.bit_generator, np.array([-np.pi]))
  397. assert_raises(TypeError, self.bit_generator, np.array([np.pi, -np.pi]))
  398. assert_raises(TypeError, self.bit_generator, np.array([0, np.pi]))
  399. assert_raises(TypeError, self.bit_generator, [np.pi])
  400. assert_raises(TypeError, self.bit_generator, [0, np.pi])
  401. def test_state_tuple(self):
  402. rs = Generator(self.bit_generator(*self.data1['seed']))
  403. bit_generator = rs.bit_generator
  404. state = bit_generator.state
  405. desired = rs.integers(2 ** 16)
  406. tup = (state['bit_generator'], state['state']['key'],
  407. state['state']['pos'])
  408. bit_generator.state = tup
  409. actual = rs.integers(2 ** 16)
  410. assert_equal(actual, desired)
  411. tup = tup + (0, 0.0)
  412. bit_generator.state = tup
  413. actual = rs.integers(2 ** 16)
  414. assert_equal(actual, desired)
  415. class TestSFC64(Base):
  416. @classmethod
  417. def setup_class(cls):
  418. cls.bit_generator = SFC64
  419. cls.bits = 64
  420. cls.dtype = np.uint64
  421. cls.data1 = cls._read_csv(
  422. join(pwd, './data/sfc64-testset-1.csv'))
  423. cls.data2 = cls._read_csv(
  424. join(pwd, './data/sfc64-testset-2.csv'))
  425. cls.seed_error_type = (ValueError, TypeError)
  426. cls.invalid_init_types = [(3.2,), ([None],), (1, None)]
  427. cls.invalid_init_values = [(-1,)]
  428. def test_legacy_pickle(self):
  429. # Pickling format was changed in 2.0.x
  430. import gzip
  431. import pickle
  432. expected_state = np.array(
  433. [
  434. 9957867060933711493,
  435. 532597980065565856,
  436. 14769588338631205282,
  437. 13
  438. ],
  439. dtype=np.uint64
  440. )
  441. base_path = os.path.split(os.path.abspath(__file__))[0]
  442. pkl_file = os.path.join(base_path, "data", "sfc64_np126.pkl.gz")
  443. with gzip.open(pkl_file) as gz:
  444. sfc = pickle.load(gz)
  445. assert isinstance(sfc, SFC64)
  446. assert_equal(sfc.state["state"]["state"], expected_state)
  447. class TestDefaultRNG:
  448. def test_seed(self):
  449. for args in [(), (None,), (1234,), ([1234, 5678],)]:
  450. rg = default_rng(*args)
  451. assert isinstance(rg.bit_generator, PCG64)
  452. def test_passthrough(self):
  453. bg = Philox()
  454. rg = default_rng(bg)
  455. assert rg.bit_generator is bg
  456. rg2 = default_rng(rg)
  457. assert rg2 is rg
  458. assert rg2.bit_generator is bg
  459. def test_coercion_RandomState_Generator(self):
  460. # use default_rng to coerce RandomState to Generator
  461. rs = RandomState(1234)
  462. rg = default_rng(rs)
  463. assert isinstance(rg.bit_generator, MT19937)
  464. assert rg.bit_generator is rs._bit_generator
  465. # RandomState with a non MT19937 bit generator
  466. _original = np.random.get_bit_generator()
  467. bg = PCG64(12342298)
  468. np.random.set_bit_generator(bg)
  469. rs = np.random.mtrand._rand
  470. rg = default_rng(rs)
  471. assert rg.bit_generator is bg
  472. # vital to get global state back to original, otherwise
  473. # other tests start to fail.
  474. np.random.set_bit_generator(_original)