BUG: Ensure integer-type stream on 32bit

Ensure integer type is stream compatible on 32 bit Fix incorrect clause end Add integer-generator tests that check long streams
numpy · seberg · May 28, 2019 · Feb 20, 2018 · Feb 22, 2018 · Feb 23, 2018
commit b42a5ca0a076b40c612014dc540ca5f9bcf10f41
diff --git a/numpy/random/common.pxd b/numpy/random/common.pxd
@@ -96,6 +96,6 @@ cdef object cont_broadcast_3(void *func, void *state, object size, object lock,
                              np.ndarray c_arr, object c_name, constraint_type c_constraint)
 
 cdef object discrete_broadcast_iii(void *func, void *state, object size, object lock,
-                                  np.ndarray a_arr, object a_name, constraint_type a_constraint,
-                                  np.ndarray b_arr, object b_name, constraint_type b_constraint,
-                                  np.ndarray c_arr, object c_name, constraint_type c_constraint)
+                                   np.ndarray a_arr, object a_name, constraint_type a_constraint,
+                                   np.ndarray b_arr, object b_name, constraint_type b_constraint,
+                                   np.ndarray c_arr, object c_name, constraint_type c_constraint)
diff --git a/numpy/random/common.pyx b/numpy/random/common.pyx
@@ -666,9 +666,9 @@ cdef object discrete_broadcast_di(void *func, void *state, object size, object l
     return randoms
 
 cdef object discrete_broadcast_iii(void *func, void *state, object size, object lock,
-                                  np.ndarray a_arr, object a_name, constraint_type a_constraint,
-                                  np.ndarray b_arr, object b_name, constraint_type b_constraint,
-                                  np.ndarray c_arr, object c_name, constraint_type c_constraint):
+                                   np.ndarray a_arr, object a_name, constraint_type a_constraint,
+                                   np.ndarray b_arr, object b_name, constraint_type b_constraint,
+                                   np.ndarray c_arr, object c_name, constraint_type c_constraint):
     cdef np.ndarray randoms
     cdef int64_t *randoms_data
     cdef np.broadcast it

diff --git a/numpy/random/legacy_distributions.pxd b/numpy/random/legacy_distributions.pxd
@@ -5,7 +5,7 @@ from libc.stdint cimport int64_t
 import numpy as np
 cimport numpy as np
 
-from .distributions cimport bitgen_t
+from .distributions cimport bitgen_t, binomial_t
 
 cdef extern from "distributions-boxmuller.h":
 
@@ -35,6 +35,12 @@ cdef extern from "distributions-boxmuller.h":
     double legacy_wald(aug_bitgen_t *aug_state, double mean, double scale) nogil
     double legacy_lognormal(aug_bitgen_t *aug_state, double mean, double sigma) nogil
     int64_t legacy_negative_binomial(aug_bitgen_t *aug_state, double n, double p) nogil
+    int64_t legacy_random_hypergeometric(bitgen_t *bitgen_state, int64_t good, int64_t bad, int64_t sample) nogil
+    int64_t legacy_random_logseries(bitgen_t *bitgen_state, double p) nogil
+    int64_t legacy_random_poisson(bitgen_t *bitgen_state, double lam) nogil
+    int64_t legacy_random_zipf(bitgen_t *bitgen_state, double a) nogil
+    int64_t legacy_random_geometric(bitgen_t *bitgen_state, double p) nogil
+    void legacy_random_multinomial(bitgen_t *bitgen_state, long n, long *mnix, double *pix, np.npy_intp d, binomial_t *binomial) nogil
     double legacy_standard_cauchy(aug_bitgen_t *state) nogil
     double legacy_beta(aug_bitgen_t *aug_state, double a, double b) nogil
     double legacy_f(aug_bitgen_t *aug_state, double dfnum, double dfden) nogil

diff --git a/numpy/random/mtrand.pyx b/numpy/random/mtrand.pyx
@@ -21,6 +21,19 @@ from .legacy_distributions cimport *
 
 np.import_array()
 
+cdef object int64_to_long(object x):
+    """
+    Convert int64 to long for legacy compatibility, which used long for integer
+    distributions
+    """
+    cdef int64_t x64
+
+    if np.isscalar(x):
+        x64 = x
+        return <long>x64
+    return x.astype('l', casting='unsafe')
+
+
 cdef class RandomState:
     """
     RandomState(bit_generator=None)
@@ -3031,53 +3044,53 @@ cdef class RandomState:
 
         # Uses a custom implementation since self._binomial is required
         cdef double _dp = 0
-        cdef int64_t _in = 0
+        cdef long _in = 0
         cdef bint is_scalar = True
         cdef np.npy_intp i, cnt
         cdef np.ndarray randoms
-        cdef np.int64_t *randoms_data
+        cdef long *randoms_data
         cdef np.broadcast it
 
         p_arr = <np.ndarray>np.PyArray_FROM_OTF(p, np.NPY_DOUBLE, np.NPY_ALIGNED)
         is_scalar = is_scalar and np.PyArray_NDIM(p_arr) == 0
-        n_arr = <np.ndarray>np.PyArray_FROM_OTF(n, np.NPY_INT64, np.NPY_ALIGNED)
+        n_arr = <np.ndarray>np.PyArray_FROM_OTF(n, np.NPY_LONG, np.NPY_ALIGNED)
         is_scalar = is_scalar and np.PyArray_NDIM(n_arr) == 0
 
         if not is_scalar:
             check_array_constraint(p_arr, 'p', CONS_BOUNDED_0_1)
             check_array_constraint(n_arr, 'n', CONS_NON_NEGATIVE)
             if size is not None:
-                randoms = <np.ndarray>np.empty(size, np.int64)
+                randoms = <np.ndarray>np.empty(size, int)
             else:
                 it = np.PyArray_MultiIterNew2(p_arr, n_arr)
-                randoms = <np.ndarray>np.empty(it.shape, np.int64)
+                randoms = <np.ndarray>np.empty(it.shape, int)
 
-            randoms_data = <np.int64_t *>np.PyArray_DATA(randoms)
+            randoms_data = <long *>np.PyArray_DATA(randoms)
             cnt = np.PyArray_SIZE(randoms)
 
             it = np.PyArray_MultiIterNew3(randoms, p_arr, n_arr)
             with self.lock, nogil:
                 for i in range(cnt):
                     _dp = (<double*>np.PyArray_MultiIter_DATA(it, 1))[0]
-                    _in = (<int64_t*>np.PyArray_MultiIter_DATA(it, 2))[0]
-                    (<int64_t*>np.PyArray_MultiIter_DATA(it, 0))[0] = random_binomial(&self._bitgen, _dp, _in, &self._binomial)
+                    _in = (<long*>np.PyArray_MultiIter_DATA(it, 2))[0]
+                    (<long*>np.PyArray_MultiIter_DATA(it, 0))[0] = random_binomial(&self._bitgen, _dp, _in, &self._binomial)
 
                     np.PyArray_MultiIter_NEXT(it)
 
             return randoms
 
         _dp = PyFloat_AsDouble(p)
-        _in = <int64_t>n
+        _in = <long>n
         check_constraint(_dp, 'p', CONS_BOUNDED_0_1)
         check_constraint(<double>_in, 'n', CONS_NON_NEGATIVE)
 
         if size is None:
             with self.lock:
                 return random_binomial(&self._bitgen, _dp, _in, &self._binomial)
 
-        randoms = <np.ndarray>np.empty(size, np.int64)
+        randoms = <np.ndarray>np.empty(size, int)
         cnt = np.PyArray_SIZE(randoms)
-        randoms_data = <np.int64_t *>np.PyArray_DATA(randoms)
+        randoms_data = <long *>np.PyArray_DATA(randoms)
 
         with self.lock, nogil:
             for i in range(cnt):
@@ -3157,10 +3170,12 @@ cdef class RandomState:
         ...    print(i, "wells drilled, probability of one success =", probability)
 
         """
-        return disc(&legacy_negative_binomial, &self._aug_state, size, self.lock, 2, 0,
-                    n, 'n', CONS_POSITIVE,
-                    p, 'p', CONS_BOUNDED_0_1,
-                    0.0, '', CONS_NONE)
+        out = disc(&legacy_negative_binomial, &self._aug_state, size, self.lock, 2, 0,
+                   n, 'n', CONS_POSITIVE,
+                   p, 'p', CONS_BOUNDED_0_1,
+                   0.0, '', CONS_NONE)
+        # Match historical output type
+        return int64_to_long(out)
 
     def poisson(self, lam=1.0, size=None):
         """
@@ -3228,10 +3243,12 @@ cdef class RandomState:
         >>> s = np.random.poisson(lam=(100., 500.), size=(100, 2))
 
         """
-        return disc(&random_poisson, &self._bitgen, size, self.lock, 1, 0,
-                    lam, 'lam', CONS_POISSON,
-                    0.0, '', CONS_NONE,
-                    0.0, '', CONS_NONE)
+        out = disc(&legacy_random_poisson, &self._bitgen, size, self.lock, 1, 0,
+                   lam, 'lam', CONS_POISSON,
+                   0.0, '', CONS_NONE,
+                   0.0, '', CONS_NONE)
+        # Match historical output type
+        return int64_to_long(out)
 
     def zipf(self, a, size=None):
         """
@@ -3307,10 +3324,12 @@ cdef class RandomState:
         >>> plt.show()
 
         """
-        return disc(&random_zipf, &self._bitgen, size, self.lock, 1, 0,
-                    a, 'a', CONS_GT_1,
-                    0.0, '', CONS_NONE,
-                    0.0, '', CONS_NONE)
+        out = disc(&legacy_random_zipf, &self._bitgen, size, self.lock, 1, 0,
+                   a, 'a', CONS_GT_1,
+                   0.0, '', CONS_NONE,
+                   0.0, '', CONS_NONE)
+        # Match historical output type
+        return int64_to_long(out)
 
     def geometric(self, p, size=None):
         """
@@ -3358,10 +3377,12 @@ cdef class RandomState:
         0.34889999999999999 #random
 
         """
-        return disc(&random_geometric, &self._bitgen, size, self.lock, 1, 0,
-                    p, 'p', CONS_BOUNDED_GT_0_1,
-                    0.0, '', CONS_NONE,
-                    0.0, '', CONS_NONE)
+        out = disc(&legacy_random_geometric, &self._bitgen, size, self.lock, 1, 0,
+                   p, 'p', CONS_BOUNDED_GT_0_1,
+                   0.0, '', CONS_NONE,
+                   0.0, '', CONS_NONE)
+        # Match historical output type
+        return int64_to_long(out)
 
     def hypergeometric(self, ngood, nbad, nsample, size=None):
         """
@@ -3458,9 +3479,10 @@ cdef class RandomState:
         cdef np.ndarray ongood, onbad, onsample
         cdef int64_t lngood, lnbad, lnsample
 
-        ongood = <np.ndarray>np.PyArray_FROM_OTF(ngood, np.NPY_INT64, np.NPY_ALIGNED)
-        onbad = <np.ndarray>np.PyArray_FROM_OTF(nbad, np.NPY_INT64, np.NPY_ALIGNED)
-        onsample = <np.ndarray>np.PyArray_FROM_OTF(nsample, np.NPY_INT64, np.NPY_ALIGNED)
+        # This cast to long is required to ensure that the values are inbounds
+        ongood = <np.ndarray>np.PyArray_FROM_OTF(ngood, np.NPY_LONG, np.NPY_ALIGNED)
+        onbad = <np.ndarray>np.PyArray_FROM_OTF(nbad, np.NPY_LONG, np.NPY_ALIGNED)
+        onsample = <np.ndarray>np.PyArray_FROM_OTF(nsample, np.NPY_LONG, np.NPY_ALIGNED)
 
         if np.PyArray_NDIM(ongood) == np.PyArray_NDIM(onbad) == np.PyArray_NDIM(onsample) == 0:
 
@@ -3470,17 +3492,25 @@ cdef class RandomState:
 
             if lngood + lnbad < lnsample:
                 raise ValueError("ngood + nbad < nsample")
-            return disc(&random_hypergeometric, &self._bitgen, size, self.lock, 0, 3,
-                        lngood, 'ngood', CONS_NON_NEGATIVE,
-                        lnbad, 'nbad', CONS_NON_NEGATIVE,
-                        lnsample, 'nsample', CONS_GTE_1)
+            out = disc(&legacy_random_hypergeometric, &self._bitgen, size, self.lock, 0, 3,
+                       lngood, 'ngood', CONS_NON_NEGATIVE,
+                       lnbad, 'nbad', CONS_NON_NEGATIVE,
+                       lnsample, 'nsample', CONS_GTE_1)
+            # Match historical output type
+            return int64_to_long(out)
 
         if np.any(np.less(np.add(ongood, onbad), onsample)):
             raise ValueError("ngood + nbad < nsample")
-        return discrete_broadcast_iii(&random_hypergeometric, &self._bitgen, size, self.lock,
-                                      ongood, 'ngood', CONS_NON_NEGATIVE,
-                                      onbad, 'nbad', CONS_NON_NEGATIVE,
-                                      onsample, 'nsample', CONS_GTE_1)
+        # Convert to int64, if necessary, to use int64 infrastructure
+        ongood = ongood.astype(np.int64)
+        onbad = onbad.astype(np.int64)
+        onbad = onbad.astype(np.int64)
+        out = discrete_broadcast_iii(&legacy_random_hypergeometric,&self._bitgen, size, self.lock,
+                                     ongood, 'ngood', CONS_NON_NEGATIVE,
+                                     onbad, 'nbad', CONS_NON_NEGATIVE,
+                                     onsample, 'nsample', CONS_GTE_1)
+        # Match historical output type
+        return int64_to_long(out)
 
     def logseries(self, p, size=None):
         """
@@ -3557,10 +3587,12 @@ cdef class RandomState:
         >>> plt.show()
 
         """
-        return disc(&random_logseries, &self._bitgen, size, self.lock, 1, 0,
-                 p, 'p', CONS_BOUNDED_0_1,
-                 0.0, '', CONS_NONE,
-                 0.0, '', CONS_NONE)
+        out = disc(&legacy_random_logseries, &self._bitgen, size, self.lock, 1, 0,
+                   p, 'p', CONS_BOUNDED_0_1,
+                   0.0, '', CONS_NONE,
+                   0.0, '', CONS_NONE)
+        # Match historical output type
+        return int64_to_long(out)
 
     # Multivariate distributions:
     def multivariate_normal(self, mean, cov, size=None, check_valid='warn',
@@ -3808,8 +3840,8 @@ cdef class RandomState:
         cdef np.npy_intp d, i, sz, offset
         cdef np.ndarray parr, mnarr
         cdef double *pix
-        cdef int64_t *mnix
-        cdef int64_t ni
+        cdef long *mnix
+        cdef long ni
 
         d = len(pvals)
         parr = <np.ndarray>np.PyArray_FROM_OTF(pvals, np.NPY_DOUBLE, np.NPY_ALIGNED)
@@ -3826,16 +3858,16 @@ cdef class RandomState:
             except:
                 shape = tuple(size) + (d,)
 
-        multin = np.zeros(shape, dtype=np.int64)
+        multin = np.zeros(shape, dtype=int)
         mnarr = <np.ndarray>multin
-        mnix = <int64_t*>np.PyArray_DATA(mnarr)
+        mnix = <long*>np.PyArray_DATA(mnarr)
         sz = np.PyArray_SIZE(mnarr)
         ni = n
         check_constraint(ni, 'n', CONS_NON_NEGATIVE)
         offset = 0
         with self.lock, nogil:
             for i in range(sz // d):
-                random_multinomial(&self._bitgen, ni, &mnix[offset], pix, d, &self._binomial)
+                legacy_random_multinomial(&self._bitgen, ni, &mnix[offset], pix, d, &self._binomial)
                 offset += d
 
         return multin

diff --git a/numpy/random/setup.py b/numpy/random/setup.py
@@ -84,9 +84,8 @@ def generate_libraries(ext, build_dir):
         if DEBUG:
             EXTRA_LINK_ARGS += ['-debug']
             EXTRA_COMPILE_ARGS += ["-Zi", "/Od"]
-        if sys.version_info < (3, 0):
-            EXTRA_INCLUDE_DIRS += [join(MOD_DIR, 'src', 'common')]
 
+    LEGACY_DEFS = [('NP_RANDOM_LEGACY', '1')]
     DSFMT_DEFS = [('DSFMT_MEXP', '19937')]
     if USE_SSE2:
         if os.name == 'nt':
@@ -180,7 +179,7 @@ def generate_libraries(ext, build_dir):
                         extra_compile_args=EXTRA_COMPILE_ARGS,
                         extra_link_args=EXTRA_LINK_ARGS,
                         depends=['mtrand.pyx'],
-                        define_macros=defs + DSFMT_DEFS,
+                        define_macros=defs + DSFMT_DEFS + LEGACY_DEFS,
                         )
     return config