svds(solver=’propack’) — SciPy v1.15.2 Manual (original) (raw)

scipy.sparse.linalg.svds(A, k=6, ncv=None, tol=0, which='LM', v0=None, maxiter=None, return_singular_vectors=True, solver='arpack', rng=None, options=None)

Partial singular value decomposition of a sparse matrix using PROPACK.

Compute the largest or smallest k singular values and corresponding singular vectors of a sparse matrix A. The order in which the singular values are returned is not guaranteed.

In the descriptions below, let M, N = A.shape.

Parameters:

Asparse matrix or LinearOperator

Matrix to decompose. If A is a LinearOperatorobject, it must define both matvec and rmatvec methods.

kint, default: 6

Number of singular values and singular vectors to compute. Must satisfy 1 <= k <= min(M, N).

ncvint, optional

Ignored.

tolfloat, optional

The desired relative accuracy for computed singular values. Zero (default) means machine precision.

which{‘LM’, ‘SM’}

Which k singular values to find: either the largest magnitude (‘LM’) or smallest magnitude (‘SM’) singular values. Note that choosingwhich='SM' will force the irl option to be set True.

v0ndarray, optional

Starting vector for iterations: must be of length A.shape[0]. If not specified, PROPACK will generate a starting vector.

maxiterint, optional

Maximum number of iterations / maximal dimension of the Krylov subspace. Default is 10 * k.

return_singular_vectors{True, False, “u”, “vh”}

Singular values are always computed and returned; this parameter controls the computation and return of singular vectors.

solver{‘arpack’, ‘propack’, ‘lobpcg’}, optional

This is the solver-specific documentation for solver='propack'.‘arpack’ and‘lobpcg’are also supported.

rngnumpy.random.Generator, optional

Pseudorandom number generator state. When rng is None, a newnumpy.random.Generator is created using entropy from the operating system. Types other than numpy.random.Generator are passed to numpy.random.default_rng to instantiate a Generator.

optionsdict, optional

A dictionary of solver-specific options. No solver-specific options are currently supported; this parameter is reserved for future use.

Returns:

undarray, shape=(M, k)

Unitary matrix having left singular vectors as columns.

sndarray, shape=(k,)

The singular values.

vhndarray, shape=(k, N)

Unitary matrix having right singular vectors as rows.

Notes

This is an interface to the Fortran library PROPACK [1]. The current default is to run with IRL mode disabled unless seeking the smallest singular values/vectors (which='SM').

References

Examples

Construct a matrix A from singular values and vectors.

import numpy as np from scipy.stats import ortho_group from scipy.sparse import csc_array, diags_array from scipy.sparse.linalg import svds rng = np.random.default_rng() orthogonal = csc_array(ortho_group.rvs(10, random_state=rng)) s = [0.0001, 0.001, 3, 4, 5] # singular values u = orthogonal[:, :5] # left singular vectors vT = orthogonal[:, 5:].T # right singular vectors A = u @ diags_array(s) @ vT

With only three singular values/vectors, the SVD approximates the original matrix.

u2, s2, vT2 = svds(A, k=3, solver='propack') A2 = u2 @ np.diag(s2) @ vT2 np.allclose(A2, A.todense(), atol=1e-3) True

With all five singular values/vectors, we can reproduce the original matrix.

u3, s3, vT3 = svds(A, k=5, solver='propack') A3 = u3 @ np.diag(s3) @ vT3 np.allclose(A3, A.todense()) True

The singular values match the expected singular values, and the singular vectors are as expected up to a difference in sign.

(np.allclose(s3, s) and ... np.allclose(np.abs(u3), np.abs(u.toarray())) and ... np.allclose(np.abs(vT3), np.abs(vT.toarray()))) True

The singular vectors are also orthogonal.

(np.allclose(u3.T @ u3, np.eye(5)) and ... np.allclose(vT3 @ vT3.T, np.eye(5))) True