Skip to content

Conversation

vlad-perevezentsev
Copy link
Collaborator

This PR suggests adding dpnp.linalg.lu_solve() for 2D arrays similar to scipy.linalg.lu_solve
Support for ND inputs will be added in the next phase.

  • Have you provided a meaningful PR description?
  • Have you added a test, reproducer or referred to an issue with a reproducer?
  • Have you tested your changes locally for CPU and GPU devices?
  • Have you made sure that new changes do not introduce compiler warnings?
  • Have you checked performance impact of proposed changes?
  • Have you added documentation for your changes, if necessary?
  • Have you added your changes to the changelog?

@vlad-perevezentsev vlad-perevezentsev added this to the 0.19.0 release milestone Sep 4, 2025
@vlad-perevezentsev vlad-perevezentsev self-assigned this Sep 4, 2025
@antonwolfy antonwolfy changed the title Implement dpnp.linalg.solve() 2D inputs Implement dpnp.linalg.lu_solve() 2D inputs Sep 4, 2025
Copy link
Contributor

github-actions bot commented Sep 4, 2025

Array API standard conformance tests for dpnp=0.19.0dev3=py313h509198e_27 ran successfully.
Passed: 1226
Failed: 1
Skipped: 9

Copy link
Contributor

github-actions bot commented Sep 4, 2025

View rendered docs @ https://intelpython.github.io/dpnp/pull/2575/index.html

@@ -966,6 +968,75 @@ def lu_factor(a, overwrite_a=False, check_finite=True):
return dpnp_lu_factor(a, overwrite_a=overwrite_a, check_finite=check_finite)


def lu_solve(lu_and_piv, b, trans=0, overwrite_b=False, check_finite=True):
"""
Solve an equation system, a x = b, given the LU factorization of `a`
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
Solve an equation system, a x = b, given the LU factorization of `a`
Solve an equation system, a x = b, given the LU factorization of `a`.


Parameters
----------
(lu, piv) : {tuple of dpnp.ndarrays or usm_ndarrays}
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Now it is rendered in a wrong way in the documentation

Suggested change
(lu, piv) : {tuple of dpnp.ndarrays or usm_ndarrays}
lu, piv : {tuple of dpnp.ndarrays or usm_ndarrays}

Parameters
----------
(lu, piv) : {tuple of dpnp.ndarrays or usm_ndarrays}
LU factorization of matrix `a` ((M, N)) together with pivot indices.
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why do we need double brackets here?

Suggested change
LU factorization of matrix `a` ((M, N)) together with pivot indices.
LU factorization of matrix `a` (M, N) together with pivot indices.

A[i, i] = A.dtype.type(off + 1.0)
return A

@pytest.mark.parametrize("shape", [(1, 1), (2, 2), (3, 3), (5, 5)])
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I failed to find a test with (M, N) shapes, only with (M, M)

_manager = dpu.SequentialOrderManager[exec_q]
dep_evs = _manager.submitted_events

# oneMKL LAPACK getrf overwrites `a`.
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
# oneMKL LAPACK getrf overwrites `a`.
# oneMKL LAPACK getrf overwrites `lu`.

Comment on lines +2585 to +2592
ht_ev, dep_ev = ti._copy_usm_ndarray_into_usm_ndarray(
src=b_usm_arr,
dst=b_h.get_array(),
sycl_queue=b.sycl_queue,
depends=_manager.submitted_events,
)
_manager.add_event_pair(ht_ev, dep_ev)
dep_ev = [dep_ev]
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Isn't we need to override dep_evs with list of copy events?

Suggested change
ht_ev, dep_ev = ti._copy_usm_ndarray_into_usm_ndarray(
src=b_usm_arr,
dst=b_h.get_array(),
sycl_queue=b.sycl_queue,
depends=_manager.submitted_events,
)
_manager.add_event_pair(ht_ev, dep_ev)
dep_ev = [dep_ev]
ht_ev, b_copy_ev = ti._copy_usm_ndarray_into_usm_ndarray(
src=b_usm_arr,
dst=b_h.get_array(),
sycl_queue=b.sycl_queue,
depends=dep_evs,
)
_manager.add_event_pair(ht_ev, b_copy_ev)
dep_evs = [lu_copy_ev, b_copy_ev]

else:
# input is suitable for in-place modification
b_h = b
dep_ev = _manager.submitted_events
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

and then here:

Suggested change
dep_ev = _manager.submitted_events
dep_evs = [lu_copy_ev]

b_h = b
dep_ev = _manager.submitted_events

# oneMKL LAPACK getrf overwrites `a`.
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
# oneMKL LAPACK getrf overwrites `a`.
# oneMKL LAPACK getrf overwrites `piv`.

Comment on lines +2598 to +2611
# oneMKL LAPACK getrf overwrites `a`.
piv_h = dpnp.empty_like(piv, order="F", usm_type=res_usm_type)

# use DPCTL tensor function to fill the сopy of the pivot array
# from the pivot array
ht_ev, piv_copy_ev = ti._copy_usm_ndarray_into_usm_ndarray(
src=piv_usm_arr,
dst=piv_h.get_array(),
sycl_queue=piv.sycl_queue,
depends=dep_evs,
)
_manager.add_event_pair(ht_ev, piv_copy_ev)
# MKL lapack uses 1-origin while SciPy uses 0-origin
piv_h += 1
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

it'd make sense to do before line#2581, to build a list of copy events dep_evs = [lu_copy_ev, piv_copy_ev, b_copy_ev]

depends=dep_evs,
)
_manager.add_event_pair(ht_ev, piv_copy_ev)
# MKL lapack uses 1-origin while SciPy uses 0-origin
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I guess SciPy also uses MKL to call getrs, so it seems unclear for me.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants