fix: validate the sampling method and reject an empty mean in multinormal_cholesky

CHANGELOG.md

@@ -12,8 +12,10 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 * Added a thread-safety section to the how-to guide for free-threaded Python [gh-159](https://github.com/IntelPython/mkl_random/pull/159)
 * Added an [ASV](https://asv.readthedocs.io/en/stable/) benchmark suite under `benchmarks/` and a `benchmark` optional dependency group [gh-184](https://github.com/IntelPython/mkl_random/pull/184)
 * Enabled support of Python 3.15 [gh-188](https://github.com/IntelPython/mkl_random/pull/188)
+* Exported the sampling method ids `ICDF`, `BOXMULLER`, `BOXMULLER2`, `POISNORM` and `PTPE` from `mkl_random`, as the `method` docstrings already documented [gh-189](https://github.com/IntelPython/mkl_random/pull/189)
 
 ### Changed
+* Unsupported `method` values (e.g. a misspelled name, a method the function does not implement, `None`, `bool`, `float`) now emit a `UserWarning` listing the supported methods before falling back to the function's default method, instead of falling back silently [gh-189](https://github.com/IntelPython/mkl_random/pull/189)
 * Sped up `normal`, `uniform`, `exponential`, `laplace`, `gumbel`, `logistic`, `rayleigh` and `lognormal` for array-valued parameters. Seeded results change for these array paths; scalar paths are unchanged [gh-171](https://github.com/IntelPython/mkl_random/pull/171)
 * Array parameters for these distributions must broadcast to the requested `size` without adding dimensions; previously accepted mismatches now raise `ValueError` [gh-171](https://github.com/IntelPython/mkl_random/pull/171)
 * `uniform` with array-valued bounds may return `high` due to floating-point rounding [gh-171](https://github.com/IntelPython/mkl_random/pull/171)
@@ -28,6 +30,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 * Reduced the memory use and sped up integer generation with `WH`, `MCG31`, `R250` and `MRG32K3A`, which lack `viRngUniformBits` support [gh-178](https://github.com/IntelPython/mkl_random/pull/178)
 
 ### Fixed
+* Fixed `multinormal_cholesky` and the other functions taking `method` routing integer-like ids such as `np.int64(0)`, `np.int64(2)` or `False` to `BoxMuller`; NumPy integers now select the method they name [gh-189](https://github.com/IntelPython/mkl_random/pull/189)
+* Fixed `multinormal_cholesky` raising `ZeroDivisionError` for an empty `mean`; it now raises `ValueError` [gh-189](https://github.com/IntelPython/mkl_random/pull/189)
 * Fixed an out-of-range integer `brng` indexing `brng_list` past its end, which returned uninitialized memory as random values; it now warns and falls back to `MT19937` [gh-177](https://github.com/IntelPython/mkl_random/pull/177)
 * Fixed `brng=0` being treated as unset, which left the state unseeded [gh-177](https://github.com/IntelPython/mkl_random/pull/177)
 * Fixed `MKLRandomState(seed, brng=None)` reading past `brng_list` and returning uninitialized memory; it now uses `MT19937` [gh-177](https://github.com/IntelPython/mkl_random/pull/177)

mkl_random/__init__.py

@@ -28,6 +28,11 @@
 from . import _init_helper
 from ._version import __version__
 from .mklrand import (
+    BOXMULLER,
+    BOXMULLER2,
+    ICDF,
+    POISNORM,
+    PTPE,
     MKLRandomState,
     RandomState,
     beta,
@@ -105,6 +110,11 @@
 __all__ = [
     "MKLRandomState",
     "RandomState",
+    "ICDF",
+    "BOXMULLER",
+    "BOXMULLER2",
+    "POISNORM",
+    "PTPE",
     "seed",
     "get_state",
     "set_state",

mkl_random/mklrand.pyx

@@ -1650,18 +1650,51 @@ _method_alias_dict_poisson = {
 }
 
 
+_method_names = {
+    ICDF: "ICDF",
+    BOXMULLER: "BoxMuller",
+    BOXMULLER2: "BoxMuller2",
+    POISNORM: "POISNORM",
+    PTPE: "PTPE",
+}
+
+
 def choose_method(method, mlist, alias_dict = None):
-    if (method not in mlist):
-        if (alias_dict is None) or (not isinstance(alias_dict, dict)):
-            # issue warning
-            return mlist[0]
-        else:
-            if method not in alias_dict.keys():
-                return mlist[0]
-            else:
-                return alias_dict[method]
-    else:
-        return method
+    """
+    Resolve ``method`` to one of the module-level method constants in
+    ``mlist``. Callers compare the result with ``is``, so the constant itself
+    is returned, never an equal-valued object such as ``np.int64(2)``.
+    Unsupported values, including ``bool`` and ``float``, which compare
+    equal to the integer constants, warn with the list of supported methods
+    and fall back to ``mlist[0]``.
+    """
+    for candidate in mlist:
+        # the constants themselves, including the default arguments
+        if candidate is method:
+            return candidate
+
+    index = None
+    if isinstance(method, str):
+        if isinstance(alias_dict, dict):
+            index = alias_dict.get(method)
+    elif not isinstance(method, (bool, np.bool_)):
+        try:
+            index = operator.index(method)
+        except TypeError:
+            pass
+    for candidate in mlist:
+        if candidate == index:
+            return candidate
+
+    names = [_method_names.get(m, m) for m in mlist]
+    supported = ", ".join(f"\"{name}\"" for name in names)
+    warnings.warn(
+        f"The sampling method {method!r} is not supported. "
+        f"Supported methods are {supported}; \"{names[0]}\" "
+        "will be used instead",
+        UserWarning
+    )
+    return mlist[0]
 
 
 _brng_dict = {
@@ -3197,7 +3230,7 @@ cdef class _MKLRandomState:
             Output shape.  If the given shape is, e.g., ``(m, n, k)``, then
             ``m * n * k`` samples are drawn.  Default is None, in which case a
             single value is returned.
-        method : 'ICDF, 'BoxMuller', 'BoxMuller2', optional
+        method : 'ICDF', 'BoxMuller', 'BoxMuller2', optional
             Sampling method used by Intel MKL. Can also be specified using
             tokens mkl_random.ICDF, mkl_random.BOXMULLER, mkl_random.BOXMULLER2
 
@@ -3270,7 +3303,7 @@ cdef class _MKLRandomState:
             Output shape.  If the given shape is, e.g., ``(m, n, k)``, then
             ``m * n * k`` samples are drawn.  Default is None, in which case a
             single value is returned.
-        method : 'ICDF, 'BoxMuller', 'BoxMuller2', optional
+        method : 'ICDF', 'BoxMuller', 'BoxMuller2', optional
             Sampling method used by Intel MKL. Can also be specified using
             tokens mkl_random.ICDF, mkl_random.BOXMULLER, mkl_random.BOXMULLER2
 
@@ -5258,7 +5291,7 @@ cdef class _MKLRandomState:
             Output shape.  If the given shape is, e.g., ``(m, n, k)``, then
             ``m * n * k`` samples are drawn.  Default is None, in which case a
             single value is returned.
-        method : 'ICDF, 'BoxMuller', optional
+        method : 'ICDF', 'BoxMuller', optional
             Sampling method used by Intel MKL. Can also be specified using
             tokens mkl_random.ICDF, mkl_random.BOXMULLER
 
@@ -5959,7 +5992,7 @@ cdef class _MKLRandomState:
             Output shape.  If the given shape is, e.g., ``(m, n, k)``, then
             ``m * n * k`` samples are drawn.  Default is None, in which case a
             single value is returned.
-        method : 'POISNORM, 'PTPE', optional
+        method : 'POISNORM', 'PTPE', optional
             Sampling method used by Intel MKL. Can also be specified using
             tokens mkl_random.POISNORM, mkl_random.PTPE
 
@@ -6693,9 +6726,6 @@ cdef class _MKLRandomState:
             Output shape.  If the given shape is, e.g., ``(m, n, k)``, then
             ``m * n * k`` samples are drawn.  Default is None, in which case a
             single value is returned.
-        method : 'ICDF, 'BoxMuller', 'BoxMuller2', optional
-            Sampling method used by Intel MKL. Can also be specified using
-            tokens mkl_random.ICDF, mkl_random.BOXMULLER, mkl_random.BOXMULLER2
 
         Returns
         -------
@@ -7308,7 +7338,7 @@ cdef class MKLRandomState(_MKLRandomState):
             generated, and packed in an `m`-by-`n`-by-`k` arrangement.  Because
             each sample is `N`-dimensional, the output shape is ``(m,n,k,N)``.
             If no shape is specified, a single (`N`-D) sample is returned.
-        method : 'ICDF, 'BoxMuller', 'BoxMuller2', optional
+        method : 'ICDF', 'BoxMuller', 'BoxMuller2', optional
             Sampling method used by Intel MKL. Can also be specified using
             tokens mkl_random.ICDF, mkl_random.BOXMULLER, mkl_random.BOXMULLER2
 
@@ -7409,6 +7439,8 @@ cdef class MKLRandomState(_MKLRandomState):
         if marr.ndim != 1:
             raise ValueError("mean must be 1 dimensional")
         dim = marr.shape[0]
+        if dim < 1:
+            raise ValueError("mean must have at least one element")
         if (tarr.ndim == 2):
             storage_mode = MATRIX
             if (tarr.shape[0] != tarr.shape[1]):

mkl_random/tests/test_random.py

@@ -1001,9 +1001,13 @@ def test_randomdist_lognormal(randomdist):
         ]
     )
     np.testing.assert_allclose(actual, desired, atol=1e-6, rtol=1e-10)
-    actual = rnd.lognormal(
-        mean=0.123456789, sigma=2.0, size=(3, 2), method="Box-Muller2"
-    )
+    # lognormal only supports ICDF and BoxMuller, so this falls back to ICDF
+    with pytest.warns(
+        UserWarning, match='Supported methods are "ICDF", "BoxMuller";'
+    ):
+        actual = rnd.lognormal(
+            mean=0.123456789, sigma=2.0, size=(3, 2), method="Box-Muller2"
+        )
     desired = np.array(
         [
             [0.2585388231094821, 0.43734953048924663],
@@ -1529,6 +1533,77 @@ def test_normal_family_array_methods(name, method, normal_method, size):
     )
 
 
+@pytest.mark.parametrize("ch", [np.zeros((0, 0)), np.zeros(0)])
+def test_multinormal_cholesky_empty_mean(ch):
+    with pytest.raises(ValueError, match="at least one element"):
+        rnd.MKLRandomState(123).multinormal_cholesky(np.array([]), ch, size=3)
+
+
+@pytest.mark.parametrize("method", ["ICDF", "BoxMuller", "BoxMuller2"])
+@pytest.mark.parametrize("to_int", [int, np.int64, np.uint8])
+def test_multinormal_cholesky_integer_method(method, to_int):
+    # Integer-like method ids must select the same method as the name.
+    ids = {"ICDF": 0, "BoxMuller": 1, "BoxMuller2": 2}
+    mean = np.array([0.1, -0.2])
+    chol_mat = np.array([[1.0, 0.0], [-0.5, 1.0]])
+    expected = rnd.MKLRandomState(123).multinormal_cholesky(
+        mean, chol_mat, size=8, method=method
+    )
+    with warnings.catch_warnings():
+        warnings.simplefilter("error")
+        actual = rnd.MKLRandomState(123).multinormal_cholesky(
+            mean, chol_mat, size=8, method=to_int(ids[method])
+        )
+    np.testing.assert_array_equal(actual, expected)
+
+
+@pytest.mark.parametrize(
+    "method", ["bogus", 7, -1, 1.0, None, False, True, [1]]
+)
+def test_unsupported_method_warns_and_uses_default(method):
+    # bool and float compare equal to the integer method ids, but are not
+    # valid method specifications.
+    mean = np.array([0.1, -0.2])
+    chol_mat = np.array([[1.0, 0.0], [-0.5, 1.0]])
+    expected = rnd.MKLRandomState(123).multinormal_cholesky(
+        mean, chol_mat, size=8, method="ICDF"
+    )
+    with pytest.warns(UserWarning, match="not supported"):
+        actual = rnd.MKLRandomState(123).multinormal_cholesky(
+            mean, chol_mat, size=8, method=method
+        )
+    np.testing.assert_array_equal(actual, expected)
+
+
+def test_unsupported_method_uses_each_functions_default():
+    expected = rnd.MKLRandomState(123).poisson(3.0, size=8, method="POISNORM")
+    with pytest.warns(UserWarning, match='"POISNORM" will be used instead'):
+        actual = rnd.MKLRandomState(123).poisson(3.0, size=8, method="bogus")
+    np.testing.assert_array_equal(actual, expected)
+
+
+@pytest.mark.parametrize(
+    "token,func,name",
+    [
+        ("ICDF", "standard_normal", "ICDF"),
+        ("BOXMULLER", "standard_normal", "BoxMuller"),
+        ("BOXMULLER2", "standard_normal", "BoxMuller2"),
+        ("POISNORM", "poisson", "POISNORM"),
+        ("PTPE", "poisson", "PTPE"),
+    ],
+)
+def test_method_tokens_are_exported(token, func, name):
+    # The docstrings refer to the method ids as mkl_random.ICDF etc.
+    assert token in rnd.__all__
+    expected = getattr(rnd.MKLRandomState(123), func)(size=8, method=name)
+    with warnings.catch_warnings():
+        warnings.simplefilter("error")
+        actual = getattr(rnd.MKLRandomState(123), func)(
+            size=8, method=getattr(rnd, token)
+        )
+    np.testing.assert_array_equal(actual, expected)
+
+
 @pytest.mark.parametrize(
     "name,draw,p",
     [