quimb.linalg.base_linalg#
Backend agnostic functions for solving matrices either fully or partially.
Functions
|
Return the smallest and largest eigenvalue of hermitian operator |
|
Pick a backend automatically for partial decompositions. |
|
Find all or some eigenpairs of an operator. |
|
Return a few eigenpairs from an operator. |
|
Find all or some eigenpairs of an operator. |
|
Return mid-spectrum eigenpairs from a hermitian operator. |
|
Alias for only finding the eigenvalues in a relative window. |
|
Alias for only finding the eigenvectors in a relative window. |
|
Matrix exponential, can be accelerated if explicitly hermitian. |
|
Compute the action of |
|
Alias for finding lowest eigenvalue only. |
|
Alias for finding lowest eigenvector only. |
|
Operator norms. |
|
Return the 2-norm of operator, |
Frobenius norm for dense matrices |
|
|
|
|
Returns the trace norm of operator |
|
Matrix square root, can be accelerated if explicitly hermitian. |
|
Compute full singular value decomposition of an operator, using numpy. |
|
Compute the partial singular value decomposition of an operator. |
Classes
|
Get a |
|
A simple class representing an unconstructed matrix. |
- class quimb.linalg.base_linalg.IdentityLinearOperator(*args, **kwargs)[source]#
Get a
LinearOperator
representation of the identity operator, scaled byfactor
.Examples
>>> I3 = IdentityLinearOperator(100, 1/3) >>> p = rand_ket(100) >>> np.allclose(I3 @ p, p / 3) True
- class quimb.linalg.base_linalg.Lazy(fn, *args, shape=None, factor=None, **kwargs)[source]#
A simple class representing an unconstructed matrix. This can be passed to, for example, MPI workers, who can then construct the matrix themselves. The main function
fn
should ideally take anownership
keyword to avoid forming every row.This is essentially like using
functools.partial
and assigning theshape
attribute.- Parameters
fn (callable) – A function that constructs an operator.
shape – Shape of the constructed operator.
args – Supplied to
fn
.kwargs – Supplied to
fn
.
- Returns
Lazy
- Return type
callable
Examples
Setup the lazy operator:
>>> H_lazy = Lazy(ham_heis, n=10, shape=(2**10, 2**10), sparse=True) >>> H_lazy <Lazy(ham_heis, shape=(1024, 1024), dtype=None)>
Build a matrix slice (usually done automatically by e.g.
eigs
):>>> H_lazy(ownership=(256, 512)) <256x1024 sparse matrix of type '<class 'numpy.float64'>' with 1664 stored elements in Compressed Sparse Row format>
- quimb.linalg.base_linalg.bound_spectrum(A, backend='auto', **kwargs)[source]#
Return the smallest and largest eigenvalue of hermitian operator
A
.
- quimb.linalg.base_linalg.choose_backend(A, k, int_eps=False, B=None)[source]#
Pick a backend automatically for partial decompositions.
- quimb.linalg.base_linalg.eig(A, *, isherm=False, k=- 1, sort=True, return_vecs=True, **kwargs)#
Find all or some eigenpairs of an operator.
- Parameters
A (operator) – The operator to decompose.
isherm (bool) – Whether the operator is assumed to be hermitian or not.
k (int, optional) – If negative, find all eigenpairs, else perform partial eigendecomposition and find
k
pairs. Seeeigensystem_partial()
.sort (bool, optional) – Whether to sort the eigenpairs in ascending eigenvalue order.
kwargs – Supplied to the backend function.
- Returns
el ((k,) array) – Eigenvalues.
ev ((d, k) array) – Corresponding eigenvectors as columns of array, such that
ev @ diag(el) @ ev.H == A
.
- quimb.linalg.base_linalg.eigensystem(A, isherm, *, k=- 1, sort=True, return_vecs=True, **kwargs)[source]#
Find all or some eigenpairs of an operator.
- Parameters
A (operator) – The operator to decompose.
isherm (bool) – Whether the operator is assumed to be hermitian or not.
k (int, optional) – If negative, find all eigenpairs, else perform partial eigendecomposition and find
k
pairs. Seeeigensystem_partial()
.sort (bool, optional) – Whether to sort the eigenpairs in ascending eigenvalue order.
kwargs – Supplied to the backend function.
- Returns
el ((k,) array) – Eigenvalues.
ev ((d, k) array) – Corresponding eigenvectors as columns of array, such that
ev @ diag(el) @ ev.H == A
.
- quimb.linalg.base_linalg.eigensystem_partial(A, k, isherm, *, B=None, which=None, return_vecs=True, sigma=None, ncv=None, tol=None, v0=None, sort=True, backend=None, fallback_to_scipy=False, **backend_opts)[source]#
Return a few eigenpairs from an operator.
- Parameters
A (sparse, dense or linear operator) – The operator to solve for.
k (int) – Number of eigenpairs to return.
isherm (bool) – Whether to use hermitian solve or not.
B (sparse, dense or linear operator, optional) – If given, the RHS operator defining a generalized eigen problem.
which ({'SA', 'LA', 'LM', 'SM', 'TR'}) – Where in spectrum to take eigenvalues from (see :func:
scipy.sparse.linalg.eigsh
)return_vecs (bool, optional) – Whether to return the eigenvectors.
sigma (float, optional) – Which part of spectrum to target, implies which=’TR’ if which is None.
ncv (int, optional) – number of lanczos vectors, can use to optimise speed
tol (None or float) – Tolerance with which to find eigenvalues.
v0 (None or 1D-array like) – An initial vector guess to iterate with.
sort (bool, optional) – Whether to explicitly sort by ascending eigenvalue order.
backend ({'AUTO', 'NUMPY', 'SCIPY',) – ‘LOBPCG’, ‘SLEPC’, ‘SLEPC-NOMPI’}, optional Which solver to use.
fallback_to_scipy (bool, optional) – If an error occurs and scipy is not being used, try using scipy.
backend_opts – Supplied to the backend solver.
- Returns
elk ((k,) array) – The
k
eigenvalues.evk ((d, k) array) – Array with
k
eigenvectors as columns ifreturn_vecs
.
- quimb.linalg.base_linalg.eigenvectors(A, isherm, *, k=- 1, sort=True, return_vecs=True, **kwargs)[source]#
Find all or some eigenpairs of an operator.
- Parameters
A (operator) – The operator to decompose.
isherm (bool) – Whether the operator is assumed to be hermitian or not.
k (int, optional) – If negative, find all eigenpairs, else perform partial eigendecomposition and find
k
pairs. Seeeigensystem_partial()
.sort (bool, optional) – Whether to sort the eigenpairs in ascending eigenvalue order.
kwargs – Supplied to the backend function.
- Returns
el ((k,) array) – Eigenvalues.
ev ((d, k) array) – Corresponding eigenvectors as columns of array, such that
ev @ diag(el) @ ev.H == A
.
- quimb.linalg.base_linalg.eigh(A, *, isherm=True, k=- 1, sort=True, return_vecs=True, **kwargs)#
Find all or some eigenpairs of an operator.
- Parameters
A (operator) – The operator to decompose.
isherm (bool) – Whether the operator is assumed to be hermitian or not.
k (int, optional) – If negative, find all eigenpairs, else perform partial eigendecomposition and find
k
pairs. Seeeigensystem_partial()
.sort (bool, optional) – Whether to sort the eigenpairs in ascending eigenvalue order.
kwargs – Supplied to the backend function.
- Returns
el ((k,) array) – Eigenvalues.
ev ((d, k) array) – Corresponding eigenvectors as columns of array, such that
ev @ diag(el) @ ev.H == A
.
- quimb.linalg.base_linalg.eigh_window(A, w_0, k, w_sz=None, backend='AUTO', return_vecs=True, offset_const=9.548453627934956e-06, **kwargs)[source]#
Return mid-spectrum eigenpairs from a hermitian operator.
- Parameters
A ((d, d) operator) – Operator to retrieve eigenpairs from.
w_0 (float [0.0, 1.0]) – Relative window centre to retrieve eigenpairs from.
k (int) – Target number of eigenpairs to retrieve.
w_sz (float, optional) – Relative maximum window width within which to keep eigenpairs.
backend (str, optional) – Which
eigh()
backend to use.return_vecs (bool, optional) – Whether to return eigenvectors as well.
offset_const (float, optional) – Small fudge factor (relative to window range) to avoid 1 / 0 issues.
- Returns
el ((k,) array) – Eigenvalues around w_0.
ev ((d, k) array) – The eigenvectors, if
return_vecs=True
.
- quimb.linalg.base_linalg.eigvals(A, *, isherm=False, k=- 1, sort=True, return_vecs=False, **kwargs)#
Find all or some eigenpairs of an operator.
- Parameters
A (operator) – The operator to decompose.
isherm (bool) – Whether the operator is assumed to be hermitian or not.
k (int, optional) – If negative, find all eigenpairs, else perform partial eigendecomposition and find
k
pairs. Seeeigensystem_partial()
.sort (bool, optional) – Whether to sort the eigenpairs in ascending eigenvalue order.
kwargs – Supplied to the backend function.
- Returns
el ((k,) array) – Eigenvalues.
ev ((d, k) array) – Corresponding eigenvectors as columns of array, such that
ev @ diag(el) @ ev.H == A
.
- quimb.linalg.base_linalg.eigvalsh(A, *, isherm=True, k=- 1, sort=True, return_vecs=False, **kwargs)#
Find all or some eigenpairs of an operator.
- Parameters
A (operator) – The operator to decompose.
isherm (bool) – Whether the operator is assumed to be hermitian or not.
k (int, optional) – If negative, find all eigenpairs, else perform partial eigendecomposition and find
k
pairs. Seeeigensystem_partial()
.sort (bool, optional) – Whether to sort the eigenpairs in ascending eigenvalue order.
kwargs – Supplied to the backend function.
- Returns
el ((k,) array) – Eigenvalues.
ev ((d, k) array) – Corresponding eigenvectors as columns of array, such that
ev @ diag(el) @ ev.H == A
.
- quimb.linalg.base_linalg.eigvalsh_window(*args, **kwargs)[source]#
Alias for only finding the eigenvalues in a relative window.
- quimb.linalg.base_linalg.eigvecs(A, *, isherm=False, k=- 1, sort=True, return_vecs=True, **kwargs)#
Find all or some eigenpairs of an operator.
- Parameters
A (operator) – The operator to decompose.
isherm (bool) – Whether the operator is assumed to be hermitian or not.
k (int, optional) – If negative, find all eigenpairs, else perform partial eigendecomposition and find
k
pairs. Seeeigensystem_partial()
.sort (bool, optional) – Whether to sort the eigenpairs in ascending eigenvalue order.
kwargs – Supplied to the backend function.
- Returns
el ((k,) array) – Eigenvalues.
ev ((d, k) array) – Corresponding eigenvectors as columns of array, such that
ev @ diag(el) @ ev.H == A
.
- quimb.linalg.base_linalg.eigvecsh(A, *, isherm=True, k=- 1, sort=True, return_vecs=True, **kwargs)#
Find all or some eigenpairs of an operator.
- Parameters
A (operator) – The operator to decompose.
isherm (bool) – Whether the operator is assumed to be hermitian or not.
k (int, optional) – If negative, find all eigenpairs, else perform partial eigendecomposition and find
k
pairs. Seeeigensystem_partial()
.sort (bool, optional) – Whether to sort the eigenpairs in ascending eigenvalue order.
kwargs – Supplied to the backend function.
- Returns
el ((k,) array) – Eigenvalues.
ev ((d, k) array) – Corresponding eigenvectors as columns of array, such that
ev @ diag(el) @ ev.H == A
.
- quimb.linalg.base_linalg.eigvecsh_window(*args, **kwargs)[source]#
Alias for only finding the eigenvectors in a relative window.
- quimb.linalg.base_linalg.expm(A, herm=False)[source]#
Matrix exponential, can be accelerated if explicitly hermitian.
- Parameters
A (dense or sparse operator) – Operator to exponentiate.
herm (bool, optional) – If True (not default), and
A
is dense, digonalize the matrix in order to perform the exponential.
- quimb.linalg.base_linalg.expm_multiply(mat, vec, backend='AUTO', **kwargs)[source]#
Compute the action of
expm(mat)
onvec
.- Parameters
mat (operator) – Operator with which to act with exponential on
vec
.vec (vector-like) – Vector to act with exponential of operator on.
backend ({'AUTO', 'SCIPY', 'SLEPC', 'SLEPC-KRYLOV', 'SLEPC-EXPOKIT'}) – Which backend to use.
kwargs – Supplied to backend function.
- Returns
Result of
expm(mat) @ vec
.- Return type
vector
- quimb.linalg.base_linalg.groundenergy(ham, **kwargs)[source]#
Alias for finding lowest eigenvalue only.
- quimb.linalg.base_linalg.groundstate(ham, **kwargs)[source]#
Alias for finding lowest eigenvector only.
- quimb.linalg.base_linalg.norm_2(A, **kwargs)[source]#
Return the 2-norm of operator,
A
, i.e. the largest singular value.
- quimb.linalg.base_linalg.norm_trace_dense(A, isherm=False)[source]#
Returns the trace norm of operator
A
, that is, the sum of the absolute eigenvalues.
- quimb.linalg.base_linalg.sqrtm(A, herm=True)[source]#
Matrix square root, can be accelerated if explicitly hermitian.
- Parameters
A (dense array) – Operator to take square root of.
herm (bool, optional) – If True (the default), and
A
is dense, digonalize the matrix in order to take the square root.
- Return type
array
- quimb.linalg.base_linalg.svd(A, return_vecs=True)[source]#
Compute full singular value decomposition of an operator, using numpy.
- Parameters
A ((m, n) array) – The operator.
return_vecs (bool, optional) – Whether to return the singular vectors.
- Returns
U ((m, k) array) – Left singular vectors (if
return_vecs=True
) as columns.s ((k,) array) – Singular values.
VH ((k, n) array) – Right singular vectors (if
return_vecs=True
) as rows.
- quimb.linalg.base_linalg.svds(A, k, ncv=None, return_vecs=True, backend='AUTO', **kwargs)[source]#
Compute the partial singular value decomposition of an operator.
- Parameters
A (dense, sparse or linear operator) – The operator to decompose.
k (int, optional) – number of singular value (triplets) to retrieve
ncv (int, optional) – Number of lanczos vectors to use performing decomposition.
return_vecs (bool, optional) – Whether to return the left and right vectors
backend ({'AUTO', 'SCIPY', 'SLEPC', 'SLEPC-NOMPI', 'NUMPY'}, optional) – Which solver to use to perform decomposition.
- Returns
Singular value(s) (and vectors) such that
Uk @ np.diag(sk) @ VHk
approximatesA
.- Return type
(Uk,) sk (, VHk)