# scipy.linalg.lu

scipy.linalg.lu(a, permute_l=False, overwrite_a=False, check_finite=True, p_indices=False)

  Compute LU decomposition of a matrix with partial pivoting.

  The decomposition satisfies:

      A = P @ L @ U

  where `P` is a permutation matrix, `L` lower triangular with unit diagonal elements, and `U` upper triangular. If permute_l is set to `True` then `L` is returned already permuted and hence satisfying `A = L @ U`.

  Parameters:
    a(M, N) array_like

      Array to decompose

    permute_l bool, optional

      Perform the multiplication P*L (Default: do not permute)

    overwrite_a bool, optional

      Whether to overwrite data in a (may improve performance)

    check_finite bool, optional

      Whether to check that the input matrix contains only finite numbers. Disabling may give a performance gain, but may result in problems (crashes, non-termination) if the inputs do contain infinities or NaNs.

    p_indices bool, optional

      If `True` the permutation information is returned as row indices. The default is `False` for backwards-compatibility reasons.

  Returns:
    (If `permute_l` is ``False``)
    p(…, M, M) ndarray

      Permutation arrays or vectors depending on p_indices

    l(…, M, K) ndarray

      Lower triangular or trapezoidal array with unit diagonal. `K = min(M, N)`

    u(…, K, N) ndarray

      Upper triangular or trapezoidal array

    (If `permute_l` is ``True``)
    pl(…, M, K) ndarray

      Permuted L matrix. `K = min(M, N)`

    u(…, K, N) ndarray

      Upper triangular or trapezoidal array

  Notes

  Permutation matrices are costly since they are nothing but row reorder of `L` and hence indices are strongly recommended to be used instead if the permutation is required. The relation in the 2D case then becomes simply `A = L[P, :] @ U`. In higher dimensions, it is better to use permute_l to avoid complicated indexing tricks.

  In 2D case, if one has the indices however, for some reason, the permutation matrix is still needed then it can be constructed by `np.eye(M)[P, :]`.

  Examples

      >>> import numpy as np
      >>> from scipy.linalg import lu
      >>> A = np.array([[2, 5, 8, 7], [5, 2, 2, 8], [7, 5, 6, 6], [5, 4, 4, 8]])
      >>> p, l, u = lu(A)
      >>> np.allclose(A, p @ l @ u)
      True
      >>> p  # Permutation matrix
      array([[0., 1., 0., 0.],  # Row index 1
             [0., 0., 0., 1.],  # Row index 3
             [1., 0., 0., 0.],  # Row index 0
             [0., 0., 1., 0.]]) # Row index 2
      >>> p, _, _ = lu(A, p_indices=True)
      >>> p
      array([1, 3, 0, 2])  # as given by row indices above
      >>> np.allclose(A, l[p, :] @ u)
      True

  We can also use nd-arrays, for example, a demonstration with 4D array:

      >>> rng = np.random.default_rng()
      >>> A = rng.uniform(low=-4, high=4, size=[3, 2, 4, 8])
      >>> p, l, u = lu(A)
      >>> p.shape, l.shape, u.shape
      ((3, 2, 4, 4), (3, 2, 4, 4), (3, 2, 4, 8))
      >>> np.allclose(A, p @ l @ u)
      True
      >>> PL, U = lu(A, permute_l=True)
      >>> np.allclose(A, PL @ U)
      True