Merge pull request #22786 from asmeurer/linalg-namedtuplesHEADmain

ENH: Add namedtuple return types to linalg functions that return tuples

1 files changed, 149 insertions, 108 deletions

diff --git a/numpy/linalg/linalg.py b/numpy/linalg/linalg.py
index 78927d1ae..b838b9397 100644
--- a/numpy/linalg/linalg.py
+++ b/numpy/linalg/linalg.py
@@ -17,6 +17,7 @@ __all__ = ['matrix_power', 'solve', 'tensorsolve', 'tensorinv', 'inv',
 import functools
 import operator
 import warnings
+from typing import NamedTuple, Any
 
 from .._utils import set_module
 from numpy.core import (
@@ -33,6 +34,28 @@ from numpy.core import overrides
 from numpy.lib.twodim_base import triu, eye
 from numpy.linalg import _umath_linalg
 
+from numpy._typing import NDArray
+
+class EigResult(NamedTuple):
+    eigenvalues: NDArray[Any]
+    eigenvectors: NDArray[Any]
+
+class EighResult(NamedTuple):
+    eigenvalues: NDArray[Any]
+    eigenvectors: NDArray[Any]
+
+class QRResult(NamedTuple):
+    Q: NDArray[Any]
+    R: NDArray[Any]
+
+class SlogdetResult(NamedTuple):
+    sign: NDArray[Any]
+    logabsdet: NDArray[Any]
+
+class SVDResult(NamedTuple):
+    U: NDArray[Any]
+    S: NDArray[Any]
+    Vh: NDArray[Any]
 
 array_function_dispatch = functools.partial(
     overrides.array_function_dispatch, module='numpy.linalg')
@@ -778,10 +801,9 @@ def qr(a, mode='reduced'):
     mode : {'reduced', 'complete', 'r', 'raw'}, optional
         If K = min(M, N), then
 
-        * 'reduced'  : returns q, r with dimensions
-                       (..., M, K), (..., K, N) (default)
-        * 'complete' : returns q, r with dimensions (..., M, M), (..., M, N)
-        * 'r'        : returns r only with dimensions (..., K, N)
+        * 'reduced'  : returns Q, R with dimensions (..., M, K), (..., K, N) (default)
+        * 'complete' : returns Q, R with dimensions (..., M, M), (..., M, N)
+        * 'r'        : returns R only with dimensions (..., K, N)
         * 'raw'      : returns h, tau with dimensions (..., N, M), (..., K,)
 
         The options 'reduced', 'complete, and 'raw' are new in numpy 1.8,
@@ -797,14 +819,17 @@ def qr(a, mode='reduced'):
 
     Returns
     -------
-    q : ndarray of float or complex, optional
+    When mode is 'reduced' or 'complete', the result will be a namedtuple with
+    the attributes `Q` and `R`.
+
+    Q : ndarray of float or complex, optional
         A matrix with orthonormal columns. When mode = 'complete' the
         result is an orthogonal/unitary matrix depending on whether or not
         a is real/complex. The determinant may be either +/- 1 in that
         case. In case the number of dimensions in the input array is
         greater than 2 then a stack of the matrices with above properties
         is returned.
-    r : ndarray of float or complex, optional
+    R : ndarray of float or complex, optional
         The upper-triangular matrix or a stack of upper-triangular
         matrices if the number of dimensions in the input array is greater
         than 2.
@@ -849,19 +874,19 @@ def qr(a, mode='reduced'):
     Examples
     --------
     >>> a = np.random.randn(9, 6)
-    >>> q, r = np.linalg.qr(a)
-    >>> np.allclose(a, np.dot(q, r))  # a does equal qr
+    >>> Q, R = np.linalg.qr(a)
+    >>> np.allclose(a, np.dot(Q, R))  # a does equal QR
     True
-    >>> r2 = np.linalg.qr(a, mode='r')
-    >>> np.allclose(r, r2)  # mode='r' returns the same r as mode='full'
+    >>> R2 = np.linalg.qr(a, mode='r')
+    >>> np.allclose(R, R2)  # mode='r' returns the same R as mode='full'
     True
     >>> a = np.random.normal(size=(3, 2, 2)) # Stack of 2 x 2 matrices as input
-    >>> q, r = np.linalg.qr(a)
-    >>> q.shape
+    >>> Q, R = np.linalg.qr(a)
+    >>> Q.shape
     (3, 2, 2)
-    >>> r.shape
+    >>> R.shape
     (3, 2, 2)
-    >>> np.allclose(a, np.matmul(q, r))
+    >>> np.allclose(a, np.matmul(Q, R))
     True
 
     Example illustrating a common use of `qr`: solving of least squares
@@ -876,8 +901,8 @@ def qr(a, mode='reduced'):
       x = array([[y0], [m]])
       b = array([[1], [0], [2], [1]])
 
-    If A = qr such that q is orthonormal (which is always possible via
-    Gram-Schmidt), then ``x = inv(r) * (q.T) * b``.  (In numpy practice,
+    If A = QR such that Q is orthonormal (which is always possible via
+    Gram-Schmidt), then ``x = inv(R) * (Q.T) * b``.  (In numpy practice,
     however, we simply use `lstsq`.)
 
     >>> A = np.array([[0, 1], [1, 1], [1, 1], [2, 1]])
@@ -887,9 +912,9 @@ def qr(a, mode='reduced'):
            [1, 1],
            [2, 1]])
     >>> b = np.array([1, 2, 2, 3])
-    >>> q, r = np.linalg.qr(A)
-    >>> p = np.dot(q.T, b)
-    >>> np.dot(np.linalg.inv(r), p)
+    >>> Q, R = np.linalg.qr(A)
+    >>> p = np.dot(Q.T, b)
+    >>> np.dot(np.linalg.inv(R), p)
     array([  1.,   1.])
 
     """
@@ -961,7 +986,7 @@ def qr(a, mode='reduced'):
     q = q.astype(result_t, copy=False)
     r = r.astype(result_t, copy=False)
 
-    return wrap(q), wrap(r)
+    return QRResult(wrap(q), wrap(r))
 
 # Eigenvalues
 
@@ -1178,7 +1203,9 @@ def eig(a):
 
     Returns
     -------
-    w : (..., M) array
+    A namedtuple with the following attributes:
+
+    eigenvalues : (..., M) array
         The eigenvalues, each repeated according to its multiplicity.
         The eigenvalues are not necessarily ordered. The resulting
         array will be of complex type, unless the imaginary part is
@@ -1186,10 +1213,10 @@ def eig(a):
         is real the resulting eigenvalues will be real (0 imaginary
         part) or occur in conjugate pairs
 
-    v : (..., M, M) array
+    eigenvectors : (..., M, M) array
         The normalized (unit "length") eigenvectors, such that the
-        column ``v[:,i]`` is the eigenvector corresponding to the
-        eigenvalue ``w[i]``.
+        column ``eigenvectors[:,i]`` is the eigenvector corresponding to the
+        eigenvalue ``eigenvalues[i]``.
 
     Raises
     ------
@@ -1219,30 +1246,30 @@ def eig(a):
     This is implemented using the ``_geev`` LAPACK routines which compute
     the eigenvalues and eigenvectors of general square arrays.
 
-    The number `w` is an eigenvalue of `a` if there exists a vector
-    `v` such that ``a @ v = w * v``. Thus, the arrays `a`, `w`, and
-    `v` satisfy the equations ``a @ v[:,i] = w[i] * v[:,i]``
-    for :math:`i \\in \\{0,...,M-1\\}`.
+    The number `w` is an eigenvalue of `a` if there exists a vector `v` such
+    that ``a @ v = w * v``. Thus, the arrays `a`, `eigenvalues`, and
+    `eigenvectors` satisfy the equations ``a @ eigenvectors[:,i] =
+    eigenvalues[i] * eigenvalues[:,i]`` for :math:`i \\in \\{0,...,M-1\\}`.
 
-    The array `v` of eigenvectors may not be of maximum rank, that is, some
-    of the columns may be linearly dependent, although round-off error may
-    obscure that fact. If the eigenvalues are all different, then theoretically
-    the eigenvectors are linearly independent and `a` can be diagonalized by
-    a similarity transformation using `v`, i.e, ``inv(v) @ a @ v`` is diagonal.
+    The array `eigenvectors` may not be of maximum rank, that is, some of the
+    columns may be linearly dependent, although round-off error may obscure
+    that fact. If the eigenvalues are all different, then theoretically the
+    eigenvectors are linearly independent and `a` can be diagonalized by a
+    similarity transformation using `eigenvectors`, i.e, ``inv(eigenvectors) @
+    a @ eigenvectors`` is diagonal.
 
     For non-Hermitian normal matrices the SciPy function `scipy.linalg.schur`
-    is preferred because the matrix `v` is guaranteed to be unitary, which is
-    not the case when using `eig`. The Schur factorization produces an
-    upper triangular matrix rather than a diagonal matrix, but for normal
-    matrices only the diagonal of the upper triangular matrix is needed, the
-    rest is roundoff error.
-
-    Finally, it is emphasized that `v` consists of the *right* (as in
-    right-hand side) eigenvectors of `a`.  A vector `y` satisfying
-    ``y.T @ a = z * y.T`` for some number `z` is called a *left*
-    eigenvector of `a`, and, in general, the left and right eigenvectors
-    of a matrix are not necessarily the (perhaps conjugate) transposes
-    of each other.
+    is preferred because the matrix `eigenvectors` is guaranteed to be
+    unitary, which is not the case when using `eig`. The Schur factorization
+    produces an upper triangular matrix rather than a diagonal matrix, but for
+    normal matrices only the diagonal of the upper triangular matrix is
+    needed, the rest is roundoff error.
+
+    Finally, it is emphasized that `eigenvectors` consists of the *right* (as
+    in right-hand side) eigenvectors of `a`. A vector `y` satisfying ``y.T @ a
+    = z * y.T`` for some number `z` is called a *left* eigenvector of `a`,
+    and, in general, the left and right eigenvectors of a matrix are not
+    necessarily the (perhaps conjugate) transposes of each other.
 
     References
     ----------
@@ -1253,41 +1280,45 @@ def eig(a):
     --------
     >>> from numpy import linalg as LA
 
-    (Almost) trivial example with real e-values and e-vectors.
+    (Almost) trivial example with real eigenvalues and eigenvectors.
 
-    >>> w, v = LA.eig(np.diag((1, 2, 3)))
-    >>> w; v
+    >>> eigenvalues, eigenvectors = LA.eig(np.diag((1, 2, 3)))
+    >>> eigenvalues
     array([1., 2., 3.])
+    >>> eigenvectors
     array([[1., 0., 0.],
            [0., 1., 0.],
            [0., 0., 1.]])
 
-    Real matrix possessing complex e-values and e-vectors; note that the
-    e-values are complex conjugates of each other.
+    Real matrix possessing complex eigenvalues and eigenvectors; note that the
+    eigenvalues are complex conjugates of each other.
 
-    >>> w, v = LA.eig(np.array([[1, -1], [1, 1]]))
-    >>> w; v
+    >>> eigenvalues, eigenvectors = LA.eig(np.array([[1, -1], [1, 1]]))
+    >>> eigenvalues
     array([1.+1.j, 1.-1.j])
+    >>> eigenvectors
     array([[0.70710678+0.j        , 0.70710678-0.j        ],
            [0.        -0.70710678j, 0.        +0.70710678j]])
 
-    Complex-valued matrix with real e-values (but complex-valued e-vectors);
+    Complex-valued matrix with real eigenvalues (but complex-valued eigenvectors);
     note that ``a.conj().T == a``, i.e., `a` is Hermitian.
 
     >>> a = np.array([[1, 1j], [-1j, 1]])
-    >>> w, v = LA.eig(a)
-    >>> w; v
+    >>> eigenvalues, eigenvectors = LA.eig(a)
+    >>> eigenvalues
     array([2.+0.j, 0.+0.j])
+    >>> eigenvectors
     array([[ 0.        +0.70710678j,  0.70710678+0.j        ], # may vary
            [ 0.70710678+0.j        , -0.        +0.70710678j]])
 
     Be careful about round-off error!
 
     >>> a = np.array([[1 + 1e-9, 0], [0, 1 - 1e-9]])
-    >>> # Theor. e-values are 1 +/- 1e-9
-    >>> w, v = LA.eig(a)
-    >>> w; v
+    >>> # Theor. eigenvalues are 1 +/- 1e-9
+    >>> eigenvalues, eigenvectors = LA.eig(a)
+    >>> eigenvalues
     array([1., 1.])
+    >>> eigenvectors
     array([[1., 0.],
            [0., 1.]])
 
@@ -1311,7 +1342,7 @@ def eig(a):
         result_t = _complexType(result_t)
 
     vt = vt.astype(result_t, copy=False)
-    return w.astype(result_t, copy=False), wrap(vt)
+    return EigResult(w.astype(result_t, copy=False), wrap(vt))
 
 
 @array_function_dispatch(_eigvalsh_dispatcher)
@@ -1339,13 +1370,15 @@ def eigh(a, UPLO='L'):
 
     Returns
     -------
-    w : (..., M) ndarray
+    A namedtuple with the following attributes:
+
+    eigenvalues : (..., M) ndarray
         The eigenvalues in ascending order, each repeated according to
         its multiplicity.
-    v : {(..., M, M) ndarray, (..., M, M) matrix}
-        The column ``v[:, i]`` is the normalized eigenvector corresponding
-        to the eigenvalue ``w[i]``.  Will return a matrix object if `a` is
-        a matrix object.
+    eigenvectors : {(..., M, M) ndarray, (..., M, M) matrix}
+        The column ``eigenvectors[:, i]`` is the normalized eigenvector
+        corresponding to the eigenvalue ``eigenvalues[i]``.  Will return a
+        matrix object if `a` is a matrix object.
 
     Raises
     ------
@@ -1372,10 +1405,10 @@ def eigh(a, UPLO='L'):
     The eigenvalues/eigenvectors are computed using LAPACK routines ``_syevd``,
     ``_heevd``.
 
-    The eigenvalues of real symmetric or complex Hermitian matrices are
-    always real. [1]_ The array `v` of (column) eigenvectors is unitary
-    and `a`, `w`, and `v` satisfy the equations
-    ``dot(a, v[:, i]) = w[i] * v[:, i]``.
+    The eigenvalues of real symmetric or complex Hermitian matrices are always
+    real. [1]_ The array `eigenvalues` of (column) eigenvectors is unitary and
+    `a`, `eigenvalues`, and `eigenvectors` satisfy the equations ``dot(a,
+    eigenvectors[:, i]) = eigenvalues[i] * eigenvectors[:, i]``.
 
     References
     ----------
@@ -1389,24 +1422,26 @@ def eigh(a, UPLO='L'):
     >>> a
     array([[ 1.+0.j, -0.-2.j],
            [ 0.+2.j,  5.+0.j]])
-    >>> w, v = LA.eigh(a)
-    >>> w; v
+    >>> eigenvalues, eigenvectors = LA.eigh(a)
+    >>> eigenvalues
     array([0.17157288, 5.82842712])
+    >>> eigenvectors
     array([[-0.92387953+0.j        , -0.38268343+0.j        ], # may vary
            [ 0.        +0.38268343j,  0.        -0.92387953j]])
 
-    >>> np.dot(a, v[:, 0]) - w[0] * v[:, 0] # verify 1st e-val/vec pair
+    >>> np.dot(a, eigenvectors[:, 0]) - eigenvalues[0] * eigenvectors[:, 0] # verify 1st eigenval/vec pair
     array([5.55111512e-17+0.0000000e+00j, 0.00000000e+00+1.2490009e-16j])
-    >>> np.dot(a, v[:, 1]) - w[1] * v[:, 1] # verify 2nd e-val/vec pair
+    >>> np.dot(a, eigenvectors[:, 1]) - eigenvalues[1] * eigenvectors[:, 1] # verify 2nd eigenval/vec pair
     array([0.+0.j, 0.+0.j])
 
     >>> A = np.matrix(a) # what happens if input is a matrix object
     >>> A
     matrix([[ 1.+0.j, -0.-2.j],
             [ 0.+2.j,  5.+0.j]])
-    >>> w, v = LA.eigh(A)
-    >>> w; v
+    >>> eigenvalues, eigenvectors = LA.eigh(A)
+    >>> eigenvalues
     array([0.17157288, 5.82842712])
+    >>> eigenvectors
     matrix([[-0.92387953+0.j        , -0.38268343+0.j        ], # may vary
             [ 0.        +0.38268343j,  0.        -0.92387953j]])
 
@@ -1430,6 +1465,7 @@ def eigh(a, UPLO='L'):
            [ 0.        +0.89442719j,  0.        -0.4472136j ]])
     array([[ 0.89442719+0.j       , -0.        +0.4472136j],
            [-0.        +0.4472136j,  0.89442719+0.j       ]])
+
     """
     UPLO = UPLO.upper()
     if UPLO not in ('L', 'U'):
@@ -1451,7 +1487,7 @@ def eigh(a, UPLO='L'):
     w, vt = gufunc(a, signature=signature, extobj=extobj)
     w = w.astype(_realType(result_t), copy=False)
     vt = vt.astype(result_t, copy=False)
-    return w, wrap(vt)
+    return EighResult(w, wrap(vt))
 
 
 # Singular value decomposition
@@ -1493,16 +1529,19 @@ def svd(a, full_matrices=True, compute_uv=True, hermitian=False):
 
     Returns
     -------
-    u : { (..., M, M), (..., M, K) } array
+    When `compute_uv` is True, the result is a namedtuple with the following
+    attribute names:
+
+    U : { (..., M, M), (..., M, K) } array
         Unitary array(s). The first ``a.ndim - 2`` dimensions have the same
         size as those of the input `a`. The size of the last two dimensions
         depends on the value of `full_matrices`. Only returned when
         `compute_uv` is True.
-    s : (..., K) array
+    S : (..., K) array
         Vector(s) with the singular values, within each vector sorted in
         descending order. The first ``a.ndim - 2`` dimensions have the same
         size as those of the input `a`.
-    vh : { (..., N, N), (..., K, N) } array
+    Vh : { (..., N, N), (..., K, N) } array
         Unitary array(s). The first ``a.ndim - 2`` dimensions have the same
         size as those of the input `a`. The size of the last two dimensions
         depends on the value of `full_matrices`. Only returned when
@@ -1555,45 +1594,45 @@ def svd(a, full_matrices=True, compute_uv=True, hermitian=False):
 
     Reconstruction based on full SVD, 2D case:
 
-    >>> u, s, vh = np.linalg.svd(a, full_matrices=True)
-    >>> u.shape, s.shape, vh.shape
+    >>> U, S, Vh = np.linalg.svd(a, full_matrices=True)
+    >>> U.shape, S.shape, Vh.shape
     ((9, 9), (6,), (6, 6))
-    >>> np.allclose(a, np.dot(u[:, :6] * s, vh))
+    >>> np.allclose(a, np.dot(U[:, :6] * S, Vh))
     True
     >>> smat = np.zeros((9, 6), dtype=complex)
-    >>> smat[:6, :6] = np.diag(s)
-    >>> np.allclose(a, np.dot(u, np.dot(smat, vh)))
+    >>> smat[:6, :6] = np.diag(S)
+    >>> np.allclose(a, np.dot(U, np.dot(smat, Vh)))
     True
 
     Reconstruction based on reduced SVD, 2D case:
 
-    >>> u, s, vh = np.linalg.svd(a, full_matrices=False)
-    >>> u.shape, s.shape, vh.shape
+    >>> U, S, Vh = np.linalg.svd(a, full_matrices=False)
+    >>> U.shape, S.shape, Vh.shape
     ((9, 6), (6,), (6, 6))
-    >>> np.allclose(a, np.dot(u * s, vh))
+    >>> np.allclose(a, np.dot(U * S, Vh))
     True
-    >>> smat = np.diag(s)
-    >>> np.allclose(a, np.dot(u, np.dot(smat, vh)))
+    >>> smat = np.diag(S)
+    >>> np.allclose(a, np.dot(U, np.dot(smat, Vh)))
     True
 
     Reconstruction based on full SVD, 4D case:
 
-    >>> u, s, vh = np.linalg.svd(b, full_matrices=True)
-    >>> u.shape, s.shape, vh.shape
+    >>> U, S, Vh = np.linalg.svd(b, full_matrices=True)
+    >>> U.shape, S.shape, Vh.shape
     ((2, 7, 8, 8), (2, 7, 3), (2, 7, 3, 3))
-    >>> np.allclose(b, np.matmul(u[..., :3] * s[..., None, :], vh))
+    >>> np.allclose(b, np.matmul(U[..., :3] * S[..., None, :], Vh))
     True
-    >>> np.allclose(b, np.matmul(u[..., :3], s[..., None] * vh))
+    >>> np.allclose(b, np.matmul(U[..., :3], S[..., None] * Vh))
     True
 
     Reconstruction based on reduced SVD, 4D case:
 
-    >>> u, s, vh = np.linalg.svd(b, full_matrices=False)
-    >>> u.shape, s.shape, vh.shape
+    >>> U, S, Vh = np.linalg.svd(b, full_matrices=False)
+    >>> U.shape, S.shape, Vh.shape
     ((2, 7, 8, 3), (2, 7, 3), (2, 7, 3, 3))
-    >>> np.allclose(b, np.matmul(u * s[..., None, :], vh))
+    >>> np.allclose(b, np.matmul(U * S[..., None, :], Vh))
     True
-    >>> np.allclose(b, np.matmul(u, s[..., None] * vh))
+    >>> np.allclose(b, np.matmul(U, S[..., None] * Vh))
     True
 
     """
@@ -1614,7 +1653,7 @@ def svd(a, full_matrices=True, compute_uv=True, hermitian=False):
             u = _nx.take_along_axis(u, sidx[..., None, :], axis=-1)
             # singular values are unsigned, move the sign into v
             vt = transpose(u * sgn[..., None, :]).conjugate()
-            return wrap(u), s, wrap(vt)
+            return SVDResult(wrap(u), s, wrap(vt))
         else:
             s = eigvalsh(a)
             s = abs(s)
@@ -1643,7 +1682,7 @@ def svd(a, full_matrices=True, compute_uv=True, hermitian=False):
         u = u.astype(result_t, copy=False)
         s = s.astype(_realType(result_t), copy=False)
         vh = vh.astype(result_t, copy=False)
-        return wrap(u), s, wrap(vh)
+        return SVDResult(wrap(u), s, wrap(vh))
     else:
         if m < n:
             gufunc = _umath_linalg.svd_m
@@ -2012,15 +2051,17 @@ def slogdet(a):
 
     Returns
     -------
+    A namedtuple with the following attributes:
+
     sign : (...) array_like
         A number representing the sign of the determinant. For a real matrix,
         this is 1, 0, or -1. For a complex matrix, this is a complex number
         with absolute value 1 (i.e., it is on the unit circle), or else 0.
-    logdet : (...) array_like
+    logabsdet : (...) array_like
         The natural log of the absolute value of the determinant.
 
-    If the determinant is zero, then `sign` will be 0 and `logdet` will be
-    -Inf. In all cases, the determinant is equal to ``sign * np.exp(logdet)``.
+    If the determinant is zero, then `sign` will be 0 and `logabsdet` will be
+    -Inf. In all cases, the determinant is equal to ``sign * np.exp(logabsdet)``.
 
     See Also
     --------
@@ -2045,10 +2086,10 @@ def slogdet(a):
     The determinant of a 2-D array ``[[a, b], [c, d]]`` is ``ad - bc``:
 
     >>> a = np.array([[1, 2], [3, 4]])
-    >>> (sign, logdet) = np.linalg.slogdet(a)
-    >>> (sign, logdet)
+    >>> (sign, logabsdet) = np.linalg.slogdet(a)
+    >>> (sign, logabsdet)
     (-1, 0.69314718055994529) # may vary
-    >>> sign * np.exp(logdet)
+    >>> sign * np.exp(logabsdet)
     -2.0
 
     Computing log-determinants for a stack of matrices:
@@ -2056,10 +2097,10 @@ def slogdet(a):
     >>> a = np.array([ [[1, 2], [3, 4]], [[1, 2], [2, 1]], [[1, 3], [3, 1]] ])
     >>> a.shape
     (3, 2, 2)
-    >>> sign, logdet = np.linalg.slogdet(a)
-    >>> (sign, logdet)
+    >>> sign, logabsdet = np.linalg.slogdet(a)
+    >>> (sign, logabsdet)
     (array([-1., -1., -1.]), array([ 0.69314718,  1.09861229,  2.07944154]))
-    >>> sign * np.exp(logdet)
+    >>> sign * np.exp(logabsdet)
     array([-2., -3., -8.])
 
     This routine succeeds where ordinary `det` does not:
@@ -2079,7 +2120,7 @@ def slogdet(a):
     sign, logdet = _umath_linalg.slogdet(a, signature=signature)
     sign = sign.astype(result_t, copy=False)
     logdet = logdet.astype(real_t, copy=False)
-    return sign, logdet
+    return SlogdetResult(sign, logdet)
 
 
 @array_function_dispatch(_unary_dispatcher)