SciPy

scipy.linalg.interpolative.interp_decomp

scipy.linalg.interpolative.interp_decomp(A, eps_or_k, rand=True)[source]

Compute ID of a matrix.

An ID of a matrix A is a factorization defined by a rank k, a column index array idx, and interpolation coefficients proj such that:

numpy.dot(A[:,idx[:k]], proj) = A[:,idx[k:]]

The original matrix can then be reconstructed as:

numpy.hstack([A[:,idx[:k]],
                            numpy.dot(A[:,idx[:k]], proj)]
                        )[:,numpy.argsort(idx)]

or via the routine reconstruct_matrix_from_id. This can equivalently be written as:

numpy.dot(A[:,idx[:k]],
                    numpy.hstack([numpy.eye(k), proj])
                  )[:,np.argsort(idx)]

in terms of the skeleton and interpolation matrices:

B = A[:,idx[:k]]

and:

P = numpy.hstack([numpy.eye(k), proj])[:,np.argsort(idx)]

respectively. See also reconstruct_interp_matrix and reconstruct_skel_matrix.

The ID can be computed to any relative precision or rank (depending on the value of eps_or_k). If a precision is specified (eps_or_k < 1), then this function has the output signature:

k, idx, proj = interp_decomp(A, eps_or_k)

Otherwise, if a rank is specified (eps_or_k >= 1), then the output signature is:

idx, proj = interp_decomp(A, eps_or_k)
Parameters:

A : numpy.ndarray or scipy.sparse.linalg.LinearOperator with rmatvec

Matrix to be factored

eps_or_k : float or int

Relative error (if eps_or_k < 1) or rank (if eps_or_k >= 1) of approximation.

rand : bool, optional

Whether to use random sampling if A is of type numpy.ndarray (randomized algorithms are always used if A is of type scipy.sparse.linalg.LinearOperator).

Returns:

k : int

Rank required to achieve specified relative precision if eps_or_k < 1.

idx : numpy.ndarray

Column index array.

proj : numpy.ndarray

Interpolation coefficients.