From 4bbb5f298130a4c3a6b2ba5deebe3716888c6ad3 Mon Sep 17 00:00:00 2001 From: Vahid Tavanashad Date: Wed, 13 Mar 2024 20:12:52 -0500 Subject: [PATCH 1/4] update dpnp.cond --- dpnp/linalg/dpnp_algo_linalg.pyx | 25 ----- dpnp/linalg/dpnp_iface_linalg.py | 57 +++++++++--- dpnp/linalg/dpnp_utils_linalg.py | 46 ++++++++- tests/skipped_tests.tbl | 13 --- tests/skipped_tests_gpu.tbl | 13 --- tests/test_linalg.py | 155 ++++++++++++++++++++++++++----- tests/test_sycl_queue.py | 24 +++++ tests/test_usm_type.py | 14 +++ 8 files changed, 258 insertions(+), 89 deletions(-) diff --git a/dpnp/linalg/dpnp_algo_linalg.pyx b/dpnp/linalg/dpnp_algo_linalg.pyx index 126e7bdbea50..367247c23685 100644 --- a/dpnp/linalg/dpnp_algo_linalg.pyx +++ b/dpnp/linalg/dpnp_algo_linalg.pyx @@ -45,7 +45,6 @@ cimport numpy cimport dpnp.dpnp_utils as utils __all__ = [ - "dpnp_cond", "dpnp_eig", "dpnp_eigvals", ] @@ -60,30 +59,6 @@ ctypedef c_dpctl.DPCTLSyclEventRef(*custom_linalg_2in_1out_func_ptr_t)(c_dpctl.D const c_dpctl.DPCTLEventVectorRef) -cpdef object dpnp_cond(object input, object p): - if p in ('f', 'fro'): - # TODO: change order='K' when support is implemented - input = dpnp.ravel(input, order='C') - sqnorm = dpnp.dot(input, input) - res = dpnp.sqrt(sqnorm) - ret = dpnp.array([res]) - elif p == dpnp.inf: - dpnp_sum_val = dpnp.sum(dpnp.abs(input), axis=1) - ret = dpnp.max(dpnp_sum_val) - elif p == -dpnp.inf: - dpnp_sum_val = dpnp.sum(dpnp.abs(input), axis=1) - ret = dpnp.min(dpnp_sum_val) - elif p == 1: - dpnp_sum_val = dpnp.sum(dpnp.abs(input), axis=0) - ret = dpnp.max(dpnp_sum_val) - elif p == -1: - dpnp_sum_val = dpnp.sum(dpnp.abs(input), axis=0) - ret = dpnp.min(dpnp_sum_val) - else: - ret = dpnp.array([input.item(0)]) - return ret - - cpdef tuple dpnp_eig(utils.dpnp_descriptor x1): cdef shape_type_c x1_shape = x1.shape diff --git a/dpnp/linalg/dpnp_iface_linalg.py b/dpnp/linalg/dpnp_iface_linalg.py index 07b2d46e0f7c..87931661d07d 100644 --- a/dpnp/linalg/dpnp_iface_linalg.py +++ b/dpnp/linalg/dpnp_iface_linalg.py @@ -48,6 +48,7 @@ check_stacked_2d, check_stacked_square, dpnp_cholesky, + dpnp_cond, dpnp_det, dpnp_eigh, dpnp_inv, @@ -144,32 +145,60 @@ def cholesky(a, upper=False): return dpnp_cholesky(a, upper=upper) -def cond(input, p=None): +def cond(x, p=None): """ Compute the condition number of a matrix. For full documentation refer to :obj:`numpy.linalg.cond`. - Limitations - ----------- - Input array is supported as :obj:`dpnp.ndarray`. - Parameter p=[None, 1, -1, 2, -2, dpnp.inf, -dpnp.inf, 'fro'] is supported. + Parameters + ---------- + x : {dpnp.ndarray, usm_ndarray} + The matrix whose condition number is sought. + p : {None, 1, -1, 2, -2, inf, -inf, "fro"}, optional + Order of the norm used in the condition number computation. + inf means dpnp's `inf` object. The default is ``None``. + + Returns + ------- + out : dpnp.ndarray + The condition number of the matrix. May be infinite. See Also -------- :obj:`dpnp.norm` : Matrix or vector norm. - """ - if not use_origin_backend(input): - if p in [None, 1, -1, 2, -2, dpnp.inf, -dpnp.inf, "fro"]: - result_obj = dpnp_cond(input, p) - result = dpnp.convert_single_elem_array_to_scalar(result_obj) + Examples + -------- + >>> import dpnp as np + >>> a = np.array([[1, 0, -1], [0, 1, 0], [1, 0, 1]]) + >>> a + array([[ 1, 0, -1], + [ 0, 1, 0], + [ 1, 0, 1]]) + >>> np.linalg.cond(a) + array(1.4142135623730951) + >>> np.linalg.cond(a, 'fro') + array(3.1622776601683795) + >>> np.linalg.cond(a, np.inf) + array(2.) + >>> np.linalg.cond(a, -np.inf) + array(1.) + >>> np.linalg.cond(a, 1) + array(2.) + >>> np.linalg.cond(a, -1) + array(1.) + >>> np.linalg.cond(a, 2) + array(1.4142135623730951) + >>> np.linalg.cond(a, -2) + array(0.70710678118654746) # may vary + >>> min(np.linalg.svd(a, compute_uv=False))*min(np.linalg.svd(np.linalg.inv(a), compute_uv=False)) + array(0.70710678118654746) # may vary - return result - else: - pass + """ - return call_origin(numpy.linalg.cond, input, p) + dpnp.check_supported_arrays_type(x) + return dpnp_cond(x, p) def det(a): diff --git a/dpnp/linalg/dpnp_utils_linalg.py b/dpnp/linalg/dpnp_utils_linalg.py index 3a1b565f02eb..ef1d27b908aa 100644 --- a/dpnp/linalg/dpnp_utils_linalg.py +++ b/dpnp/linalg/dpnp_utils_linalg.py @@ -38,6 +38,7 @@ "check_stacked_2d", "check_stacked_square", "dpnp_cholesky", + "dpnp_cond", "dpnp_det", "dpnp_eigh", "dpnp_inv", @@ -199,6 +200,11 @@ def _common_inexact_type(default_dtype, *dtypes): return dpnp.result_type(*inexact_dtypes) +def _is_empty_2d(arr): + # check size first for efficiency + return arr.size == 0 and prod(arr.shape[-2:]) == 0 + + def _lu_factor(a, res_type): """ Compute pivoted LU decomposition. @@ -841,6 +847,40 @@ def dpnp_cholesky(a, upper): return a_h +def dpnp_cond(x, p=None): + """Compute the condition number of a matrix.""" + + if _is_empty_2d(x): + raise dpnp.linalg.LinAlgError("cond is not defined on empty arrays") + if p is None or p == 2 or p == -2: + s = dpnp.linalg.svd(x, compute_uv=False) + with numpy.errstate(all="ignore"): + if p == -2: + r = s[..., -1] / s[..., 0] + else: + r = s[..., 0] / s[..., -1] + else: + # Call inv(x) ignoring errors. The result array will + # contain nans in the entries where inversion failed. + check_stacked_2d(x) + check_stacked_square(x) + result_t = _common_type(x) + with numpy.errstate(all="ignore"): + invx = dpnp.linalg.inv(x) + r = dpnp.linalg.norm(x, p, axis=(-2, -1)) * dpnp.linalg.norm( + invx, p, axis=(-2, -1) + ) + r = r.astype(result_t, copy=False) + + # Convert nans to infs unless the original array had nan entries + nan_mask = dpnp.isnan(r) + if nan_mask.any(): + nan_mask &= ~dpnp.isnan(x).any(axis=(-2, -1)) + r[nan_mask] = dpnp.inf + + return r + + def dpnp_det(a): """ dpnp_det(a) @@ -1222,18 +1262,18 @@ def dpnp_multi_dot(n, arrays, out=None): """Compute the dot product of two or more arrays in a single function call.""" if not arrays[0].ndim in [1, 2]: - raise numpy.linalg.LinAlgError( + raise dpnp.linalg.LinAlgError( f"{arrays[0].ndim}-dimensional array given. First array must be 1-D or 2-D." ) if not arrays[-1].ndim in [1, 2]: - raise numpy.linalg.LinAlgError( + raise dpnp.linalg.LinAlgError( f"{arrays[-1].ndim}-dimensional array given. Last array must be 1-D or 2-D." ) for arr in arrays[1:-1]: if arr.ndim != 2: - raise numpy.linalg.LinAlgError( + raise dpnp.linalg.LinAlgError( f"{arr.ndim}-dimensional array given. Inner arrays must be 2-D." ) diff --git a/tests/skipped_tests.tbl b/tests/skipped_tests.tbl index e19be4d2751f..4a2f8c22a73e 100644 --- a/tests/skipped_tests.tbl +++ b/tests/skipped_tests.tbl @@ -37,19 +37,6 @@ tests/third_party/cupy/fft_tests/test_fft.py::TestFftn_param_23_{axes=None, norm tests/third_party/intel/test_zero_copy_test1.py::test_dpnp_interaction_with_dpctl_memory -tests/test_linalg.py::test_cond[-1-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] -tests/test_linalg.py::test_cond[1-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] -tests/test_linalg.py::test_cond[-2-[[1, 0, -1], [0, 1, 0], [1, 0, 1]]] -tests/test_linalg.py::test_cond[2-[[1, 0, -1], [0, 1, 0], [1, 0, 1]]] -tests/test_linalg.py::test_cond[-2-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] -tests/test_linalg.py::test_cond[2-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] -tests/test_linalg.py::test_cond["fro"-[[1, 0, -1], [0, 1, 0], [1, 0, 1]]] -tests/test_linalg.py::test_cond["fro"-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] -tests/test_linalg.py::test_cond[None-[[1, 0, -1], [0, 1, 0], [1, 0, 1]]] -tests/test_linalg.py::test_cond[None-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] -tests/test_linalg.py::test_cond[-numpy.inf-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] -tests/test_linalg.py::test_cond[numpy.inf-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] - tests/test_linalg.py::test_matrix_rank[None-[0, 1]-float64] tests/test_linalg.py::test_matrix_rank[None-[0, 1]-float32] tests/test_linalg.py::test_matrix_rank[None-[0, 1]-int64] diff --git a/tests/skipped_tests_gpu.tbl b/tests/skipped_tests_gpu.tbl index e2b95476da0a..08ad452a137a 100644 --- a/tests/skipped_tests_gpu.tbl +++ b/tests/skipped_tests_gpu.tbl @@ -158,19 +158,6 @@ tests/third_party/cupy/random_tests/test_distributions.py::TestDistributionsMult tests/third_party/intel/test_zero_copy_test1.py::test_dpnp_interaction_with_dpctl_memory -tests/test_linalg.py::test_cond[-1-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] -tests/test_linalg.py::test_cond[1-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] -tests/test_linalg.py::test_cond[-2-[[1, 0, -1], [0, 1, 0], [1, 0, 1]]] -tests/test_linalg.py::test_cond[2-[[1, 0, -1], [0, 1, 0], [1, 0, 1]]] -tests/test_linalg.py::test_cond[-2-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] -tests/test_linalg.py::test_cond[2-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] -tests/test_linalg.py::test_cond["fro"-[[1, 0, -1], [0, 1, 0], [1, 0, 1]]] -tests/test_linalg.py::test_cond["fro"-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] -tests/test_linalg.py::test_cond[None-[[1, 0, -1], [0, 1, 0], [1, 0, 1]]] -tests/test_linalg.py::test_cond[None-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] -tests/test_linalg.py::test_cond[-numpy.inf-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] -tests/test_linalg.py::test_cond[numpy.inf-[[1, 2, 3], [4, 5, 6], [7, 8, 9]]] - tests/test_linalg.py::test_matrix_rank[None-[0, 1]-float64] tests/test_linalg.py::test_matrix_rank[None-[0, 1]-float32] tests/test_linalg.py::test_matrix_rank[None-[0, 1]-int64] diff --git a/tests/test_linalg.py b/tests/test_linalg.py index 6b9800a8d9b3..5e4109c77060 100644 --- a/tests/test_linalg.py +++ b/tests/test_linalg.py @@ -209,25 +209,139 @@ def test_cholesky_errors(self): assert_raises(inp.linalg.LinAlgError, inp.linalg.cholesky, a_dp) -@pytest.mark.parametrize( - "arr", - [[[1, 0, -1], [0, 1, 0], [1, 0, 1]], [[1, 2, 3], [4, 5, 6], [7, 8, 9]]], - ids=[ - "[[1, 0, -1], [0, 1, 0], [1, 0, 1]]", - "[[1, 2, 3], [4, 5, 6], [7, 8, 9]]", - ], -) -@pytest.mark.parametrize( - "p", - [None, 1, -1, 2, -2, numpy.inf, -numpy.inf, "fro"], - ids=["None", "1", "-1", "2", "-2", "numpy.inf", "-numpy.inf", '"fro"'], -) -def test_cond(arr, p): - a = numpy.array(arr) - ia = inp.array(a) - result = inp.linalg.cond(ia, p) - expected = numpy.linalg.cond(a, p) - assert_array_equal(expected, result) +class TestCond: + def setup_method(self): + numpy.random.seed(42) + + @pytest.mark.parametrize( + "shape", [(0, 4, 4), (4, 0, 3, 3)], ids=["(0, 5, 3)", "(4, 0, 2, 3)"] + ) + @pytest.mark.parametrize( + "p", + [None, -inp.Inf, -2, -1, 1, 2, inp.Inf, "fro"], + ids=["None", "-dpnp.Inf", "-2", "-1", "1", "2", "dpnp.Inf", "fro"], + ) + def test_cond_empty(self, shape, p): + a = numpy.empty(shape) + ia = inp.array(a) + + result = inp.linalg.cond(ia, p=p) + expected = numpy.linalg.cond(a, p=p) + assert_dtype_allclose(result, expected) + + @pytest.mark.parametrize( + "dtype", get_all_dtypes(no_bool=True, no_complex=True) + ) + @pytest.mark.parametrize( + "shape", [(4, 4), (2, 4, 3, 3)], ids=["(4, 4)", "(2, 4, 3, 3)"] + ) + @pytest.mark.parametrize( + "p", + [None, -inp.Inf, -2, -1, 1, 2, inp.Inf, "fro"], + ids=["None", "-dpnp.Inf", "-2", "-1", "1", "2", "dpnp.Inf", "fro"], + ) + def test_cond(self, dtype, shape, p): + a = numpy.array( + numpy.random.uniform(-5, 5, numpy.prod(shape)), dtype=dtype + ).reshape(shape) + ia = inp.array(a) + + result = inp.linalg.cond(ia, p=p) + expected = numpy.linalg.cond(a, p=p) + assert_dtype_allclose(result, expected) + + @pytest.mark.parametrize( + "p", + [None, -inp.Inf, -2, -1, 1, 2, inp.Inf, "fro"], + ids=["None", "-dpnp.Inf", "-2", "-1", "1", "2", "dpnp.Inf", "fro"], + ) + def test_cond_bool(self, p): + a = numpy.array([[True, True], [True, False]]) + ia = inp.array(a) + + result = inp.linalg.cond(ia, p=p) + expected = numpy.linalg.cond(a, p=p) + assert_dtype_allclose(result, expected) + + @pytest.mark.parametrize("dtype", get_complex_dtypes()) + @pytest.mark.parametrize( + "shape", [(4, 4), (2, 4, 3, 3)], ids=["(4, 4)", "(2, 4, 3, 3)"] + ) + @pytest.mark.parametrize( + "p", + [None, -inp.Inf, -2, -1, 1, 2, inp.Inf, "fro"], + ids=["None", "-dpnp.Inf", "-2", "-1", "1", "2", "dpnp.Inf", "fro"], + ) + def test_cond_complex(self, dtype, shape, p): + x1 = numpy.random.uniform(-5, 5, numpy.prod(shape)) + x2 = numpy.random.uniform(-5, 5, numpy.prod(shape)) + a = numpy.array(x1 + 1j * x2, dtype=dtype).reshape(shape) + ia = inp.array(a) + + result = inp.linalg.cond(ia, p=p) + expected = numpy.linalg.cond(a, p=p) + assert_dtype_allclose(result, expected) + + @pytest.mark.parametrize( + "p", + [-inp.Inf, -1, 1, inp.Inf, "fro"], + ids=["-dpnp.Inf", "-1", "1", "dpnp.Inf", "fro"], + ) + def test_cond_nan_input(self, p): + a = numpy.array(numpy.random.uniform(-10, 10, 9)).reshape(3, 3) + a[1, 1] = numpy.nan + ia = inp.array(a) + + result = inp.linalg.cond(ia, p=p) + expected = numpy.linalg.cond(a, p=p) + assert_dtype_allclose(result, expected) + + @pytest.mark.parametrize( + "p", + [None, -inp.Inf, -2, -1, 1, 2, inp.Inf, "fro"], + ids=["None", "-dpnp.Inf", "-2", "-1", "1", "2", "dpnp.Inf", "fro"], + ) + def test_cond_nan(self, p): + a = numpy.array(numpy.random.uniform(-5, 5, 16)).reshape(2, 2, 2, 2) + a[0, 0] = 0 + a[1, 1] = 0 + ia = inp.array(a) + + result = inp.linalg.cond(ia, p=p) + expected = numpy.linalg.cond(a, p=p) + assert_dtype_allclose(result, expected) + + @pytest.mark.parametrize( + "p", + [None, -inp.Inf, -2, -1, 1, 2, inp.Inf, "fro"], + ids=["None", "-dpnp.Inf", "-2", "-1", "1", "2", "dpnp.Inf", "fro"], + ) + @pytest.mark.parametrize( + "stride", + [(-2, -3, 2, -2), (-2, 4, -4, -4), (2, 3, 4, 4), (-1, 3, 3, -3)], + ids=[ + "(-2, -3, 2, -2)", + "(-2, 4, -4, -4)", + "(2, 3, 4, 4)", + "(-1, 3, 3, -3)", + ], + ) + def test_cond_strided(self, p, stride): + A = numpy.random.rand(6, 8, 10, 10) + B = inp.asarray(A) + slices = tuple(slice(None, None, stride[i]) for i in range(A.ndim)) + a = A[slices] + b = B[slices] + + result = inp.linalg.cond(b, p=p) + expected = numpy.linalg.cond(a, p=p) + assert_dtype_allclose(result, expected) + + def test_cond_error(self): + # cond is not defined on empty arrays + ia = inp.empty((2, 0)) + with pytest.raises(ValueError): + inp.linalg.cond(ia, p=1) class TestDet: @@ -1110,8 +1224,7 @@ def test_norm_strided_ND(self, axis, stride): assert_dtype_allclose(result, expected) def test_norm_error(self): - a = numpy.arange(120).reshape(2, 3, 4, 5) - ia = inp.array(a) + ia = inp.arange(120).reshape(2, 3, 4, 5) # Duplicate axes given with pytest.raises(ValueError): diff --git a/tests/test_sycl_queue.py b/tests/test_sycl_queue.py index 1043d4444e93..66891ce417f4 100644 --- a/tests/test_sycl_queue.py +++ b/tests/test_sycl_queue.py @@ -1142,6 +1142,30 @@ def test_cholesky(data, is_empty, device): assert_sycl_queue_equal(result_queue, expected_queue) +@pytest.mark.parametrize( + "device", + valid_devices, + ids=[device.filter_string for device in valid_devices], +) +@pytest.mark.parametrize( + "p", + [None, -dpnp.Inf, -2, -1, 1, 2, dpnp.Inf, "fro"], + ids=["None", "-dpnp.Inf", "-2", "-1", "1", "2", "dpnp.Inf", "fro"], +) +def test_cond(device, p): + numpy.random.seed(42) + a = numpy.array(numpy.random.uniform(-5, 5, 16)).reshape(4, 4) + ia = dpnp.array(a, device=device) + + result = dpnp.linalg.cond(ia, p=p) + expected = numpy.linalg.cond(a, p=p) + assert_dtype_allclose(result, expected) + + expected_queue = ia.get_array().sycl_queue + result_queue = result.get_array().sycl_queue + assert_sycl_queue_equal(result_queue, expected_queue) + + @pytest.mark.parametrize( "device", valid_devices, diff --git a/tests/test_usm_type.py b/tests/test_usm_type.py index 2eccb30ba283..0373ab547257 100644 --- a/tests/test_usm_type.py +++ b/tests/test_usm_type.py @@ -680,6 +680,20 @@ def test_concat_stack(func, data1, data2, usm_type_x, usm_type_y): assert z.usm_type == du.get_coerced_usm_type([usm_type_x, usm_type_y]) +@pytest.mark.parametrize("usm_type", list_of_usm_types, ids=list_of_usm_types) +@pytest.mark.parametrize( + "p", + [None, -dp.Inf, -2, -1, 1, 2, dp.Inf, "fro"], + ids=["None", "-dpnp.Inf", "-2", "-1", "1", "2", "dpnp.Inf", "fro"], +) +def test_cond(usm_type, p): + ia = dp.arange(32, usm_type=usm_type).reshape(2, 4, 4) + + result = dp.linalg.cond(ia, p=p) + assert ia.usm_type == usm_type + assert result.usm_type == usm_type + + @pytest.mark.parametrize("usm_type", list_of_usm_types, ids=list_of_usm_types) def test_multi_dot(usm_type): numpy_array_list = [] From 85563356ca4ebeea924643849100305fbb8db1d3 Mon Sep 17 00:00:00 2001 From: Vahid Tavanashad Date: Wed, 3 Apr 2024 08:33:30 -0500 Subject: [PATCH 2/4] address comments --- dpnp/linalg/dpnp_iface_linalg.py | 36 +++++++++++-------- dpnp/linalg/dpnp_utils_linalg.py | 26 ++++++-------- tests/test_linalg.py | 2 +- tests/test_product.py | 12 +++---- .../cupy/linalg_tests/test_norms.py | 3 -- 5 files changed, 39 insertions(+), 40 deletions(-) diff --git a/dpnp/linalg/dpnp_iface_linalg.py b/dpnp/linalg/dpnp_iface_linalg.py index 87931661d07d..60cb1fd6e07a 100644 --- a/dpnp/linalg/dpnp_iface_linalg.py +++ b/dpnp/linalg/dpnp_iface_linalg.py @@ -157,7 +157,8 @@ def cond(x, p=None): The matrix whose condition number is sought. p : {None, 1, -1, 2, -2, inf, -inf, "fro"}, optional Order of the norm used in the condition number computation. - inf means dpnp's `inf` object. The default is ``None``. + inf means the `dpnp.inf` object, and the Frobenius norm is + the root-of-sum-of-squares norm. The default is ``None``. Returns ------- @@ -166,7 +167,7 @@ def cond(x, p=None): See Also -------- - :obj:`dpnp.norm` : Matrix or vector norm. + :obj:`dpnp.linalg.norm` : Matrix or vector norm. Examples -------- @@ -177,9 +178,9 @@ def cond(x, p=None): [ 0, 1, 0], [ 1, 0, 1]]) >>> np.linalg.cond(a) - array(1.4142135623730951) + array(1.41421356) >>> np.linalg.cond(a, 'fro') - array(3.1622776601683795) + array(3.16227766) >>> np.linalg.cond(a, np.inf) array(2.) >>> np.linalg.cond(a, -np.inf) @@ -189,11 +190,11 @@ def cond(x, p=None): >>> np.linalg.cond(a, -1) array(1.) >>> np.linalg.cond(a, 2) - array(1.4142135623730951) + array(1.41421356) >>> np.linalg.cond(a, -2) - array(0.70710678118654746) # may vary + array(0.70710678) # may vary >>> min(np.linalg.svd(a, compute_uv=False))*min(np.linalg.svd(np.linalg.inv(a), compute_uv=False)) - array(0.70710678118654746) # may vary + array(0.70710678) # may vary """ @@ -375,6 +376,11 @@ def inv(a): out : (..., M, M) dpnp.ndarray (Multiplicative) inverse of the matrix a. + See Also + -------- + :obj:`dpnp.linalg.cond` : Compute the condition number of a matrix. + :obj:`dpnp.linalg.svd` : Compute the singular value decomposition. + Examples -------- >>> import dpnp as np @@ -676,11 +682,11 @@ def norm(x, ord=None, axis=None, keepdims=False): [ 2, 3, 4]]) >>> np.linalg.norm(a) - array(7.745966692414834) + array(7.74596669) >>> np.linalg.norm(b) - array(7.745966692414834) + array(7.74596669) >>> np.linalg.norm(b, 'fro') - array(7.745966692414834) + array(7.74596669) >>> np.linalg.norm(a, np.inf) array(4.) >>> np.linalg.norm(b, np.inf) @@ -699,16 +705,16 @@ def norm(x, ord=None, axis=None, keepdims=False): >>> np.linalg.norm(b, -1) array(6.) >>> np.linalg.norm(a, 2) - array(7.745966692414834) + array(7.74596669) >>> np.linalg.norm(b, 2) - array(7.3484692283495345) + array(7.34846923) >>> np.linalg.norm(a, -2) array(0.) >>> np.linalg.norm(b, -2) - array(1.8570331885190563e-016) # may vary + array(4.35106603e-18) # may vary >>> np.linalg.norm(a, 3) - array(5.8480354764257312) # may vary + array(5.84803548) # may vary >>> np.linalg.norm(a, -3) array(0.) @@ -729,7 +735,7 @@ def norm(x, ord=None, axis=None, keepdims=False): >>> np.linalg.norm(m, axis=(1,2)) array([ 3.74165739, 11.22497216]) >>> np.linalg.norm(m[0, :, :]), np.linalg.norm(m[1, :, :]) - (array(3.7416573867739413), array(11.224972160321824)) + (array(3.74165739), array(11.22497216)) """ diff --git a/dpnp/linalg/dpnp_utils_linalg.py b/dpnp/linalg/dpnp_utils_linalg.py index ef1d27b908aa..df35f05d4de3 100644 --- a/dpnp/linalg/dpnp_utils_linalg.py +++ b/dpnp/linalg/dpnp_utils_linalg.py @@ -202,7 +202,7 @@ def _common_inexact_type(default_dtype, *dtypes): def _is_empty_2d(arr): # check size first for efficiency - return arr.size == 0 and prod(arr.shape[-2:]) == 0 + return arr.size == 0 and numpy.prod(arr.shape[-2:]) == 0 def _lu_factor(a, res_type): @@ -854,22 +854,18 @@ def dpnp_cond(x, p=None): raise dpnp.linalg.LinAlgError("cond is not defined on empty arrays") if p is None or p == 2 or p == -2: s = dpnp.linalg.svd(x, compute_uv=False) - with numpy.errstate(all="ignore"): - if p == -2: - r = s[..., -1] / s[..., 0] - else: - r = s[..., 0] / s[..., -1] + if p == -2: + r = s[..., -1] / s[..., 0] + else: + r = s[..., 0] / s[..., -1] else: - # Call inv(x) ignoring errors. The result array will - # contain nans in the entries where inversion failed. - check_stacked_2d(x) - check_stacked_square(x) result_t = _common_type(x) - with numpy.errstate(all="ignore"): - invx = dpnp.linalg.inv(x) - r = dpnp.linalg.norm(x, p, axis=(-2, -1)) * dpnp.linalg.norm( - invx, p, axis=(-2, -1) - ) + # The result array will contain nans in the entries + # where inversion failed + invx = dpnp.linalg.inv(x) + r = dpnp.linalg.norm(x, p, axis=(-2, -1)) * dpnp.linalg.norm( + invx, p, axis=(-2, -1) + ) r = r.astype(result_t, copy=False) # Convert nans to infs unless the original array had nan entries diff --git a/tests/test_linalg.py b/tests/test_linalg.py index 5e4109c77060..4e1f4cdb04d7 100644 --- a/tests/test_linalg.py +++ b/tests/test_linalg.py @@ -211,7 +211,7 @@ def test_cholesky_errors(self): class TestCond: def setup_method(self): - numpy.random.seed(42) + numpy.random.seed(70) @pytest.mark.parametrize( "shape", [(0, 4, 4), (4, 0, 3, 3)], ids=["(0, 5, 3)", "(4, 0, 2, 3)"] diff --git a/tests/test_product.py b/tests/test_product.py index 5661a33b74de..4b8ab53b0f28 100644 --- a/tests/test_product.py +++ b/tests/test_product.py @@ -929,23 +929,23 @@ def test_multi_dot_error(self): a = dpnp.ones((5, 8, 10)) b = dpnp.ones((10, 5)) - c = dpnp.ones((8, 15)) + c = dpnp.ones((5, 15)) # First array must be 1-D or 2-D - with pytest.raises(numpy.linalg.LinAlgError): + with pytest.raises(dpnp.linalg.LinAlgError): dpnp.linalg.multi_dot([a, b, c]) - a = dpnp.ones((5, 8)) + a = dpnp.ones((5, 10)) b = dpnp.ones((10, 5)) - c = dpnp.ones((8, 15, 6)) + c = dpnp.ones((5, 15, 6)) # Last array must be 1-D or 2-D - with pytest.raises(numpy.linalg.LinAlgError): + with pytest.raises(dpnp.linalg.LinAlgError): dpnp.linalg.multi_dot([a, b, c]) a = dpnp.ones((5, 10)) b = dpnp.ones((10, 5, 8)) c = dpnp.ones((8, 15)) # Inner array must be 2-D - with pytest.raises(numpy.linalg.LinAlgError): + with pytest.raises(dpnp.linalg.LinAlgError): dpnp.linalg.multi_dot([a, b, c]) a = dpnp.ones((5, 10)) diff --git a/tests/third_party/cupy/linalg_tests/test_norms.py b/tests/third_party/cupy/linalg_tests/test_norms.py index e866d1eb11ed..7af691fb8c27 100644 --- a/tests/third_party/cupy/linalg_tests/test_norms.py +++ b/tests/third_party/cupy/linalg_tests/test_norms.py @@ -87,9 +87,6 @@ def test_matrix_rank(self, xp, dtype): return y -# TODO: Remove the use of fixture for all tests in this file -# when dpnp.prod() will support complex dtypes on Gen9 -@pytest.mark.usefixtures("allow_fall_back_on_numpy") class TestDet(unittest.TestCase): @testing.for_dtypes("fdFD") @testing.numpy_cupy_allclose(rtol=1e-3, atol=1e-4) From e9f0643f37f04a86808b4fe606dc2944ebf866ab Mon Sep 17 00:00:00 2001 From: Vahid Tavanashad Date: Thu, 4 Apr 2024 10:32:39 -0500 Subject: [PATCH 3/4] update pinv function --- dpnp/linalg/dpnp_utils_linalg.py | 10 +++------- tests/test_linalg.py | 2 +- 2 files changed, 4 insertions(+), 8 deletions(-) diff --git a/dpnp/linalg/dpnp_utils_linalg.py b/dpnp/linalg/dpnp_utils_linalg.py index 4172eca0b098..0e9cc5154752 100644 --- a/dpnp/linalg/dpnp_utils_linalg.py +++ b/dpnp/linalg/dpnp_utils_linalg.py @@ -1437,13 +1437,9 @@ def dpnp_pinv(a, rcond=1e-15, hermitian=False): """ - if a.size == 0: + if _is_empty_2d(a): m, n = a.shape[-2:] - if m == 0 or n == 0: - res_type = a.dtype - else: - res_type = _common_type(a) - return dpnp.empty_like(a, shape=(a.shape[:-2] + (n, m)), dtype=res_type) + return dpnp.empty_like(a, shape=(a.shape[:-2] + (n, m)), dtype=a.dtype) if dpnp.is_supported_array_type(rcond): # Check that `a` and `rcond` are allocated on the same device @@ -1453,7 +1449,7 @@ def dpnp_pinv(a, rcond=1e-15, hermitian=False): # Allocate dpnp.ndarray if rcond is a scalar rcond = dpnp.array(rcond, usm_type=a.usm_type, sycl_queue=a.sycl_queue) - u, s, vt = dpnp_svd(a.conj(), full_matrices=False, hermitian=hermitian) + u, s, vt = dpnp_svd(dpnp.conj(a), full_matrices=False, hermitian=hermitian) # discard small singular values cutoff = rcond * dpnp.max(s, axis=-1) diff --git a/tests/test_linalg.py b/tests/test_linalg.py index 815fb388bb55..9a53b315763f 100644 --- a/tests/test_linalg.py +++ b/tests/test_linalg.py @@ -336,7 +336,7 @@ def test_cond_strided(self, p, stride): result = inp.linalg.cond(b, p=p) expected = numpy.linalg.cond(a, p=p) - assert_dtype_allclose(result, expected) + assert_dtype_allclose(result, expected, factor=24) def test_cond_error(self): # cond is not defined on empty arrays From 2844d6713be9ee568aadb4001af1befeb3791804 Mon Sep 17 00:00:00 2001 From: Anton <100830759+antonwolfy@users.noreply.github.com> Date: Sun, 7 Apr 2024 15:30:20 +0200 Subject: [PATCH 4/4] Update dpnp/linalg/dpnp_utils_linalg.py --- dpnp/linalg/dpnp_utils_linalg.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/dpnp/linalg/dpnp_utils_linalg.py b/dpnp/linalg/dpnp_utils_linalg.py index 0e9cc5154752..ed9d0b88472d 100644 --- a/dpnp/linalg/dpnp_utils_linalg.py +++ b/dpnp/linalg/dpnp_utils_linalg.py @@ -1439,7 +1439,7 @@ def dpnp_pinv(a, rcond=1e-15, hermitian=False): if _is_empty_2d(a): m, n = a.shape[-2:] - return dpnp.empty_like(a, shape=(a.shape[:-2] + (n, m)), dtype=a.dtype) + return dpnp.empty_like(a, shape=(a.shape[:-2] + (n, m))) if dpnp.is_supported_array_type(rcond): # Check that `a` and `rcond` are allocated on the same device