Fix dpnp.linspace returning NaN for equal infinite endpoints

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
 

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,

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.
 

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;

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,

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"):