add explicit shape comparison

Author: Vahid Tavanashad

dpnp/dpnp_iface_linearalgebra.py

@@ -69,16 +69,31 @@
 
 
 # TODO: implement a specific scalar-array kernel
-def _call_multiply(a, b, out=None):
-    """Call multiply function for special cases of scalar-array dots."""
+def _call_multiply(a, b, out=None, outer_calc=False):
+    """
+    Call multiply function for special cases of scalar-array dots.
+
+    if `sc` is an scalar and `a` is an array of type float32, we have
+    dpnp.multiply(a, sc).dtype == dpnp.float32 and
+    numpy.multiply(a, sc).dtype == dpnp.float32.
+
+    However, for scalar-array dots such as dot function we have
+    dpnp.dot(a, sc).dtype == dpnp.float32 while
+    numpy.dot(a, sc).dtype == dpnp.float64.
+
+    We need to adjust the behavior of the multiply function when it is
+    being used for special cases of scalar-array dots.
+
+    """
 
     sc, arr = (a, b) if dpnp.isscalar(a) else (b, a)
     sc_dtype = map_dtype_to_device(type(sc), arr.sycl_device)
     res_dtype = dpnp.result_type(sc_dtype, arr)
+    multiply_func = dpnp.multiply.outer if outer_calc else dpnp.multiply
     if out is not None and out.dtype == arr.dtype:
-        res = dpnp.multiply(a, b, out=out)
+        res = multiply_func(a, b, out=out)
     else:
-        res = dpnp.multiply(a, b, dtype=res_dtype)
+        res = multiply_func(a, b, dtype=res_dtype)
     return dpnp.get_result_array(res, out, casting="no")
 
 
@@ -1109,16 +1124,15 @@ def outer(a, b, out=None):
 
     dpnp.check_supported_arrays_type(a, b, scalar_type=True, all_scalars=False)
     if dpnp.isscalar(a):
-        x1 = a
         x2 = dpnp.ravel(b)[None, :]
+        result = _call_multiply(a, x2, out=out, outer_calc=True)
     elif dpnp.isscalar(b):
         x1 = dpnp.ravel(a)[:, None]
-        x2 = b
+        result = _call_multiply(x1, b, out=out, outer_calc=True)
     else:
-        x1 = dpnp.ravel(a)
-        x2 = dpnp.ravel(b)
+        result = dpnp.multiply.outer(dpnp.ravel(a), dpnp.ravel(b), out=out)
 
-    return dpnp.multiply.outer(x1, x2, out=out)
+    return result
 
 
 def tensordot(a, b, axes=2):
@@ -1288,13 +1302,13 @@ def vdot(a, b):
         if b.size != 1:
             raise ValueError("The second array should be of size one.")
         a_conj = numpy.conj(a)
-        return _call_multiply(a_conj, b)
+        return dpnp.squeeze(_call_multiply(a_conj, b))
 
     if dpnp.isscalar(b):
         if a.size != 1:
             raise ValueError("The first array should be of size one.")
         a_conj = dpnp.conj(a)
-        return _call_multiply(a_conj, b)
+        return dpnp.squeeze(_call_multiply(a_conj, b))
 
     if a.ndim == 1 and b.ndim == 1:
         return dpnp_dot(a, b, out=None, conjugate=True)

dpnp/dpnp_utils/dpnp_utils_linearalgebra.py

@@ -1108,7 +1108,7 @@ def dpnp_multiplication(
                 result = dpnp.moveaxis(result, (-2, -1), axes_res)
             elif len(axes_res) == 1:
                 result = dpnp.moveaxis(result, (-1,), axes_res)
-            return dpnp.ascontiguousarray(result)
+            return result
 
         return dpnp.asarray(result, order=order)
 

dpnp/tests/helper.py

@@ -10,13 +10,28 @@
 from . import config
 
 
+def _assert_dtype(a_dt, b_dt, check_only_type_kind=False):
+    if check_only_type_kind:
+        assert a_dt.kind == b_dt.kind, f"{a_dt.kind} != {b_dt.kind}"
+    else:
+        assert a_dt == b_dt, f"{a_dt} != {b_dt}"
+
+
+def _assert_shape(a, b):
+    if hasattr(b, "shape"):
+        assert a.shape == b.shape, f"{a.shape} != {b.shape}"
+    else:
+        # numpy output is scalar, then dpnp is 0-D array
+        assert a.shape == (), f"{a.shape} != ()"
+
+
 def assert_dtype_allclose(
     dpnp_arr,
     numpy_arr,
     check_type=True,
     check_only_type_kind=False,
     factor=8,
-    relative_factor=None,
+    check_shape=True,
 ):
     """
     Assert DPNP and NumPy array based on maximum dtype resolution of input arrays
@@ -37,10 +52,13 @@ def assert_dtype_allclose(
     for all data types supported by DPNP when set to True.
     It is effective only when 'check_type' is also set to True.
     The parameter `factor` scales the resolution used for comparing the arrays.
+    The parameter `check_shape`, when True (default), asserts the shape of input arrays is the same.
 
     """
 
-    list_64bit_types = [numpy.float64, numpy.complex128]
+    if check_shape:
+        _assert_shape(dpnp_arr, numpy_arr)
+
     is_inexact = lambda x: hasattr(x, "dtype") and dpnp.issubdtype(
         x.dtype, dpnp.inexact
     )
@@ -57,34 +75,32 @@ def assert_dtype_allclose(
             else -dpnp.inf
         )
         tol = factor * max(tol_dpnp, tol_numpy)
-        assert_allclose(dpnp_arr.asnumpy(), numpy_arr, atol=tol, rtol=tol)
+        assert_allclose(dpnp_arr, numpy_arr, atol=tol, rtol=tol, strict=False)
         if check_type:
+            list_64bit_types = [numpy.float64, numpy.complex128]
             numpy_arr_dtype = numpy_arr.dtype
             dpnp_arr_dtype = dpnp_arr.dtype
             dpnp_arr_dev = dpnp_arr.sycl_device
 
             if check_only_type_kind:
-                assert dpnp_arr_dtype.kind == numpy_arr_dtype.kind
+                _assert_dtype(dpnp_arr_dtype, numpy_arr_dtype, True)
             else:
                 is_np_arr_f2 = numpy_arr_dtype == numpy.float16
 
                 if is_np_arr_f2:
                     if has_support_aspect16(dpnp_arr_dev):
-                        assert dpnp_arr_dtype == numpy_arr_dtype
+                        _assert_dtype(dpnp_arr_dtype, numpy_arr_dtype)
                 elif (
                     numpy_arr_dtype not in list_64bit_types
                     or has_support_aspect64(dpnp_arr_dev)
                 ):
-                    assert dpnp_arr_dtype == numpy_arr_dtype
+                    _assert_dtype(dpnp_arr_dtype, numpy_arr_dtype)
                 else:
-                    assert dpnp_arr_dtype.kind == numpy_arr_dtype.kind
+                    _assert_dtype(dpnp_arr_dtype, numpy_arr_dtype, True)
     else:
-        assert_array_equal(dpnp_arr.asnumpy(), numpy_arr)
+        assert_array_equal(dpnp_arr, numpy_arr, strict=False)
         if check_type and hasattr(numpy_arr, "dtype"):
-            if check_only_type_kind:
-                assert dpnp_arr.dtype.kind == numpy_arr.dtype.kind
-            else:
-                assert dpnp_arr.dtype == numpy_arr.dtype
+            _assert_dtype(dpnp_arr.dtype, numpy_arr.dtype, check_only_type_kind)
 
 
 def generate_random_numpy_array(

dpnp/tests/test_arraycreation.py

@@ -952,15 +952,13 @@ def test_ascontiguousarray1(data):
     result = dpnp.ascontiguousarray(data)
     expected = numpy.ascontiguousarray(data)
     assert_dtype_allclose(result, expected)
-    assert result.shape == expected.shape
 
 
 @pytest.mark.parametrize("data", [(), 1, (2, 3), [4]])
 def test_ascontiguousarray2(data):
     result = dpnp.ascontiguousarray(dpnp.array(data))
     expected = numpy.ascontiguousarray(numpy.array(data))
     assert_dtype_allclose(result, expected)
-    assert result.shape == expected.shape
 
 
 @pytest.mark.parametrize(
@@ -970,15 +968,13 @@ def test_asfortranarray1(data):
     result = dpnp.asfortranarray(data)
     expected = numpy.asfortranarray(data)
     assert_dtype_allclose(result, expected)
-    assert result.shape == expected.shape
 
 
 @pytest.mark.parametrize("data", [(), 1, (2, 3), [4]])
 def test_asfortranarray2(data):
     result = dpnp.asfortranarray(dpnp.array(data))
     expected = numpy.asfortranarray(numpy.array(data))
     assert_dtype_allclose(result, expected)
-    assert result.shape == expected.shape
 
 
 def test_meshgrid_raise_error():

dpnp/tests/test_arraypad.py

@@ -42,7 +42,6 @@ def test_basic(self, mode):
         result = dpnp.pad(a_dp, (25, 20), mode=mode)
         if mode == "empty":
             # omit uninitialized "empty" boundary from the comparison
-            assert result.shape == expected.shape
             assert_equal(result[25:-20], expected[25:-20])
         else:
             assert_array_equal(result, expected)
@@ -70,7 +69,6 @@ def test_non_contiguous_array(self, mode):
         result = dpnp.pad(a_dp, (2, 3), mode=mode)
         if mode == "empty":
             # omit uninitialized "empty" boundary from the comparison
-            assert result.shape == expected.shape
             assert_equal(result[2:-3, 2:-3], expected[2:-3, 2:-3])
         else:
             assert_array_equal(result, expected)
@@ -287,10 +285,10 @@ def test_linear_ramp_end_values(self):
         """Ensure that end values are exact."""
         a_dp = dpnp.ones(10).reshape(2, 5)
         a = dpnp.pad(a_dp, (223, 123), mode="linear_ramp")
-        assert_equal(a[:, 0], 0.0)
-        assert_equal(a[:, -1], 0.0)
-        assert_equal(a[0, :], 0.0)
-        assert_equal(a[-1, :], 0.0)
+        assert_equal(a[:, 0], 0.0, strict=False)
+        assert_equal(a[:, -1], 0.0, strict=False)
+        assert_equal(a[0, :], 0.0, strict=False)
+        assert_equal(a[-1, :], 0.0, strict=False)
 
     @pytest.mark.parametrize(
         "dtype", [numpy.uint32, numpy.uint64] + get_all_dtypes(no_none=True)
@@ -426,7 +424,6 @@ def test_empty(self):
         expected = numpy.pad(a_np, [(2, 3), (3, 1)], "empty")
         result = dpnp.pad(a_dp, [(2, 3), (3, 1)], "empty")
         # omit uninitialized "empty" boundary from the comparison
-        assert result.shape == expected.shape
         assert_equal(result[2:-3, 3:-1], expected[2:-3, 3:-1])
 
     # Check how padding behaves on arrays with an empty dimension.

dpnp/tests/test_fill.py

@@ -38,7 +38,7 @@ def test_fill_strided_array():
     expected = dpnp.tile(dpnp.asarray([0, 1], dtype=a.dtype), 50)
 
     b.fill(1)
-    assert_array_equal(b, 1)
+    assert_array_equal(b, 1, strict=False)
     assert_array_equal(a, expected)
 
 
@@ -51,7 +51,7 @@ def test_fill_strided_2d_array(order):
     expected[::-2, ::2] = 1
 
     b.fill(1)
-    assert_array_equal(b, 1)
+    assert_array_equal(b, 1, strict=False)
     assert_array_equal(a, expected)
 
 
@@ -60,27 +60,27 @@ def test_fill_memset(order):
     a = dpnp.ones((10, 10), dtype="i4", order=order)
     a.fill(0)
 
-    assert_array_equal(a, 0)
+    assert_array_equal(a, 0, strict=False)
 
 
 def test_fill_float_complex_to_int():
     a = dpnp.ones((10, 10), dtype="i4")
 
     a.fill(complex(2, 0))
-    assert_array_equal(a, 2)
+    assert_array_equal(a, 2, strict=False)
 
     a.fill(float(3))
-    assert_array_equal(a, 3)
+    assert_array_equal(a, 3, strict=False)
 
 
 def test_fill_complex_to_float():
     a = dpnp.ones((10, 10), dtype="f4")
 
     a.fill(complex(2, 0))
-    assert_array_equal(a, 2)
+    assert_array_equal(a, 2, strict=False)
 
 
 def test_fill_bool():
     a = dpnp.full(5, fill_value=7, dtype="i4")
     a.fill(True)
-    assert_array_equal(a, 1)
+    assert_array_equal(a, 1, strict=False)

dpnp/tests/test_histogram.py

@@ -733,7 +733,8 @@ def test_bins(self, bins):
         expected_hist, expected_edges = numpy.histogramdd(v, bins)
         result_hist, result_edges = dpnp.histogramdd(iv, bins_dpnp)
         assert_allclose(result_hist, expected_hist)
-        assert_allclose(result_edges, expected_edges)
+        for x, y in zip(result_edges, expected_edges):
+            assert_allclose(x, y)
 
     def test_no_side_effects(self):
         v = dpnp.array([[1.3, 2.5, 2.3]])
@@ -752,7 +753,8 @@ def test_01d(self, data):
         result_hist, result_edges = dpnp.histogramdd(ia)
 
         assert_allclose(result_hist, expected_hist)
-        assert_allclose(result_edges, expected_edges)
+        for x, y in zip(result_edges, expected_edges):
+            assert_allclose(x, y)
 
     def test_3d(self):
         a = dpnp.ones((10, 10, 10))
@@ -822,15 +824,21 @@ def test_nan_values(self):
         )
         result_hist, result_edges = dpnp.histogramdd(ione_nan, bins=[[0, 1]])
         assert_allclose(result_hist, expected_hist)
-        assert_allclose(result_edges, expected_edges)
+        # dpnp returns both result_hist and result_edges as float64 while
+        # numpy returns result_hist as float64 but result_edges as int64
+        for x, y in zip(result_edges, expected_edges):
+            assert_allclose(x, y, strict=False)
 
         # NaN is not counted
         expected_hist, expected_edges = numpy.histogramdd(
             all_nan, bins=[[0, 1]]
         )
         result_hist, result_edges = dpnp.histogramdd(iall_nan, bins=[[0, 1]])
         assert_allclose(result_hist, expected_hist)
-        assert_allclose(result_edges, expected_edges)
+        # dpnp returns both result_hist and result_edges as float64 while
+        # numpy returns result_hist as float64 but result_edges as int64
+        for x, y in zip(result_edges, expected_edges):
+            assert_allclose(x, y, strict=False)
 
     def test_bins_another_sycl_queue(self):
         v = dpnp.arange(7, 12, sycl_queue=dpctl.SyclQueue())
@@ -866,7 +874,10 @@ def test_different_bins_amount(self, bins_count):
         expected_hist, expected_edges = numpy.histogramdd(v, bins=[bins_count])
         result_hist, result_edges = dpnp.histogramdd(iv, bins=[bins_count])
         assert_array_equal(result_hist, expected_hist)
-        assert_allclose(result_edges, expected_edges)
+        # for float32 input, dpnp returns both result_hist and result_edges
+        # as float64 while numpy returns result_hist as float64 but
+        for x, y in zip(result_edges, expected_edges):
+            assert_allclose(x, y, strict=False)
 
 
 class TestHistogram2d:
@@ -1045,8 +1056,10 @@ def test_nan_values(self):
             ione_nan, ione_nan, bins=[[0, 1]] * 2
         )
         assert_allclose(result_hist, expected_hist)
-        assert_allclose(result_edges_x, expected_edges_x)
-        assert_allclose(result_edges_y, expected_edges_y)
+        # dpnp returns both result_hist and result_edges as float64 while
+        # numpy returns result_hist as float64 but result_edges as int64
+        assert_allclose(result_edges_x, expected_edges_x, strict=False)
+        assert_allclose(result_edges_y, expected_edges_y, strict=False)
 
         # NaN is not counted
         expected_hist, expected_edges_x, expected_edges_y = numpy.histogram2d(
@@ -1056,8 +1069,10 @@ def test_nan_values(self):
             iall_nan, iall_nan, bins=[[0, 1]] * 2
         )
         assert_allclose(result_hist, expected_hist)
-        assert_allclose(result_edges_x, expected_edges_x)
-        assert_allclose(result_edges_y, expected_edges_y)
+        # dpnp returns both result_hist and result_edges as float64 while
+        # numpy returns result_hist as float64 but result_edges as int64
+        assert_allclose(result_edges_x, expected_edges_x, strict=False)
+        assert_allclose(result_edges_y, expected_edges_y, strict=False)
 
     def test_bins_another_sycl_queue(self):
         x = y = dpnp.arange(7, 12, sycl_queue=dpctl.SyclQueue())
@@ -1107,5 +1122,11 @@ def test_different_bins_amount(self, bins_count):
             ix, iy, bins=bins_count
         )
         assert_array_equal(result_hist, expected_hist)
-        assert_allclose(result_edges_x, expected_edges_x, rtol=1e-6)
-        assert_allclose(result_edges_y, expected_edges_y, rtol=1e-6)
+        # dpnp returns both result_hist and result_edges as float64 while
+        # numpy returns result_hist as float64 but result_edges as float32
+        assert_allclose(
+            result_edges_x, expected_edges_x, rtol=1e-6, strict=False
+        )
+        assert_allclose(
+            result_edges_y, expected_edges_y, rtol=1e-6, strict=False
+        )

dpnp/tests/test_indexing.py

@@ -292,7 +292,7 @@ def test_indexing_array_negative_strides(self):
 
         slices = (slice(None), dpnp.array([0, 1, 2, 3]))
         arr[slices] = 10
-        assert_array_equal(arr, 10.0)
+        assert_equal(arr, 10.0, strict=False)
 
 
 class TestIx: