diff --git a/CHANGELOG.md b/CHANGELOG.md index 97fd31e5a50..eb3f0df9799 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -89,6 +89,7 @@ This release is compatible with NumPy 2.5. * Fixed `dpnp.ndarray.view` ignoring the USM element offset of a sliced array, which also caused `dpnp.einsum` to silently return wrong results for a single sliced operand with no summed index [#3037](https://github.com/IntelPython/dpnp/pull/3037) * Fixed `dpnp.all` and `dpnp.any` aborting when reducing over an empty axis (e.g. an array with a zero-length dimension) [#3021](https://github.com/IntelPython/dpnp/pull/3021) * Released the GIL before the blocking OneMKL DFT calls in the FFT extension [#3040](https://github.com/IntelPython/dpnp/pull/3040) +* Fixed `dpnp.linspace` returning `nan` for equal infinite endpoints [#3043](https://github.com/IntelPython/dpnp/pull/3043) ### Security diff --git a/dpnp/dpnp_algo/dpnp_arraycreation.py b/dpnp/dpnp_algo/dpnp_arraycreation.py index ff83209f4d2..faaecfcfea2 100644 --- a/dpnp/dpnp_algo/dpnp_arraycreation.py +++ b/dpnp/dpnp_algo/dpnp_arraycreation.py @@ -191,23 +191,37 @@ def dpnp_linspace( step_num = (num - 1) if endpoint else num if dpnp.isscalar(start) and dpnp.isscalar(stop): - # Call linspace() function for scalars. - usm_res = dpt.linspace( - start, - stop, - num, - dtype=dt, - usm_type=_usm_type, - sycl_queue=sycl_queue_normalized, - endpoint=endpoint, - ) + if start == stop: + # equal endpoints => constant array + zero step + usm_res = dpt.full( + num, + start, + dtype=dt, + usm_type=_usm_type, + sycl_queue=sycl_queue_normalized, + ) - # calculate the used step to return - if retstep is True: - if step_num > 0: - step = (stop - start) / step_num - else: - step = dpnp.nan + # calculate the used step to return + if retstep is True: + step = dt.type(0) if step_num > 0 else dpnp.nan + else: + # Call linspace() function for scalars. + usm_res = dpt.linspace( + start, + stop, + num, + dtype=dt, + usm_type=_usm_type, + sycl_queue=sycl_queue_normalized, + endpoint=endpoint, + ) + + # calculate the used step to return + if retstep is True: + if step_num > 0: + step = (stop - start) / step_num + else: + step = dpnp.nan else: usm_start = dpt.asarray( start, @@ -219,7 +233,9 @@ def dpnp_linspace( stop, dtype=dt, usm_type=_usm_type, sycl_queue=sycl_queue_normalized ) - delta = usm_stop - usm_start + # zero the delta where endpoints coincide, else `inf - inf = NaN` + # propagates (NaN != NaN untouched) + delta = dpt.where((usm_stop == usm_start), 0, (usm_stop - usm_start)) usm_res = dpt.arange( 0, diff --git a/dpnp/dpnp_iface_arraycreation.py b/dpnp/dpnp_iface_arraycreation.py index 98d9945be33..2090fe848d1 100644 --- a/dpnp/dpnp_iface_arraycreation.py +++ b/dpnp/dpnp_iface_arraycreation.py @@ -2864,7 +2864,7 @@ def linspace( There are `num` equally spaced samples in the closed interval [`start`, `stop`] or the half-open interval [`start`, `stop`) (depending on whether `endpoint` is ``True`` or ``False``). - step : float, optional + step : dpnp.ndarray, optional Only returned if `retstep` is ``True``. Size of spacing between samples. diff --git a/dpnp/tensor/libtensor/include/kernels/constructors.hpp b/dpnp/tensor/libtensor/include/kernels/constructors.hpp index 6ea42d17a9e..ac03b15004d 100644 --- a/dpnp/tensor/libtensor/include/kernels/constructors.hpp +++ b/dpnp/tensor/libtensor/include/kernels/constructors.hpp @@ -152,6 +152,13 @@ class LinearSequenceAffineFunctor void operator()(sycl::id<1> wiid) const { auto i = wiid.get(0); + + // equal endpoints => constant, and avoids `inf * 0 = NaN` + if (start_v == end_v) { + p[i] = start_v; + return; + } + wTy wc = wTy(i) / n; wTy w = wTy(n - i) / n; using dpnp::tensor::type_utils::is_complex; diff --git a/dpnp/tests/tensor/test_usm_ndarray_ctor.py b/dpnp/tests/tensor/test_usm_ndarray_ctor.py index 8a791524546..d76c15ea9b4 100644 --- a/dpnp/tests/tensor/test_usm_ndarray_ctor.py +++ b/dpnp/tests/tensor/test_usm_ndarray_ctor.py @@ -1573,6 +1573,18 @@ def test_linspace_int(): assert np.array_equal(dpt.asnumpy(X), Xnp) +@pytest.mark.parametrize("dtype", ["f2", "f4", "f8", "c8", "c16"]) +@pytest.mark.parametrize("endpoint", [True, False]) +def test_linspace_inf_equal_endpoints(dtype, endpoint): + q = get_queue_or_skip() + skip_if_dtype_not_supported(dtype, q) + val = complex(np.inf, np.inf) if dpt.dtype(dtype).kind == "c" else np.inf + X = dpt.linspace( + val, val, num=5, endpoint=endpoint, dtype=dtype, sycl_queue=q + ) + assert np.array_equal(dpt.asnumpy(X), np.full(5, val, dtype=dtype)) + + @pytest.mark.parametrize( "dt", _all_dtypes, diff --git a/dpnp/tests/test_arraycreation.py b/dpnp/tests/test_arraycreation.py index 1dd8e905577..f0776e4f1ad 100644 --- a/dpnp/tests/test_arraycreation.py +++ b/dpnp/tests/test_arraycreation.py @@ -21,6 +21,7 @@ get_array, get_float_dtypes, has_support_aspect64, + numpy_version, ) from .third_party.cupy import testing @@ -250,6 +251,59 @@ def test_axis(self, axis): func = lambda xp: xp.linspace([2, 3], [20, 15], num=10, axis=axis) assert_allclose(func(dpnp), func(numpy)) + @pytest.mark.parametrize("val", [numpy.inf, -numpy.inf, numpy.inf + 1j]) + @pytest.mark.parametrize("num", [1, 5]) + @pytest.mark.parametrize("endpoint", [True, False]) + def test_inf_equal_endpoints_scalar(self, val, num, endpoint): + result, step = dpnp.linspace( + val, val, num, endpoint=endpoint, retstep=True + ) + if numpy_version() >= "2.6.0": + expected, exp_step = numpy.linspace( + val, val, num, endpoint=endpoint, retstep=True + ) + assert_dtype_allclose(step, exp_step) + else: + expected = numpy.full(num, val) + step_val = step.asnumpy() + if (num - endpoint) > 0: + assert step_val == 0 + else: + assert numpy.isnan(step_val) + assert_dtype_allclose(result, expected) + + def test_inf_equal_endpoints_array(self): + start = numpy.array([numpy.inf, -numpy.inf, 1.0]) + stop = numpy.array([numpy.inf, -numpy.inf, 1.0]) + + result = dpnp.linspace(start, stop, num=4) + if numpy_version() >= "2.6.0": + expected = numpy.linspace(start, stop, num=4) + else: + expected = numpy.full((4, 3), [numpy.inf, -numpy.inf, 1.0]) + assert_dtype_allclose(result, expected) + + def test_inf_mixed_endpoints_array(self): + start = numpy.array([numpy.inf, numpy.inf]) + stop = numpy.array([numpy.inf, 2.0]) + + result = dpnp.linspace(start, stop, num=3) + if numpy_version() >= "2.6.0": + expected = numpy.linspace(start, stop, num=3) + assert_dtype_allclose(result, expected) + else: + # mixed infinities still yield NaN interior; equal column stays inf + res = result.asnumpy() + assert res[0, 0] == numpy.inf and res[-1, 0] == numpy.inf + assert numpy.isnan(res[1, 1]) + assert res[-1, 1] == 2.0 + + @pytest.mark.parametrize("num", [1, 5]) + def test_nan_endpoints(self, num): + result = dpnp.linspace(numpy.nan, numpy.nan, num) + expected = numpy.linspace(numpy.nan, numpy.nan, num) + assert_dtype_allclose(result, expected) + @pytest.mark.parametrize("xp", [dpnp, numpy]) def test_negative_num(self, xp): with pytest.raises(ValueError, match="must be non-negative"):