Merge pull request #282 from bashtage/edge-cases

MAINT: Port output check from NumPy

doc/source/change-log.rst

@@ -16,8 +16,10 @@ Change Log
   cannot update NumPy.
 
 
-Since v1.20.0
-=============
+v1.20.1
+=======
+- Fixed a bug that affects :func:`~randomgen.generator.Generator.standard_gamma` when
+  used with ``out`` and a Fortran contiguous array.
 - Added :func:`~randomgen.generator.ExtendedGenerator.multivariate_complex_normal`.
 - Added :func:`~randomgen.generator.ExtendedGenerator.standard_wishart` and
   :func:`~randomgen.generator.ExtendedGenerator.wishart` variate generators.

doc/source/names_wordlist.txt

@@ -115,4 +115,11 @@ Wooldridge
 Zhenyu
 ecuyer
 Overton
-Horner
+Horner
+Feiveson
+Jour
+Uhlig
+Dıaz
+Garcıa
+Jáimez
+Mardia

doc/source/spelling_wordlist.txt

@@ -279,3 +279,8 @@ Intrinsics
 precomputed
 args
 kwargs
+Wishart
+wishart
+Fortran
+loc
+trivariate

randomgen/common.pyx

@@ -533,23 +533,44 @@ cdef validate_output_shape(iter_shape, np.ndarray output):
         )
 
 
-cdef check_output(object out, object dtype, object size):
+cdef check_output(object out, object dtype, object size, bint require_c_array):
+    """
+    Check user-supplied output array properties and shape
+    
+    Parameters
+    ----------
+    out : {ndarray, None}
+        The array to check.  If None, returns immediately.
+    dtype : dtype
+        The required dtype of out.
+    size : {None, int, tuple[int]}
+        The size passed.  If out is an ndarray, verifies that the shape of out
+        matches size.
+    require_c_array : bool
+        Whether out must be a C-array.  If False, out can be either C- or F-
+        ordered.  If True, must be C-ordered. In either case, must be
+        contiguous, writable, aligned and in native byte-order.
+    """
     if out is None:
         return
     cdef np.ndarray out_array = <np.ndarray>out
-    if not (np.PyArray_CHKFLAGS(out_array, api.NPY_ARRAY_CARRAY) or
-            np.PyArray_CHKFLAGS(out_array, api.NPY_ARRAY_FARRAY)):
-        raise ValueError("Supplied output array is not contiguous, writable or aligned.")
+    if not (np.PyArray_ISCARRAY(out_array) or
+            (np.PyArray_ISFARRAY(out_array) and not require_c_array)):
+        req = "C-" if require_c_array else ""
+        raise ValueError(
+            f'Supplied output array must be {req}contiguous, writable, '
+            f'aligned, and in machine byte-order.'
+        )
     if out_array.dtype != dtype:
-        raise TypeError("Supplied output array has the wrong type. "
-                        "Expected {0}, got {0}".format(dtype, out_array.dtype))
+        raise TypeError('Supplied output array has the wrong type. '
+                        'Expected {0}, got {1}'.format(np.dtype(dtype), out_array.dtype))
     if size is not None:
         try:
             tup_size = tuple(size)
         except TypeError:
             tup_size = tuple([size])
         if tup_size != out.shape:
-            raise ValueError("size must match out.shape when used together")
+            raise ValueError('size must match out.shape when used together')
 
 
 cdef object double_fill(void *func, bitgen_t *state, object size, object lock, object out):
@@ -565,7 +586,7 @@ cdef object double_fill(void *func, bitgen_t *state, object size, object lock, o
             return out_val
 
     if out is not None:
-        check_output(out, np.float64, size)
+        check_output(out, np.float64, size, False)
         out_array = <np.ndarray>out
     else:
         out_array = <np.ndarray>np.empty(size, np.double)
@@ -587,7 +608,7 @@ cdef object float_fill(void *func, bitgen_t *state, object size, object lock, ob
             return random_func(state)
 
     if out is not None:
-        check_output(out, np.float32, size)
+        check_output(out, np.float32, size, False)
         out_array = <np.ndarray>out
     else:
         out_array = <np.ndarray>np.empty(size, np.float32)
@@ -610,7 +631,7 @@ cdef object float_fill_from_double(void *func, bitgen_t *state, object size, obj
             return <float>random_func(state)
 
     if out is not None:
-        check_output(out, np.float32, size)
+        check_output(out, np.float32, size, False)
         out_array = <np.ndarray>out
     else:
         out_array = <np.ndarray>np.empty(size, np.float32)
@@ -826,7 +847,7 @@ cdef object cont(void *func, void *state, object size, object lock, int narg,
     cdef np.ndarray a_arr, b_arr, c_arr
     cdef double _a = 0.0, _b = 0.0, _c = 0.0
     cdef bint is_scalar = True
-    check_output(out, np.float64, size)
+    check_output(out, np.float64, size, narg > 0)
     if narg > 0:
         a_arr = <np.ndarray>np.PyArray_FROM_OTF(a, np.NPY_DOUBLE, api.NPY_ARRAY_ALIGNED)
         is_scalar = is_scalar and np.PyArray_NDIM(a_arr) == 0
@@ -1276,7 +1297,7 @@ cdef object cont_f(void *func, bitgen_t *state, object size, object lock,
     cdef float _a
     cdef bint is_scalar = True
     cdef int requirements = api.NPY_ARRAY_ALIGNED | api.NPY_ARRAY_FORCECAST
-    check_output(out, np.float32, size)
+    check_output(out, np.float32, size, True)
     a_arr = <np.ndarray>np.PyArray_FROMANY(a, np.NPY_FLOAT32, 0, 0, requirements)
     is_scalar = np.PyArray_NDIM(a_arr) == 0
 

randomgen/generator.pyx

@@ -5355,7 +5355,7 @@ cdef class ExtendedGenerator:
 
         Notes
         -----
-        Uses the method of Odell and Fieveson [1]_ when `df` >= `dim`.
+        Uses the method of Odell and Feiveson [1]_ when `df` >= `dim`.
         Otherwise variates are directly generated as the inner product
         of `df` by `dim` arrays of standard normal random variates.
 
@@ -5405,12 +5405,12 @@ cdef class ExtendedGenerator:
         """
         wishart(df, scale, size=None, *, check_valid="warn", tol=None, rank=None, method="svd")
 
-        Draw samples from the Wishart and psuedo-Wishart distributions.
+        Draw samples from the Wishart and pseudo-Wishart distributions.
 
         Parameters
         ----------
         df : {int, array_like[int]}
-            Degree-of-freedom values. In array-like must boradcast with all
+            Degree-of-freedom values. In array-like must broadcast with all
             but the final two dimensions of ``shape``.
         scale : array_like
             Shape matrix of the distribution. It must be symmetric and
@@ -5459,7 +5459,7 @@ cdef class ExtendedGenerator:
 
         Notes
         -----
-        Uses the method of Odell and Fieveson [1]_ when `df` >= `dim`.
+        Uses the method of Odell and Feiveson [1]_ when `df` >= `dim`.
         Otherwise variates are directly generated as the inner product
         of `df` by `dim` arrays of standard normal random variates.
 
@@ -5491,10 +5491,10 @@ cdef class ExtendedGenerator:
         shape_arr = <np.ndarray>np.asarray(scale, dtype=np.float64, order="C")
         shape_nd = np.PyArray_NDIM(shape_arr)
         msg = (
-            "scale must have at least 2 dimensions. The final two dimensions "
-            "must be the same so that scale's shape is (...,N,N)."
+            "scale must be non-empty and have at least 2 dimensions. The final "
+            "two dimensions must be the same so that scale's shape is (...,N,N)."
         )
-        if shape_nd < 2:
+        if shape_nd < 2 or shape_arr.size == 0:
             raise ValueError(msg)
         dim = np.shape(shape_arr)[shape_nd-1]
         rank_val = dim if rank is None else int(rank)
@@ -5720,8 +5720,8 @@ and the trailing dimensions must match exactly so that
                                                np.NPY_ARRAY_ALIGNED |
                                                np.NPY_ARRAY_C_CONTIGUOUS)
         ldim = np.PyArray_NDIM(larr)
-        if ldim < 1:
-            raise ValueError("loc must be at least 1-dimensional")
+        if ldim < 1 or larr.size == 0:
+            raise ValueError("loc must be non-empty and at least 1-dimensional")
         dim = np.PyArray_DIMS(larr)[ldim - 1]
 
         if gamma is None:
@@ -5734,10 +5734,16 @@ and the trailing dimensions must match exactly so that
 
         gdim = np.PyArray_NDIM(garr)
         gshape = np.PyArray_DIMS(garr)
-        if gdim < 2 or gshape[gdim - 2] != gshape[gdim - 1] or gshape[gdim - 1] != dim:
+        if (
+                gdim < 2 or
+                gshape[gdim - 2] != gshape[gdim - 1] or
+                gshape[gdim - 1] != dim or
+                garr.size == 0
+        ):
             raise ValueError(
-                "gamma must be at least 2-dimensional and the final two dimensions "
-                f"must match the final dimension of loc, {dim}."
+                "gamma must be non-empty with at least 2-dimensional and the "
+                "final two dimensions must match the final dimension of loc,"
+                f" {dim}."
             )
         if relation is None:
             rarr = <np.ndarray>np.zeros((dim,dim), dtype=complex)
@@ -5748,10 +5754,16 @@ and the trailing dimensions must match exactly so that
                                                    np.NPY_ARRAY_C_CONTIGUOUS)
         rdim = np.PyArray_NDIM(rarr)
         rshape = np.PyArray_DIMS(rarr)
-        if rdim < 2 or rshape[rdim - 2] != rshape[rdim - 1] or rshape[rdim - 1] != dim:
+        if (
+                rdim < 2 or
+                rshape[rdim - 2] != rshape[rdim - 1] or
+                rshape[rdim - 1] != dim or
+                rarr.size == 0
+        ):
             raise ValueError(
-                "relation must be at least 2-dimensional and the final two dimensions "
-                f"must match the final dimension of loc, {dim}."
+                "relation must be non-empty with at least 2-dimensional and the "
+                "final two dimensions must match the final dimension of loc,"
+                f" {dim}."
             )
         can_bcast, cov_shape = broadcast_shape(np.shape(garr), np.shape(rarr), False)
         if not can_bcast:

randomgen/rdrand.pyx

@@ -178,7 +178,7 @@ cdef class RDRAND(BitGenerator):
     state structure, or use PyErr_Occurred to see if an error occurred
     during generation.
 
-    To see the exception you will generatr, you can run this invalid code
+    To see the exception you will generate, you can run this invalid code
 
     >>> from numpy.random import Generator
     >>> from randomgen import RDRAND

randomgen/tests/test_extended_generator.py

@@ -436,11 +436,11 @@ def test_missing_scipy_exception():
 
 def test_wishart_exceptions():
     eg = ExtendedGenerator()
-    with pytest.raises(ValueError, match="scale must have at"):
+    with pytest.raises(ValueError, match="scale must be non-empty"):
         eg.wishart(10, [10])
-    with pytest.raises(ValueError, match="scale must have at"):
+    with pytest.raises(ValueError, match="scale must be non-empty"):
         eg.wishart(10, 10)
-    with pytest.raises(ValueError, match="scale must have at"):
+    with pytest.raises(ValueError, match="scale must be non-empty"):
         eg.wishart(10, np.array([[1, 2]]))
     with pytest.raises(ValueError, match="At least one"):
         eg.wishart([], np.eye(2))
@@ -596,3 +596,23 @@ def test_mv_complex_normal_exceptions(extended_gen):
         extended_gen.multivariate_complex_normal(
             [0.0, 0.0], np.ones((4, 1, 3, 2, 2)), np.ones((1, 1, 2, 2, 2))
         )
+
+
+def test_wishart_edge(extended_gen):
+    with pytest.raises(ValueError, match="scale must be non-empty"):
+        extended_gen.wishart(5, np.empty((0, 0)))
+    with pytest.raises(ValueError, match="scale must be non-empty"):
+        extended_gen.wishart(5, np.empty((0, 2, 2)))
+    with pytest.raises(ValueError, match="scale must be non-empty"):
+        extended_gen.wishart(5, [[]])
+    with pytest.raises(ValueError, match="At least one value is required"):
+        extended_gen.wishart(np.empty((0, 2, 3)), np.eye(2))
+
+
+def test_mv_complex_normal_edge(extended_gen):
+    with pytest.raises(ValueError, match="loc must be non-empty and at least"):
+        extended_gen.multivariate_complex_normal(np.empty((0, 2)))
+    with pytest.raises(ValueError, match="gamma must be non-empty"):
+        extended_gen.multivariate_complex_normal([0, 0], np.empty((0, 2, 2)))
+    with pytest.raises(ValueError, match="relation must be non-empty"):
+        extended_gen.multivariate_complex_normal([0, 0], np.eye(2), np.empty((0, 2, 2)))