quimb.tensor.decomp =================== .. py:module:: quimb.tensor.decomp .. autoapi-nested-parse:: Functions for decomposing and projecting matrices. Attributes ---------- .. autoapisummary:: quimb.tensor.decomp._CUTOFF_MODE_MAP quimb.tensor.decomp._similarity_compress_fns quimb.tensor.decomp._ISOMETRIZE_METHODS Functions --------- .. autoapisummary:: quimb.tensor.decomp.map_cutoff_mode quimb.tensor.decomp.rdmul quimb.tensor.decomp.rddiv quimb.tensor.decomp.ldmul quimb.tensor.decomp.lddiv quimb.tensor.decomp.dag_numba quimb.tensor.decomp.rdmul_numba quimb.tensor.decomp.rddiv_numba quimb.tensor.decomp.ldmul_numba quimb.tensor.decomp.lddiv_numba quimb.tensor.decomp.sgn quimb.tensor.decomp.sgn_numba quimb.tensor.decomp.sgn_tf quimb.tensor.decomp._trim_and_renorm_svd_result quimb.tensor.decomp.svd_truncated quimb.tensor.decomp._compute_number_svals_to_keep_numba quimb.tensor.decomp._compute_svals_renorm_factor_numba quimb.tensor.decomp._trim_and_renorm_svd_result_numba quimb.tensor.decomp.svd_truncated_numba quimb.tensor.decomp.svd_truncated_numpy quimb.tensor.decomp.svd_truncated_lazy quimb.tensor.decomp.lu_truncated quimb.tensor.decomp.svdvals quimb.tensor.decomp._svd_via_eig_truncated_numba quimb.tensor.decomp.svd_via_eig_truncated quimb.tensor.decomp.svdvals_eig quimb.tensor.decomp.eigh_truncated quimb.tensor.decomp.eigh_truncated_numba quimb.tensor.decomp._choose_k quimb.tensor.decomp.svds quimb.tensor.decomp.isvd quimb.tensor.decomp._rsvd_numpy quimb.tensor.decomp.rsvd quimb.tensor.decomp.eigsh quimb.tensor.decomp.qr_stabilized quimb.tensor.decomp.qr_stabilized_numba quimb.tensor.decomp.qr_stabilized_lazy quimb.tensor.decomp.lq_stabilized quimb.tensor.decomp.lq_stabilized_numba quimb.tensor.decomp._cholesky_numba quimb.tensor.decomp.cholesky quimb.tensor.decomp.polar_right quimb.tensor.decomp.polar_right_numba quimb.tensor.decomp.polar_left quimb.tensor.decomp.polar_left_numba quimb.tensor.decomp._similarity_compress_eig quimb.tensor.decomp._similarity_compress_eig_numba quimb.tensor.decomp._similarity_compress_eigh quimb.tensor.decomp._similarity_compress_eigh_numba quimb.tensor.decomp._similarity_compress_svd quimb.tensor.decomp._similarity_compress_svd_numba quimb.tensor.decomp._similarity_compress_biorthog quimb.tensor.decomp._similarity_compress_biorthog_numba quimb.tensor.decomp.similarity_compress quimb.tensor.decomp.isometrize_qr quimb.tensor.decomp.isometrize_svd quimb.tensor.decomp.isometrize_exp quimb.tensor.decomp.isometrize_cayley quimb.tensor.decomp.isometrize_modified_gram_schmidt quimb.tensor.decomp.isometrize_householder quimb.tensor.decomp.isometrize_torch_householder quimb.tensor.decomp.isometrize quimb.tensor.decomp.squared_op_to_reduced_factor quimb.tensor.decomp.squared_op_to_reduced_factor_numba quimb.tensor.decomp.compute_oblique_projectors Module Contents --------------- .. py:data:: _CUTOFF_MODE_MAP .. py:function:: map_cutoff_mode(cutoff_mode) Map mode to an integer for compatibility with numba. .. py:function:: rdmul(x, d) Right-multiplication a matrix by a vector representing a diagonal. .. py:function:: rddiv(x, d) Right-multiplication of a matrix by a vector representing an inverse diagonal. .. py:function:: ldmul(d, x) Left-multiplication a matrix by a vector representing a diagonal. .. py:function:: lddiv(d, x) Left-multiplication of a matrix by a vector representing an inverse diagonal. .. py:function:: dag_numba(x) .. py:function:: rdmul_numba(x, d) .. py:function:: rddiv_numba(x, d) .. py:function:: ldmul_numba(d, x) .. py:function:: lddiv_numba(d, x) .. py:function:: sgn(x) Get the 'sign' of ``x``, such that ``x / sgn(x)`` is real and non-negative. .. py:function:: sgn_numba(x) .. py:function:: sgn_tf(x) .. py:function:: _trim_and_renorm_svd_result(U, s, VH, cutoff, cutoff_mode, max_bond, absorb, renorm) Give full SVD decomposion result ``U``, ``s``, ``VH``, optionally trim, renormalize, and absorb the singular values. See ``svd_truncated`` for details. .. py:function:: svd_truncated(x, cutoff=-1.0, cutoff_mode=4, max_bond=-1, absorb=0, renorm=0, backend=None) Truncated svd or raw array ``x``. :param cutoff: Singular value cutoff threshold, if ``cutoff <= 0.0``, then only ``max_bond`` is used. :type cutoff: float, optional :param cutoff_mode: How to perform the trim: - 1: ['abs'], trim values below ``cutoff`` - 2: ['rel'], trim values below ``s[0] * cutoff`` - 3: ['sum2'], trim s.t. ``sum(s_trim**2) < cutoff``. - 4: ['rsum2'], trim s.t. ``sum(s_trim**2) < sum(s**2) * cutoff``. - 5: ['sum1'], trim s.t. ``sum(s_trim**1) < cutoff``. - 6: ['rsum1'], trim s.t. ``sum(s_trim**1) < sum(s**1) * cutoff``. :type cutoff_mode: {1, 2, 3, 4, 5, 6}, optional :param max_bond: An explicit maximum bond dimension, use -1 for none. :type max_bond: int, optional :param absorb: How to absorb the singular values. -1: left, 0: both, 1: right and None: don't absorb (return). :type absorb: {-1, 0, 1, None}, optional :param renorm: Whether to renormalize the singular values (depends on `cutoff_mode`). :type renorm: {0, 1}, optional .. py:function:: _compute_number_svals_to_keep_numba(s, cutoff, cutoff_mode) Find the number of singular values to keep of ``s`` given ``cutoff`` and ``cutoff_mode``. .. py:function:: _compute_svals_renorm_factor_numba(s, n_chi, renorm) Find the normalization constant for ``s`` such that the new sum squared of the ``n_chi`` largest values equals the sum squared of all the old ones. .. py:function:: _trim_and_renorm_svd_result_numba(U, s, VH, cutoff, cutoff_mode, max_bond, absorb, renorm) Accelerate version of ``_trim_and_renorm_svd_result``. .. py:function:: svd_truncated_numba(x, cutoff=-1.0, cutoff_mode=4, max_bond=-1, absorb=0, renorm=0) Accelerated version of ``svd_truncated`` for numpy arrays. .. py:function:: svd_truncated_numpy(x, cutoff=-1.0, cutoff_mode=4, max_bond=-1, absorb=0, renorm=0) Numpy version of ``svd_truncated``, trying the accelerated version first, then falling back to the more stable scipy version. .. py:function:: svd_truncated_lazy(x, cutoff=-1.0, cutoff_mode=4, max_bond=-1, absorb=0, renorm=0) .. py:function:: lu_truncated(x, cutoff=-1.0, cutoff_mode=4, max_bond=-1, absorb=0, renorm=0, backend=None) .. py:function:: svdvals(x) SVD-decomposition, but return singular values only. .. py:function:: _svd_via_eig_truncated_numba(x, cutoff=-1.0, cutoff_mode=4, max_bond=-1, absorb=0, renorm=0) SVD-split via eigen-decomposition. .. py:function:: svd_via_eig_truncated(x, cutoff=-1.0, cutoff_mode=4, max_bond=-1, absorb=0, renorm=0) .. py:function:: svdvals_eig(x) SVD-decomposition via eigen, but return singular values only. .. py:function:: eigh_truncated(x, cutoff=-1.0, cutoff_mode=4, max_bond=-1, absorb=0, renorm=0, backend=None) .. py:function:: eigh_truncated_numba(x, cutoff=-1.0, cutoff_mode=4, max_bond=-1, absorb=0, renorm=0) SVD-decomposition, using hermitian eigen-decomposition, only works if ``x`` is hermitian. .. py:function:: _choose_k(x, cutoff, max_bond) Choose the number of singular values to target. .. py:function:: svds(x, cutoff=0.0, cutoff_mode=4, max_bond=-1, absorb=0, renorm=0) SVD-decomposition using iterative methods. Allows the computation of only a certain number of singular values, e.g. max_bond, from the get-go, and is thus more efficient. Can also supply ``scipy.sparse.linalg.LinearOperator``. .. py:function:: isvd(x, cutoff=0.0, cutoff_mode=4, max_bond=-1, absorb=0, renorm=0) SVD-decomposition using interpolative matrix random methods. Allows the computation of only a certain number of singular values, e.g. max_bond, from the get-go, and is thus more efficient. Can also supply ``scipy.sparse.linalg.LinearOperator``. .. py:function:: _rsvd_numpy(x, cutoff=0.0, cutoff_mode=4, max_bond=-1, absorb=0, renorm=0) .. py:function:: rsvd(x, cutoff=0.0, cutoff_mode=4, max_bond=-1, absorb=0, renorm=0) SVD-decomposition using randomized methods (due to Halko). Allows the computation of only a certain number of singular values, e.g. max_bond, from the get-go, and is thus more efficient. Can also supply ``scipy.sparse.linalg.LinearOperator``. .. py:function:: eigsh(x, cutoff=0.0, cutoff_mode=4, max_bond=-1, absorb=0, renorm=0) SVD-decomposition using iterative hermitian eigen decomp, thus assuming that ``x`` is hermitian. Allows the computation of only a certain number of singular values, e.g. max_bond, from the get-go, and is thus more efficient. Can also supply ``scipy.sparse.linalg.LinearOperator``. .. py:function:: qr_stabilized(x, backend=None) QR-decomposition, with stabilized R factor. .. py:function:: qr_stabilized_numba(x) .. py:function:: qr_stabilized_lazy(x) .. py:function:: lq_stabilized(x, backend=None) .. py:function:: lq_stabilized_numba(x) .. py:function:: _cholesky_numba(x, cutoff=-1, cutoff_mode=4, max_bond=-1, absorb=0) SVD-decomposition, using cholesky decomposition, only works if ``x`` is positive definite. .. py:function:: cholesky(x, cutoff=-1, cutoff_mode=4, max_bond=-1, absorb=0) .. py:function:: polar_right(x) Polar decomposition of ``x``. .. py:function:: polar_right_numba(x) .. py:function:: polar_left(x) Polar decomposition of ``x``. .. py:function:: polar_left_numba(x) .. py:function:: _similarity_compress_eig(X, max_bond, renorm) .. py:function:: _similarity_compress_eig_numba(X, max_bond, renorm) .. py:function:: _similarity_compress_eigh(X, max_bond, renorm) .. py:function:: _similarity_compress_eigh_numba(X, max_bond, renorm) .. py:function:: _similarity_compress_svd(X, max_bond, renorm, asymm) .. py:function:: _similarity_compress_svd_numba(X, max_bond, renorm, asymm) .. py:function:: _similarity_compress_biorthog(X, max_bond, renorm) .. py:function:: _similarity_compress_biorthog_numba(X, max_bond, renorm) .. py:data:: _similarity_compress_fns .. py:function:: similarity_compress(X, max_bond, renorm=False, method='eigh') .. py:function:: isometrize_qr(x, backend=None) Perform isometrization using the QR decomposition. .. py:function:: isometrize_svd(x, backend=None) Perform isometrization using the SVD decomposition. .. py:function:: isometrize_exp(x, backend) Perform isometrization using anti-symmetric matrix exponentiation. .. math:: U_A = \exp \left( X - X^\dagger \right) If ``x`` is rectangular it is completed with zeros first. .. py:function:: isometrize_cayley(x, backend) Perform isometrization using an anti-symmetric Cayley transform. .. math:: U_A = (I + \dfrac{A}{2})(I - \dfrac{A}{2})^{-1} where :math:`A = X - X^\dagger`. If ``x`` is rectangular it is completed with zeros first. .. py:function:: isometrize_modified_gram_schmidt(A, backend=None) Perform isometrization explicitly using the modified Gram Schmidt procedure (this is slow but a useful reference). .. py:function:: isometrize_householder(X, backend=None) .. py:function:: isometrize_torch_householder(x) Isometrize ``x`` using the Householder reflection method, as implemented by the ``torch_householder`` package. .. py:data:: _ISOMETRIZE_METHODS .. py:function:: isometrize(x, method='qr') Generate an isometric (or unitary if square) / orthogonal matrix from array ``x``. :param x: The matrix to project into isometrix form. :type x: array :param method: The method used to generate the isometry. The options are: - "qr": use the Q factor of the QR decomposition of ``x`` with the constraint that the diagonal of ``R`` is positive. - "svd": uses ``U @ VH`` of the SVD decomposition of ``x``. This is useful for finding the 'closest' isometric matrix to ``x``, such as when it has been expanded with noise etc. But is less stable for differentiation / optimization. - "exp": use the matrix exponential of ``x - dag(x)``, first completing ``x`` with zeros if it is rectangular. This is a good parametrization for optimization, but more expensive for non-square ``x``. - "cayley": use the Cayley transform of ``x - dag(x)``, first completing ``x`` with zeros if it is rectangular. This is a good parametrization for optimization (one the few compatible with `HIPS/autograd` e.g.), but more expensive for non-square ``x``. - "householder": use the Householder reflection method directly. This requires that the backend implements "linalg.householder_product". - "torch_householder": use the Householder reflection method directly, using the ``torch_householder`` package. This requires that the package is installed and that the backend is ``"torch"``. This is generally the best parametrizing method for "torch" if available. - "mgs": use a python implementation of the modified Gram Schmidt method directly. This is slow if not compiled but a useful reference. Not all backends support all methods or differentiating through all methods. :type method: str, optional :returns: **Q** -- The isometrization / orthogonalization of ``x``. :rtype: array .. py:function:: squared_op_to_reduced_factor(x2, dl, dr, right=True) Given the square, ``x2``, of an operator ``x``, compute either the left or right reduced factor matrix of the unsquared operator ``x`` with original shape ``(dl, dr)``. .. py:function:: squared_op_to_reduced_factor_numba(x2, dl, dr, right=True) .. py:function:: compute_oblique_projectors(Rl, Rr, max_bond, cutoff, absorb='both', cutoff_mode=4, **compress_opts) Compute the oblique projectors for two reduced factor matrices that describe a gauge on a bond. Concretely, assuming that ``Rl`` and ``Rr`` are the reduced factor matrices for local operator ``A``, such that: .. math:: A = Q_L R_L R_R Q_R with ``Q_L`` and ``Q_R`` isometric matrices, then the optimal inner truncation is given by: .. math:: A' = Q_L P_L P_R' Q_R :param Rl: The left reduced factor matrix. :type Rl: array :param Rr: The right reduced factor matrix. :type Rr: array :returns: * **Pl** (*array*) -- The left oblique projector. * **Pr** (*array*) -- The right oblique projector.