quimb.linalg.base_linalg ======================== .. py:module:: quimb.linalg.base_linalg .. autoapi-nested-parse:: Backend agnostic functions for solving matrices either fully or partially. Attributes ---------- .. autoapisummary:: quimb.linalg.base_linalg.eigs_slepc quimb.linalg.base_linalg._EIGS_METHODS quimb.linalg.base_linalg.eig quimb.linalg.base_linalg.eigh quimb.linalg.base_linalg.eigvals quimb.linalg.base_linalg.eigvalsh quimb.linalg.base_linalg.eigvecs quimb.linalg.base_linalg.eigvecsh quimb.linalg.base_linalg._SVDS_METHODS quimb.linalg.base_linalg._EXPM_MULTIPLY_METHODS Classes ------- .. autoapisummary:: quimb.linalg.base_linalg.IdentityLinearOperator quimb.linalg.base_linalg.Lazy Functions --------- .. autoapisummary:: quimb.linalg.base_linalg.choose_backend quimb.linalg.base_linalg.eigensystem_partial quimb.linalg.base_linalg.eigensystem quimb.linalg.base_linalg.eigenvectors quimb.linalg.base_linalg.groundstate quimb.linalg.base_linalg.groundenergy quimb.linalg.base_linalg.bound_spectrum quimb.linalg.base_linalg._rel_window_to_abs_window quimb.linalg.base_linalg.eigh_window quimb.linalg.base_linalg.eigvalsh_window quimb.linalg.base_linalg.eigvecsh_window quimb.linalg.base_linalg.svd quimb.linalg.base_linalg.svds quimb.linalg.base_linalg.norm_2 quimb.linalg.base_linalg.norm_fro_dense quimb.linalg.base_linalg.norm_fro_sparse quimb.linalg.base_linalg.norm_trace_dense quimb.linalg.base_linalg.norm quimb.linalg.base_linalg.expm quimb.linalg.base_linalg.expm_multiply quimb.linalg.base_linalg.sqrtm Module Contents --------------- .. py:data:: eigs_slepc .. py:function:: choose_backend(A, k, int_eps=False, B=None) Pick a backend automatically for partial decompositions. .. py:data:: _EIGS_METHODS .. py:function:: 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) Return a few eigenpairs from an operator. :param A: The operator to solve for. :type A: sparse, dense or linear operator :param k: Number of eigenpairs to return. :type k: int :param isherm: Whether to use hermitian solve or not. :type isherm: bool :param B: If given, the RHS operator defining a generalized eigen problem. :type B: sparse, dense or linear operator, optional :param which: Where in spectrum to take eigenvalues from (see :func:``scipy.sparse.linalg.eigsh``) :type which: {'SA', 'LA', 'LM', 'SM', 'TR'} :param return_vecs: Whether to return the eigenvectors. :type return_vecs: bool, optional :param sigma: Which part of spectrum to target, implies which='TR' if which is None. :type sigma: float, optional :param ncv: number of lanczos vectors, can use to optimise speed :type ncv: int, optional :param tol: Tolerance with which to find eigenvalues. :type tol: None or float :param v0: An initial vector guess to iterate with. :type v0: None or 1D-array like :param sort: Whether to explicitly sort by ascending eigenvalue order. :type sort: bool, optional :param backend: 'LOBPCG', 'SLEPC', 'SLEPC-NOMPI'}, optional Which solver to use. :type backend: {'AUTO', 'NUMPY', 'SCIPY', :param fallback_to_scipy: If an error occurs and scipy is not being used, try using scipy. :type fallback_to_scipy: bool, optional :param backend_opts: Supplied to the backend solver. :returns: * **elk** (*(k,) array*) -- The ``k`` eigenvalues. * **evk** (*(d, k) array*) -- Array with ``k`` eigenvectors as columns if ``return_vecs``. .. py:function:: eigensystem(A, isherm, *, k=-1, sort=True, return_vecs=True, **kwargs) Find all or some eigenpairs of an operator. :param A: The operator to decompose. :type A: operator :param isherm: Whether the operator is assumed to be hermitian or not. :type isherm: bool :param k: If negative, find all eigenpairs, else perform partial eigendecomposition and find ``k`` pairs. See :func:`~quimb.linalg.base_linalg.eigensystem_partial`. :type k: int, optional :param sort: Whether to sort the eigenpairs in ascending eigenvalue order. :type sort: bool, optional :param 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``. .. py:data:: eig .. py:data:: eigh .. py:data:: eigvals .. py:data:: eigvalsh .. py:function:: eigenvectors(A, isherm, *, sort=True, **kwargs) .. py:data:: eigvecs .. py:data:: eigvecsh .. py:function:: groundstate(ham, **kwargs) Alias for finding lowest eigenvector only. .. py:function:: groundenergy(ham, **kwargs) Alias for finding lowest eigenvalue only. .. py:function:: bound_spectrum(A, backend='auto', **kwargs) Return the smallest and largest eigenvalue of hermitian operator ``A``. .. py:function:: _rel_window_to_abs_window(el_min, el_max, w_0, w_sz=None) Convert min/max eigenvalues and relative window to absolute values. :param el_min: Smallest eigenvalue. :type el_min: float :param el_max: Largest eigenvalue. :type el_max: float :param w_0: Relative window centre. :type w_0: float [0.0 - 1.0] :param w_sz: Relative window width. :type w_sz: float, optional :returns: Absolute value of centre of window, lower and upper intervals if a window size is specified. :rtype: l_0[, l_min, l_max] .. py:function:: eigh_window(A, w_0, k, w_sz=None, backend='AUTO', return_vecs=True, offset_const=1 / 104729, **kwargs) Return mid-spectrum eigenpairs from a hermitian operator. :param A: Operator to retrieve eigenpairs from. :type A: (d, d) operator :param w_0: Relative window centre to retrieve eigenpairs from. :type w_0: float [0.0, 1.0] :param k: Target number of eigenpairs to retrieve. :type k: int :param w_sz: Relative maximum window width within which to keep eigenpairs. :type w_sz: float, optional :param backend: Which :func:`~quimb.eigh` backend to use. :type backend: str, optional :param return_vecs: Whether to return eigenvectors as well. :type return_vecs: bool, optional :param offset_const: Small fudge factor (relative to window range) to avoid 1 / 0 issues. :type offset_const: float, optional :returns: * **el** (*(k,) array*) -- Eigenvalues around w_0. * **ev** (*(d, k) array*) -- The eigenvectors, if ``return_vecs=True``. .. py:function:: eigvalsh_window(*args, **kwargs) Alias for only finding the eigenvalues in a relative window. .. py:function:: eigvecsh_window(*args, **kwargs) Alias for only finding the eigenvectors in a relative window. .. py:function:: svd(A, return_vecs=True) Compute full singular value decomposition of an operator, using numpy. :param A: The operator. :type A: (m, n) array :param return_vecs: Whether to return the singular vectors. :type return_vecs: bool, optional :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. .. py:data:: _SVDS_METHODS .. py:function:: svds(A, k, ncv=None, return_vecs=True, backend='AUTO', **kwargs) Compute the partial singular value decomposition of an operator. :param A: The operator to decompose. :type A: dense, sparse or linear operator :param k: number of singular value (triplets) to retrieve :type k: int, optional :param ncv: Number of lanczos vectors to use performing decomposition. :type ncv: int, optional :param return_vecs: Whether to return the left and right vectors :type return_vecs: bool, optional :param backend: Which solver to use to perform decomposition. :type backend: {'AUTO', 'SCIPY', 'SLEPC', 'SLEPC-NOMPI', 'NUMPY'}, optional :returns: Singular value(s) (and vectors) such that ``Uk @ np.diag(sk) @ VHk`` approximates ``A``. :rtype: (Uk,) sk (, VHk) .. py:function:: norm_2(A, **kwargs) Return the 2-norm of operator, ``A``, i.e. the largest singular value. .. py:function:: norm_fro_dense(A) Frobenius norm for dense matrices .. py:function:: norm_fro_sparse(A) .. py:function:: norm_trace_dense(A, isherm=False) Returns the trace norm of operator ``A``, that is, the sum of the absolute eigenvalues. .. py:function:: norm(A, ntype=2, **kwargs) Operator norms. :param A: The operator to find norm of. :type A: operator :param ntype: Norm to calculate, if any of: - {2, '2', 'spectral'}: largest singular value - {'f', 'fro'}: frobenius norm - {'t', 'nuc', 'tr', 'trace'}: sum of singular values :type ntype: str :returns: **x** -- The operator norm. :rtype: float .. py:function:: expm(A, herm=False) Matrix exponential, can be accelerated if explicitly hermitian. :param A: Operator to exponentiate. :type A: dense or sparse operator :param herm: If True (not default), and ``A`` is dense, digonalize the matrix in order to perform the exponential. :type herm: bool, optional .. py:data:: _EXPM_MULTIPLY_METHODS .. py:function:: expm_multiply(mat, vec, backend='AUTO', **kwargs) Compute the action of ``expm(mat)`` on ``vec``. :param mat: Operator with which to act with exponential on ``vec``. :type mat: operator :param vec: Vector to act with exponential of operator on. :type vec: vector-like :param backend: Which backend to use. :type backend: {'AUTO', 'SCIPY', 'SLEPC', 'SLEPC-KRYLOV', 'SLEPC-EXPOKIT'} :param kwargs: Supplied to backend function. :returns: Result of ``expm(mat) @ vec``. :rtype: vector .. py:function:: sqrtm(A, herm=True) Matrix square root, can be accelerated if explicitly hermitian. :param A: Operator to take square root of. :type A: dense array :param herm: If True (the default), and ``A`` is dense, digonalize the matrix in order to take the square root. :type herm: bool, optional :rtype: array .. py:class:: IdentityLinearOperator(size, factor=1) Bases: :py:obj:`scipy.sparse.linalg.LinearOperator` Get a ``LinearOperator`` representation of the identity operator, scaled by ``factor``. :param size: The size of the identity. :type size: int :param factor: The coefficient of the identity. :type factor: float .. rubric:: Examples >>> I3 = IdentityLinearOperator(100, 1/3) >>> p = rand_ket(100) >>> np.allclose(I3 @ p, p / 3) True .. py:attribute:: factor :value: 1 .. py:method:: _matvec(vec) Default matrix-vector multiplication handler. If self is a linear operator of shape (M, N), then this method will be called on a shape (N,) or (N, 1) ndarray, and should return a shape (M,) or (M, 1) ndarray. This default implementation falls back on _matmat, so defining that will define matrix-vector multiplication as well. .. py:method:: _rmatvec(vec) Default implementation of _rmatvec; defers to adjoint. .. py:method:: _matmat(mat) Default matrix-matrix multiplication handler. Falls back on the user-defined _matvec method, so defining that will define matrix multiplication (though in a very suboptimal way). .. py:class:: Lazy(fn, *args, shape=None, factor=None, **kwargs) 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 an ``ownership`` keyword to avoid forming every row. This is essentially like using ``functools.partial`` and assigning the ``shape`` attribute. :param fn: A function that constructs an operator. :type fn: callable :param shape: Shape of the constructed operator. :param args: Supplied to ``fn``. :param kwargs: Supplied to ``fn``. :returns: **Lazy** :rtype: callable .. rubric:: Examples Setup the lazy operator: >>> H_lazy = Lazy(ham_heis, n=10, shape=(2**10, 2**10), sparse=True) >>> H_lazy Build a matrix slice (usually done automatically by e.g. ``eigs``): >>> H_lazy(ownership=(256, 512)) <256x1024 sparse matrix of type '' with 1664 stored elements in Compressed Sparse Row format> .. py:attribute:: fn .. py:attribute:: args :value: () .. py:attribute:: kwargs .. py:attribute:: shape :value: None .. py:attribute:: factor :value: None .. py:attribute:: dtype :value: None .. py:method:: __imul__(x) .. py:method:: __mul__(x) .. py:method:: __rmul__(x) .. py:method:: __call__(**kwargs) .. py:method:: __repr__()