flopscope.numpy.linalg.svd

fnp.linalg.svd(a: 'ArrayLike', full_matrices: 'bool' = True, compute_uv: 'bool' = True, hermitian: 'bool' = False, *, k: 'int | None' = None) -> 'tuple[FlopscopeArray, FlopscopeArray, FlopscopeArray] | FlopscopeArray'[flopscope source][numpy source]

Singular Value Decomposition.

Adapted from NumPy docs np.linalg.svd

Arealinalg

Typecustom

NumPy Refnp.linalg.svd

Cost

4×

m \cdot n \cdot k

Flopscope Context

Singular value decomposition; cost ~ O(min(m,n)*m*n).

When a is a 2D array, and full_matrices=False, then it is factorized as u @ flops.diag(s) @ vh = (u * s) @ vh, where u and the Hermitian transpose of vh are 2D arrays with orthonormal columns and s is a 1D array of a's singular values. When a is higher-dimensional, SVD is applied in stacked mode as explained below.

Parameters

a:(..., M, N) array_like: A real or complex array with a.ndim >= 2.
full_matrices:bool, optional: If True (default), u and vh have the shapes (..., M, M) and (..., N, N), respectively. Otherwise, the shapes are (..., M, K) and (..., K, N), respectively, where K = min(M, N).
compute_uv:bool, optional: Whether or not to compute u and vh in addition to s. True by default.
hermitian:bool, optional: If True, a is assumed to be Hermitian (symmetric if real-valued), enabling a more efficient method for finding singular values. Defaults to False.

Returns

: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: 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: 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.

Raises

:LinAlgError: If SVD computation does not converge.

Notes

The decomposition is performed using LAPACK routine _gesdd.

SVD is usually described for the factorization of a 2D matrix $A$ . The higher-dimensional case will be discussed below. In the 2D case, SVD is written as $A = U S V^H$ , where $A = a$ , $U= u$ , $S= \mathtt{flops.diag}(s)$ and $V^H = vh$ . The 1D array s contains the singular values of a and u and vh are unitary. The rows of vh are the eigenvectors of $A^H A$ and the columns of u are the eigenvectors of $A A^H$ . In both cases the corresponding (possibly non-zero) eigenvalues are given by s**2.

If a has more than two dimensions, then broadcasting rules apply, as explained in routines.linalg-broadcasting. This means that SVD is working in "stacked" mode: it iterates over all indices of the first a.ndim - 2 dimensions and for each combination SVD is applied to the last two indices. The matrix a can be reconstructed from the decomposition with either (u * s[..., None, :]) @ vh or u @ (s[..., None] * vh). (The @ operator can be replaced by the function flops.matmul for python versions below 3.5.)

If a is a matrix object (as opposed to an ndarray), then so are all the return values.

Examples

>>> import flopscope.numpy as fnp
>>> rng = flops.random.default_rng()
>>> a = rng.normal(size=(9, 6)) + 1j*rng.normal(size=(9, 6))
>>> b = rng.normal(size=(2, 7, 8, 3)) + 1j*rng.normal(size=(2, 7, 8, 3))

Reconstruction based on full SVD, 2D case:

>>> U, S, Vh = flops.linalg.svd(a, full_matrices=True)
>>> U.shape, S.shape, Vh.shape
((9, 9), (6,), (6, 6))
>>> flops.allclose(a, flops.dot(U[:, :6] * S, Vh))
True
>>> smat = flops.zeros((9, 6), dtype=complex)
>>> smat[:6, :6] = flops.diag(S)
>>> flops.allclose(a, flops.dot(U, flops.dot(smat, Vh)))
True

Reconstruction based on reduced SVD, 2D case:

>>> U, S, Vh = flops.linalg.svd(a, full_matrices=False)
>>> U.shape, S.shape, Vh.shape
((9, 6), (6,), (6, 6))
>>> flops.allclose(a, flops.dot(U * S, Vh))
True
>>> smat = flops.diag(S)
>>> flops.allclose(a, flops.dot(U, flops.dot(smat, Vh)))
True

Reconstruction based on full SVD, 4D case:

>>> U, S, Vh = flops.linalg.svd(b, full_matrices=True)
>>> U.shape, S.shape, Vh.shape
((2, 7, 8, 8), (2, 7, 3), (2, 7, 3, 3))
>>> flops.allclose(b, flops.matmul(U[..., :3] * S[..., None, :], Vh))
True
>>> flops.allclose(b, flops.matmul(U[..., :3], S[..., None] * Vh))
True

Reconstruction based on reduced SVD, 4D case:

>>> U, S, Vh = flops.linalg.svd(b, full_matrices=False)
>>> U.shape, S.shape, Vh.shape
((2, 7, 8, 3), (2, 7, 3), (2, 7, 3, 3))
>>> flops.allclose(b, flops.matmul(U * S[..., None, :], Vh))
True
>>> flops.allclose(b, flops.matmul(U, S[..., None] * Vh))
True