Improve performance of choices, add randrange step

Roughly double the performance of choices and choose_weighted, add a
step parameter to randrange to match the standard library.

.appveyor.yml

@@ -1,4 +1,4 @@
-version: 0.2.1-{build}
+version: 0.3.0-{build}
 pull_requests:
   do_not_increment_build_number: true
 environment:

README.rst

@@ -9,16 +9,18 @@ Performance Comparison
 
 All times are seconds to complete one million iterations (except for shuffle which is 100,000 iterations).
 The benchmarking script can be found `here
-<https://github.com/Noctem/cyrandom/blob/master/benchmark.py>`_.
+<https://github.com/Noctem/cyrandom/blob/master/test/benchmark.py>`_.
 
 =========== ======== ========= ========
 function    stdlib   cyrandom  speedup
 =========== ======== ========= ========
-randrange   1.6768   0.12038   13.93x
-randint     1.9949   0.11875   16.8x
-choice      1.0404   0.05610   18.55x
+randrange   1.7863   0.12805   13.95x
+randint     2.1198   0.12258   17.29x
+choice      1.0969   0.05815   18.86x
 shuffle     44.806   0.98338   45.56x
-choices     3.3368   0.96160   3.47x
+choices     3.3938   0.49235   6.89x
 uniform     0.28467  0.08858   3.21x
 triangular  0.72341  0.10253   7.06x
 =========== ======== ========= ========
+
+Be aware that (for performance reasons) there is less input validation in some of these functions than in the standard library, so ensure that your arguments are valid.

cyrandom/cyrandom.pxd

@@ -1,15 +1,17 @@
+from libc.stdint cimport int32_t, uint16_t, uint32_t, uint64_t
+
 cpdef double random() nogil
 
-cdef unsigned short bit_length(unsigned long n) nogil
+cdef uint16_t bit_length(uint32_t n) nogil
 
-cdef unsigned long getrandbits(unsigned short k) nogil
+cdef uint32_t getrandbits(uint16_t k) nogil
 
-cpdef long randrange(long start, long stop) nogil
+cpdef int32_t randrange(int32_t start, int32_t stop, uint32_t step=?) nogil
 
-cpdef long randint(long a, long b) nogil
+cpdef int32_t randint(int32_t a, int32_t b) nogil
 
 cpdef double uniform(double a, double b) nogil
 
 cpdef double triangular(double low=?, double high=?, double mode=?) nogil
 
-cpdef long triangular_int(long low, long high, long mode) nogil
+cpdef int32_t triangular_int(int32_t low, int32_t high, int32_t mode) nogil

cyrandom/cyrandom.pyx

@@ -2,9 +2,9 @@
 # cython: language_level=3
 
 from libc.math cimport sqrt
+from libc.stdint cimport int32_t, uint16_t, uint32_t
 
-from bisect import bisect as _bisect
-from itertools import accumulate as _accumulate
+from cpython.mem cimport PyMem_Malloc, PyMem_Free
 
 from ._seed cimport random_seed
 from ._mersenne cimport genrand_int32, genrand_res53
@@ -19,89 +19,142 @@ cpdef double random() nogil:
     return genrand_res53()
 
 
-cdef unsigned short bit_length(unsigned long n) nogil:
-    cdef unsigned short length = 0
+cdef uint16_t bit_length(uint32_t n) nogil:
+    cdef uint16_t length = 0
     while n != 0:
         length += 1
         n >>= 1
     return length
 
 
-cdef unsigned long getrandbits(unsigned short k) nogil:
+cdef uint32_t getrandbits(uint16_t k) nogil:
     return genrand_int32() >> (32 - k)
 
 
-cpdef long randrange(long start, long stop) nogil:
+cpdef int32_t randrange(int32_t start, int32_t stop, uint32_t step=1) nogil:
     """Choose a random item from range(start, stop).
 
     This fixes the problem with randint() which includes the
     endpoint; in Python this is usually not what you want.
     """
-    cdef unsigned long width = stop - start
-    return start + _randbelow(width)
+    cdef uint32_t width = stop - start
+    if step == 1:
+        return start + _randbelow(width)
 
+    cdef uint32_t n
 
-cpdef long randint(long a, long b) nogil:
+    if step > 0:
+        n = (width + step - 1) // step
+    else:
+        n = (width + step + 1) // step
+
+    return start + _randbelow(n)
+
+
+cpdef int32_t randint(int32_t a, int32_t b) nogil:
     """Return random integer in range [a, b], including both end points.
     """
-    return randrange(a, b+1)
+    cdef uint32_t width = b - a + 1
+    return a + _randbelow(width)
 
 
-cdef unsigned long _randbelow(unsigned long n) nogil:
+cdef uint32_t _randbelow(uint32_t n) nogil:
     """Return a random int in the range [0,n).  Raises ValueError if n==0.
     """
     _getrandbits = getrandbits
-    cdef unsigned short k = bit_length(n)  # don't use (n-1) here because n can be 1
-    cdef unsigned long r = _getrandbits(k)  # 0 <= r < 2**k
+    cdef uint16_t k = bit_length(n)  # don't use (n-1) here because n can be 1
+    cdef uint32_t r = _getrandbits(k)  # 0 <= r < 2**k
     while r >= n:
         r = _getrandbits(k)
     return r
 
 
 def choice(seq):
     """Choose a random element from a non-empty sequence."""
-    cdef unsigned int i = _randbelow(len(seq))
+    cdef uint32_t i = _randbelow(len(seq))
     return seq[i]
 
 
 def shuffle(list seq):
     """Shuffle list x in place, and return None.
     """
-    cdef unsigned int i, j
+    cdef uint32_t i, j
     randbelow = _randbelow
     for i in reversed(range(1, len(seq))):
         j = randbelow(i+1)
         seq[i], seq[j] = seq[j], seq[i]
 
 
-def choices(population, tuple weights=None, *, tuple cum_weights=None, unsigned int k=1):
+def choices(population, tuple weights=None, *, tuple cum_weights=None, Py_ssize_t k=1):
     """Return a k sized list of population elements chosen with replacement.
 
     If the relative weights or cumulative weights are not specified,
     the selections are made with equal probability.
 
     """
-    random = genrand_res53
-    cdef unsigned int i, total
+    cdef Py_ssize_t i, total
+    cdef uint16_t num_weights
+    cdef uint16_t* _cum_weights
+    cdef list output = []
+
     if cum_weights is None:
         if weights is None:
             total = len(population)
-            return [population[int(random() * total)] for i in range(k)]
-        cum_weights = tuple(_accumulate(weights))
-    total = cum_weights[-1]
-    bisect = _bisect
-    return [population[bisect(cum_weights, random() * total)] for i in range(k)]
+            for i in range(k):
+                output[i] = population[int(genrand_res53() * total)]
+            return output
+
+        num_weights = len(weights)
+        _cum_weights = <uint16_t*> PyMem_Malloc(num_weights * 2)
+        total = 0
+        for i in range(num_weights):
+            total += weights[i]
+            _cum_weights[i] = total
+    else:
+        num_weights = len(cum_weights)
+        _cum_weights = <uint16_t*> PyMem_Malloc(num_weights * 2)
+        for i in range(num_weights):
+            _cum_weights[i] = cum_weights[i]
+        total = _cum_weights[num_weights - 1]
+
+    cdef uint16_t lo, mid, hi
+    cdef double x
+    for i in range(k):
+        lo = 0
+        hi = num_weights
+        x = genrand_res53() * total
+        while lo < hi:
+            mid = (lo + hi) // 2
+            if x < _cum_weights[mid]:
+                hi = mid
+            else:
+                lo = mid + 1
+        output.append(population[lo])
+
+    PyMem_Free(_cum_weights)
+    return output
 
 
 def choose_weighted(tuple population, tuple cum_weights):
     """Return an item from population according to provided weights.
     """
-    if len(cum_weights) != len(population):
+    cdef uint16_t hi = len(population)
+
+    if len(cum_weights) != hi:
         raise ValueError('The number of weights does not match the population')
-    cdef unsigned int total = cum_weights[-1]
-    bisect = _bisect
-    random = genrand_res53
-    return population[bisect(cum_weights, random() * total)]
+
+    cdef uint32_t total = cum_weights[-1]
+    cdef uint16_t lo, mid
+    cdef double x
+    lo = 0
+    x = genrand_res53() * total
+    while lo < hi:
+        mid = (lo + hi) // 2
+        if x < cum_weights[mid]:
+            hi = mid
+        else:
+            lo = mid + 1
+    return population[lo]
 
 
 cpdef double uniform(double a, double b) nogil:
@@ -127,11 +180,11 @@ cpdef double triangular(double low=0.0, double high=1.0, double mode=0.5) nogil:
     return low + (high - low) * sqrt(u * c)
 
 
-cpdef long triangular_int(long low, long high, long mode) nogil:
+cpdef int32_t triangular_int(int32_t low, int32_t high, int32_t mode) nogil:
     cdef double c, u = genrand_res53()
     c = (mode - low) / (high - low)
     if u > c:
         u = 1.0 - u
         c = 1.0 - c
         low, high = high, low
-    return <long>(low + (high - low) * sqrt(u * c))
+    return <int32_t>(low + (high - low) * sqrt(u * c))

setup.py

@@ -17,7 +17,7 @@
 
 setup(
     name="cyrandom",
-    version='0.2.1',
+    version='0.3.0',
     description='Fast random number generation.',
     long_description="A fast cython replacement for the standard library's random module.",
     url='https://github.com/Noctem/cyrandom',