Skip to content

Commit 39bda75

Browse files
committed
REF: Refactor bounded integers
Refactor to improve DRY using a common generator for the scalar cases
1 parent f58ed17 commit 39bda75

File tree

2 files changed

+101
-141
lines changed

2 files changed

+101
-141
lines changed

randomstate/bounded_integers.pxi.in

Lines changed: 99 additions & 139 deletions
Original file line numberDiff line numberDiff line change
@@ -12,6 +12,7 @@ _randint_type = {'bool': (0, 2),
1212
ctypedef np.npy_bool bool_t
1313

1414
cdef inline uint64_t _gen_mask(uint64_t max_val) nogil:
15+
"""Mask generator for use in bounded random numbers"""
1516
# Smallest bit mask >= max
1617
cdef uint64_t mask = max_val
1718
mask |= mask >> 1
@@ -21,53 +22,23 @@ cdef inline uint64_t _gen_mask(uint64_t max_val) nogil:
2122
mask |= mask >> 16
2223
mask |= mask >> 32
2324
return mask
25+
26+
2427
{{
2528
py:
26-
bc_ctypes = (('uint32', 'uint32', 'uint64', 'NPY_UINT64', 0, 0, 0, '0X100000000ULL'),
29+
type_info = (('uint32', 'uint32', 'uint64', 'NPY_UINT64', 0, 0, 0, '0X100000000ULL'),
2730
('uint16', 'uint16', 'uint32', 'NPY_UINT32', 1, 16, 0, '0X10000UL'),
2831
('uint8', 'uint8', 'uint16', 'NPY_UINT16', 3, 8, 0, '0X100UL'),
2932
('bool','bool', 'uint8', 'NPY_UINT8', 31, 1, 0, '0x2UL'),
3033
('int32', 'uint32', 'uint64', 'NPY_INT64', 0, 0, '-0x80000000LL', '0x80000000LL'),
3134
('int16', 'uint16', 'uint32', 'NPY_INT32', 1, 16, '-0x8000LL', '0x8000LL' ),
3235
('int8', 'uint8', 'uint16', 'NPY_INT16', 3, 8, '-0x80LL', '0x80LL' ),
3336
)}}
34-
{{for nptype, utype, nptype_up, npctype, remaining, bitshift, lb, ub in bc_ctypes}}
37+
{{for nptype, utype, nptype_up, npctype, remaining, bitshift, lb, ub in type_info}}
3538
{{ py: otype = nptype + '_' if nptype == 'bool' else nptype }}
36-
cdef object _rand_{{nptype}}(object low, object high, object size, aug_state *state, object lock):
37-
"""
38-
_rand_{{nptype}}(low, high, size, *state, lock)
39-
40-
Return random np.{{nptype}} integers between `low` and `high`, inclusive.
41-
42-
Return random integers from the "discrete uniform" distribution in the
43-
closed interval [`low`, `high`). If `high` is None (the default),
44-
then results are from [0, `low`). On entry the arguments are presumed
45-
to have been validated for size and order for the np.{{nptype}} type.
46-
47-
Parameters
48-
----------
49-
low : int or array-like
50-
Lowest (signed) integer to be drawn from the distribution (unless
51-
``high=None``, in which case this parameter is the *highest* such
52-
integer).
53-
high : int or array-like
54-
If provided, the largest (signed) integer to be drawn from the
55-
distribution (see above for behavior if ``high=None``).
56-
size : int or tuple of ints
57-
Output shape. If the given shape is, e.g., ``(m, n, k)``, then
58-
``m * n * k`` samples are drawn. Default is None, in which case a
59-
single value is returned.
60-
state : augmented random state
61-
State to use in the core random number generators
62-
lock : threading.Lock
63-
Lock to prevent multiple using a single RandomState simultaneously
6439

65-
Returns
66-
-------
67-
out : python scalar or ndarray of np.{{nptype}}
68-
`size`-shaped array of random integers from the appropriate
69-
distribution, or a single such random int if `size` not provided.
70-
"""
40+
cdef object _rand_{{nptype}}_broadcast(np.ndarray low, np.ndarray high, object size, aug_state *state, object lock):
41+
"""Array path for smaller integer types"""
7142
cdef {{utype}}_t rng, last_rng, off, val, mask, out_val
7243
cdef uint32_t buf
7344
cdef {{utype}}_t *out_data
@@ -77,40 +48,6 @@ cdef object _rand_{{nptype}}(object low, object high, object size, aug_state *st
7748
cdef np.broadcast it
7849
cdef int buf_rem = 0
7950

80-
if size is not None:
81-
if (np.prod(size) == 0):
82-
return np.empty(size, dtype=np.{{nptype}})
83-
84-
low = np.array(low, copy=False)
85-
high = np.array(high, copy=False)
86-
low_ndim = np.PyArray_NDIM(<np.ndarray>low)
87-
high_ndim = np.PyArray_NDIM(<np.ndarray>high)
88-
if ((low_ndim == 0 or (low_ndim==1 and low.size==1 and size is not None)) and
89-
(high_ndim == 0 or (high_ndim==1 and high.size==1 and size is not None))):
90-
low = int(low)
91-
high = int(high)
92-
93-
if low < {{lb}}:
94-
raise ValueError("low is out of bounds for {{nptype}}")
95-
if high > {{ub}}:
96-
raise ValueError("high is out of bounds for {{nptype}}")
97-
if low >= high:
98-
raise ValueError("low >= high")
99-
100-
high -= 1
101-
rng = <{{utype}}_t>(high - low)
102-
off = <{{utype}}_t>(<{{nptype}}_t>low)
103-
if size is None:
104-
with lock:
105-
random_bounded_{{utype}}_fill(state, off, rng, 1, &out_val)
106-
return np.{{otype}}(<{{nptype}}_t>out_val)
107-
else:
108-
out_arr = <np.ndarray>np.empty(size, np.{{nptype}})
109-
cnt = np.PyArray_SIZE(out_arr)
110-
out_data = <{{utype}}_t *>np.PyArray_DATA(out_arr)
111-
with lock, nogil:
112-
random_bounded_{{utype}}_fill(state, off, rng, cnt, out_data)
113-
return out_arr
11451

11552
# Array path
11653
low_arr = <np.ndarray>low
@@ -149,16 +86,91 @@ cdef object _rand_{{nptype}}(object low, object high, object size, aug_state *st
14986
out_data[i] = random_buffered_bounded_{{utype}}(state, off, rng, mask, &buf_rem, &buf)
15087

15188
np.PyArray_MultiIter_NEXT(it)
152-
15389
return out_arr
15490
{{endfor}}
91+
15592
{{
15693
py:
157-
big_bc_ctypes = (('uint64', 'uint64', 'NPY_UINT64', '0x0ULL', '0xFFFFFFFFFFFFFFFFULL'),
94+
big_type_info = (('uint64', 'uint64', 'NPY_UINT64', '0x0ULL', '0xFFFFFFFFFFFFFFFFULL'),
15895
('int64', 'uint64', 'NPY_INT64', '-0x8000000000000000LL', '0x7FFFFFFFFFFFFFFFLL' )
15996
)}}
160-
{{for nptype, utype, npctype, lb, ub in big_bc_ctypes}}
97+
{{for nptype, utype, npctype, lb, ub in big_type_info}}
16198
{{ py: otype = nptype}}
99+
cdef object _rand_{{nptype}}_broadcast(object low, object high, object size, aug_state *state, object lock):
100+
"""Array path for 64-bit integer types"""
101+
cdef np.ndarray low_arr, high_arr, out_arr, highm1_arr
102+
cdef np.npy_intp i, cnt
103+
cdef np.broadcast it
104+
cdef object closed_upper
105+
cdef uint64_t *out_data
106+
cdef {{nptype}}_t *highm1_data
107+
cdef {{nptype}}_t low_v, high_v
108+
cdef uint64_t rng, last_rng, val, mask, off, out_val
109+
110+
low_arr = <np.ndarray>low
111+
high_arr = <np.ndarray>high
112+
113+
if np.any(np.less(low_arr, {{lb}})):
114+
raise ValueError('low is out of bounds for {{nptype}}')
115+
116+
highm1_arr = <np.ndarray>np.empty_like(high_arr, dtype=np.{{nptype}})
117+
highm1_data = <{{nptype}}_t *>np.PyArray_DATA(highm1_arr)
118+
cnt = np.PyArray_SIZE(high_arr)
119+
flat = high_arr.flat
120+
for i in range(cnt):
121+
closed_upper = int(flat[i]) - 1
122+
if closed_upper > {{ub}}:
123+
raise ValueError('high is out of bounds for {{nptype}}')
124+
if closed_upper < {{lb}}:
125+
raise ValueError('low >= high')
126+
highm1_data[i] = <{{nptype}}_t>closed_upper
127+
128+
if np.any(np.greater(low_arr, highm1_arr)):
129+
raise ValueError('low >= high')
130+
131+
high_arr = highm1_arr
132+
low_arr = <np.ndarray>np.PyArray_FROM_OTF(low, np.{{npctype}}, np.NPY_ALIGNED | np.NPY_FORCECAST)
133+
134+
if size is not None:
135+
out_arr = <np.ndarray>np.empty(size, np.{{nptype}})
136+
else:
137+
it = np.PyArray_MultiIterNew2(low_arr, high_arr)
138+
out_arr = <np.ndarray>np.empty(it.shape, np.{{nptype}})
139+
140+
it = np.PyArray_MultiIterNew3(low_arr, high_arr, out_arr)
141+
out_data = <uint64_t *>np.PyArray_DATA(out_arr)
142+
n = np.PyArray_SIZE(out_arr)
143+
mask = last_rng = 0
144+
with lock, nogil:
145+
for i in range(n):
146+
low_v = (<{{nptype}}_t*>np.PyArray_MultiIter_DATA(it, 0))[0]
147+
high_v = (<{{nptype}}_t*>np.PyArray_MultiIter_DATA(it, 1))[0]
148+
rng = <{{utype}}_t>(high_v - low_v) # No -1 here since implemented above
149+
off = <{{utype}}_t>(<{{nptype}}_t>low_v)
150+
151+
if rng != last_rng:
152+
mask = _gen_mask(rng)
153+
out_data[i] = random_bounded_uint64(state, off, rng, mask)
154+
155+
np.PyArray_MultiIter_NEXT(it)
156+
157+
return out_arr
158+
{{endfor}}
159+
160+
{{
161+
py:
162+
type_info = (('uint64', 'uint64', '0x0ULL', '0xFFFFFFFFFFFFFFFFULL'),
163+
('uint32', 'uint32', '0x0UL', '0XFFFFFFFFUL'),
164+
('uint16', 'uint16', '0x0UL', '0XFFFFUL'),
165+
('uint8', 'uint8', '0x0UL', '0XFFUL'),
166+
('bool', 'bool', '0x0UL', '0x1UL'),
167+
('int64', 'uint64', '-0x8000000000000000LL', '0x7FFFFFFFFFFFFFFFL'),
168+
('int32', 'uint32', '-0x80000000L', '0x7FFFFFFFL'),
169+
('int16', 'uint16', '-0x8000L', '0x7FFFL' ),
170+
('int8', 'uint8', '-0x80L', '0x7FL' )
171+
)}}
172+
{{for nptype, utype, lb, ub in type_info}}
173+
{{ py: otype = nptype + '_' if nptype == 'bool' else nptype }}
162174
cdef object _rand_{{nptype}}(object low, object high, object size, aug_state *state, object lock):
163175
"""
164176
_rand_{{nptype}}(low, high, size, *state, lock)
@@ -194,34 +206,30 @@ cdef object _rand_{{nptype}}(object low, object high, object size, aug_state *st
194206
`size`-shaped array of random integers from the appropriate
195207
distribution, or a single such random int if `size` not provided.
196208
"""
197-
cdef np.ndarray low_arr, high_arr, out_arr, highm1_arr
209+
cdef np.ndarray out_arr, low_arr, high_arr
210+
cdef {{utype}}_t rng, off, out_val
211+
cdef {{utype}}_t *out_data
198212
cdef np.npy_intp i, cnt
199-
cdef np.broadcast it
200-
cdef object closed_upper
201-
cdef uint64_t *out_data
202-
cdef {{nptype}}_t *highm1_data
203-
cdef {{nptype}}_t low_v, high_v
204-
cdef uint64_t rng, last_rng, val, mask, off, out_val
205213

206214
if size is not None:
207215
if (np.prod(size) == 0):
208216
return np.empty(size, dtype=np.{{nptype}})
209217

210-
low = np.array(low, copy=False)
211-
high = np.array(high, copy=False)
212-
low_ndim = np.PyArray_NDIM(<np.ndarray>low)
213-
high_ndim = np.PyArray_NDIM(<np.ndarray>high)
214-
if ((low_ndim == 0 or (low_ndim==1 and low.size==1 and size is not None)) and
215-
(high_ndim == 0 or (high_ndim==1 and high.size==1 and size is not None))):
216-
low = int(low)
217-
high = int(high)
218-
high -= 1 # Use a closed interval
218+
low_arr = <np.ndarray>np.array(low, copy=False)
219+
high_arr = <np.ndarray>np.array(high, copy=False)
220+
low_ndim = np.PyArray_NDIM(low_arr)
221+
high_ndim = np.PyArray_NDIM(high_arr)
222+
if ((low_ndim == 0 or (low_ndim==1 and low_arr.size==1 and size is not None)) and
223+
(high_ndim == 0 or (high_ndim==1 and high_arr.size==1 and size is not None))):
224+
low = int(low_arr)
225+
high = int(high_arr)
226+
high -= 1
219227

220228
if low < {{lb}}:
221229
raise ValueError("low is out of bounds for {{nptype}}")
222230
if high > {{ub}}:
223231
raise ValueError("high is out of bounds for {{nptype}}")
224-
if low > high:
232+
if low > high: # -1 already subtracted, closed interval
225233
raise ValueError("low >= high")
226234

227235
rng = <{{utype}}_t>(high - low)
@@ -237,53 +245,5 @@ cdef object _rand_{{nptype}}(object low, object high, object size, aug_state *st
237245
with lock, nogil:
238246
random_bounded_{{utype}}_fill(state, off, rng, cnt, out_data)
239247
return out_arr
240-
241-
low_arr = <np.ndarray>low
242-
high_arr = <np.ndarray>high
243-
244-
if np.any(np.less(low_arr, {{lb}})):
245-
raise ValueError('low is out of bounds for {{nptype}}')
246-
247-
highm1_arr = <np.ndarray>np.empty_like(high_arr, dtype=np.{{nptype}})
248-
highm1_data = <{{nptype}}_t *>np.PyArray_DATA(highm1_arr)
249-
cnt = np.PyArray_SIZE(high_arr)
250-
flat = high_arr.flat
251-
for i in range(cnt):
252-
closed_upper = int(flat[i]) - 1
253-
if closed_upper > {{ub}}:
254-
raise ValueError('high is out of bounds for {{nptype}}')
255-
if closed_upper < {{lb}}:
256-
raise ValueError('low >= high')
257-
highm1_data[i] = <{{nptype}}_t>closed_upper
258-
259-
if np.any(np.greater(low_arr, highm1_arr)):
260-
raise ValueError('low >= high')
261-
262-
high_arr = highm1_arr
263-
low_arr = <np.ndarray>np.PyArray_FROM_OTF(low, np.{{npctype}}, np.NPY_ALIGNED | np.NPY_FORCECAST)
264-
265-
if size is not None:
266-
out_arr = <np.ndarray>np.empty(size, np.{{nptype}})
267-
else:
268-
it = np.PyArray_MultiIterNew2(low_arr, high_arr)
269-
out_arr = <np.ndarray>np.empty(it.shape, np.{{nptype}})
270-
271-
it = np.PyArray_MultiIterNew3(low_arr, high_arr, out_arr)
272-
out_data = <uint64_t *>np.PyArray_DATA(out_arr)
273-
n = np.PyArray_SIZE(out_arr)
274-
mask = last_rng = 0
275-
with lock, nogil:
276-
for i in range(n):
277-
low_v = (<{{nptype}}_t*>np.PyArray_MultiIter_DATA(it, 0))[0]
278-
high_v = (<{{nptype}}_t*>np.PyArray_MultiIter_DATA(it, 1))[0]
279-
rng = <{{utype}}_t>(high_v - low_v) # No -1 here since implemented above
280-
off = <{{utype}}_t>(<{{nptype}}_t>low_v)
281-
282-
if rng != last_rng:
283-
mask = _gen_mask(rng)
284-
out_data[i] = random_bounded_uint64(state, off, rng, mask)
285-
286-
np.PyArray_MultiIter_NEXT(it)
287-
288-
return out_arr
248+
return _rand_{{nptype}}_broadcast(low_arr, high_arr, size, state, lock)
289249
{{endfor}}

randomstate/tests/test_numpy_mt19937_regressions.py

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1,14 +1,14 @@
11
from __future__ import division, absolute_import, print_function
22

33
import sys
4-
from numpy.testing import (TestCase, run_module_suite, assert_,
4+
from numpy.testing import (run_module_suite, assert_,
55
assert_array_equal, assert_raises)
66
from numpy.compat import long
77
import numpy as np
88
from randomstate.prng.mt19937 import mt19937
99

1010

11-
class TestRegression(TestCase):
11+
class TestRegression(object):
1212

1313
def test_VonMises_range(self):
1414
# Make sure generated random variables are in [-pi, pi].

0 commit comments

Comments
 (0)