symmray.linalg ============== .. py:module:: symmray.linalg .. autoapi-nested-parse:: Interface to linear algebra submodule functions. Functions --------- .. autoapisummary:: symmray.linalg.eigh symmray.linalg.eigh_truncated symmray.linalg.norm symmray.linalg.qr symmray.linalg.qr_stabilized symmray.linalg.lq symmray.linalg.lq_stabilized symmray.linalg.solve symmray.linalg.svd symmray.linalg.svd_truncated symmray.linalg.svd_rand_truncated symmray.linalg.svd_via_eig_truncated symmray.linalg.cholesky symmray.linalg.cholesky_regularized symmray.linalg.lq_via_cholesky symmray.linalg.qr_via_cholesky Module Contents --------------- .. py:function:: eigh(x: symmray.array_common.ArrayCommon, *args, **kwargs) Hermitian eigen-decomposition of an assumed hermitian symmray array. :returns: * **w** (*VectorCommon*) -- The eigenvalues as a vector. * **u** (*ArrayCommon*) -- The array of eigenvectors. .. py:function:: eigh_truncated(x: symmray.array_common.ArrayCommon, *args, **kwargs) Truncated hermitian eigen-decomposition of an assumed hermitian symmray array. :param cutoff: Absolute eigenvalue cutoff threshold. :type cutoff: float, optional :param cutoff_mode: How to perform the truncation: - 1 or 'abs': trim values below ``cutoff`` - 2 or 'rel': trim values below ``s[0] * cutoff`` - 3 or 'sum2': trim s.t. ``sum(s_trim**2) < cutoff``. - 4 or 'rsum2': trim s.t. ``sum(s_trim**2) < sum(s**2) * cutoff``. - 5 or 'sum1': trim s.t. ``sum(s_trim**1) < cutoff``. - 6 or 'rsum1': trim s.t. ``sum(s_trim**1) < sum(s**1) * cutoff``. :type cutoff_mode: int or str, optional :param max_bond: An explicit maximum bond dimension, use -1 for none. :type max_bond: int :param absorb: How to absorb the eigenvalues. - -1 or 'left': absorb into the left factor (U). - 0 or 'both': absorb the square root into both factors. - 1 or 'right': absorb into the right factor (VH). - None: do not absorb, return eigenvalues as a BlockVector. :type absorb: {-1, 0, 1, None} :param renorm: Whether to renormalize the eigenvalues (depends on `cutoff_mode`). :type renorm: {0, 1} :returns: * **u** (*ArrayCommon*) -- The abelian array of left eigenvectors. * **w** (*VectorCommon or None*) -- The vector of eigenvalues, or None if absorbed. * **uh** (*ArrayCommon*) -- The abelian array of right eigenvectors. .. py:function:: norm(x: symmray.array_common.ArrayCommon, *args, **kwargs) .. py:function:: qr(x: symmray.array_common.ArrayCommon, *args, **kwargs) QR decomposition of a symmray array. .. py:function:: qr_stabilized(x: symmray.array_common.ArrayCommon, *args, **kwargs) Stabilized QR decomposition of a symmray array, returning results in an SVD-like ``(Q, None, R)`` format for compatibility with tensor network split drivers. .. py:function:: lq(x: symmray.array_common.ArrayCommon, *args, **kwargs) LQ decomposition of a symmray array. .. py:function:: lq_stabilized(x: symmray.array_common.ArrayCommon, *args, **kwargs) Stabilized LQ decomposition of a symmray array, returning results in an SVD-like ``(L, None, Q)`` format for compatibility with tensor network split drivers. .. py:function:: solve(x: symmray.array_common.ArrayCommon, *args, **kwargs) .. py:function:: svd(x: symmray.array_common.ArrayCommon, *args, **kwargs) -> tuple[symmray.array_common.ArrayCommon, symmray.vector_common.VectorCommon, symmray.array_common.ArrayCommon] Singular value decomposition of a symmray array, returning all singular values and vectors without truncation. .. py:function:: svd_truncated(x: symmray.array_common.ArrayCommon, *args, **kwargs) Truncated singular value decomposition of a symmray array. :param cutoff: Singular value cutoff threshold. :type cutoff: float, optional :param cutoff_mode: How to perform the truncation: - 1 or 'abs': trim values below ``cutoff`` - 2 or 'rel': trim values below ``s[0] * cutoff`` - 3 or 'sum2': trim s.t. ``sum(s_trim**2) < cutoff``. - 4 or 'rsum2': trim s.t. ``sum(s_trim**2) < sum(s**2) * cutoff``. - 5 or 'sum1': trim s.t. ``sum(s_trim**1) < cutoff``. - 6 or 'rsum1': trim s.t. ``sum(s_trim**1) < sum(s**1) * cutoff``. :type cutoff_mode: int or str, optional :param max_bond: An explicit maximum bond dimension, use -1 for none. :type max_bond: int :param absorb: How to absorb the singular values. - -1 or 'left': absorb into the left factor (U). - 0 or 'both': absorb the square root into both factors. - 1 or 'right': absorb into the right factor (VH). - None: do not absorb, return singular values as a BlockVector. :type absorb: {-1, 0, 1, None} :param renorm: Whether to renormalize the singular values (depends on `cutoff_mode`). :type renorm: {0, 1} :returns: * **u** (*ArrayCommon*) -- The abelian array of left singular vectors. * **s** (*VectorCommon or None*) -- The vector of singular values, or None if absorbed. * **vh** (*ArrayCommon*) -- The abelian array of right singular vectors. .. py:function:: svd_rand_truncated(x: symmray.array_common.ArrayCommon, *args, **kwargs) Truncated singular value decomposition of a symmray array, using randomized projection. This is efficient for low-rank approximations when a target ``max_bond`` is known. :param max_bond: Target rank / maximum bond dimension. :type max_bond: int :param absorb: How to absorb the singular values. - -1 or 'left': absorb into the left factor (U). - 0 or 'both': absorb the square root into both factors. - 1 or 'right': absorb into the right factor (VH). - None: do not absorb, return singular values. :type absorb: {-1, 0, 1, None} :param oversample: Extra sketch dimensions for accuracy. Default is 10. :type oversample: int, optional :param num_iterations: Number of power iterations for accuracy. Default is 2. :type num_iterations: int, optional :param seed: Random seed or generator for reproducibility. :type seed: int, Generator or None, optional :returns: * **u** (*ArrayCommon or None*) -- The array of left singular vectors. * **s** (*VectorCommon or None*) -- The singular values, or None if absorbed. * **vh** (*ArrayCommon or None*) -- The array of right singular vectors. .. py:function:: svd_via_eig_truncated(x: symmray.array_common.ArrayCommon, *args, **kwargs) Truncated singular value decomposition of a symmray array, using eigen-decomposition of the gram (xdag @ x or x @ xdag) matrix. This can be faster, but also incur a loss of precision due to the squaring. :param cutoff: Singular value cutoff threshold. :type cutoff: float, optional :param cutoff_mode: How to perform the truncation: - 1 or 'abs': trim values below ``cutoff`` - 2 or 'rel': trim values below ``s[0] * cutoff`` - 3 or 'sum2': trim s.t. ``sum(s_trim**2) < cutoff``. - 4 or 'rsum2': trim s.t. ``sum(s_trim**2) < sum(s**2) * cutoff``. - 5 or 'sum1': trim s.t. ``sum(s_trim**1) < cutoff``. - 6 or 'rsum1': trim s.t. ``sum(s_trim**1) < sum(s**1) * cutoff``. :type cutoff_mode: int or str, optional :param max_bond: An explicit maximum bond dimension, use -1 for none. :type max_bond: int :param absorb: How to absorb the singular values. - -1 or 'left': absorb into the left factor (U). - 0 or 'both': absorb the square root into both factors. - 1 or 'right': absorb into the right factor (VH). - None: do not absorb, return singular values as a BlockVector. :type absorb: {-1, 0, 1, None} :param renorm: Whether to renormalize the singular values (depends on `cutoff_mode`). :type renorm: {0, 1} :returns: * **u** (*ArrayCommon*) -- The abelian array of left singular vectors. * **s** (*VectorCommon or None*) -- The vector of singular values, or None if absorbed. * **vh** (*ArrayCommon*) -- The abelian array of right singular vectors. .. py:function:: cholesky(x: symmray.array_common.ArrayCommon, *args, **kwargs) Cholesky decomposition of an assumed positive-definite symmray array. :param x: The 2D block-symmetric array to decompose. :type x: ArrayCommon :param upper: Whether to return the upper triangular Cholesky factor. Default is False, returning the lower triangular factor. :type upper: bool, optional :returns: **l_or_r** -- The Cholesky factor. Lower triangular if ``upper=False``, upper triangular if ``upper=True``. :rtype: ArrayCommon .. py:function:: cholesky_regularized(x: symmray.array_common.ArrayCommon, *args, **kwargs) Cholesky decomposition with optional diagonal regularization, returning results in an SVD-like ``(left, None, right)`` format for compatibility with tensor network split drivers. :param x: The 2D block-symmetric array to decompose. Must be positive (semi-)definite. :type x: ArrayCommon :param absorb: How to return the factors: - ``0`` (``'both'``): return ``(L, None, L^H)``. - ``-12`` (``'lsqrt'``): return ``(L, None, None)``. - ``12`` (``'rsqrt'``): return ``(None, None, L^H)``. :type absorb: {-12, 0, 12}, optional :param shift: Diagonal regularization shift. If True or negative, auto-compute from dtype machine epsilon. The shift is always applied as a relative shift scaled by the trace of each block. Default is True. :type shift: float, optional :returns: * **left** (*ArrayCommon or None*) -- The lower Cholesky factor, or None. * **s** (*None*) -- Always None (no singular values). * **right** (*ArrayCommon or None*) -- The conjugate transpose of the Cholesky factor, or None. .. py:function:: lq_via_cholesky(x: symmray.array_common.ArrayCommon, *args, **kwargs) LQ decomposition via Cholesky factorization of ``x @ x^H``. Computes ``x = L @ Q`` where ``L`` is lower triangular and ``Q`` is isometric. :param x: The 2D block-symmetric array to decompose. :type x: ArrayCommon :param absorb: How to return the factors: - ``-1`` (``'left'``): return ``(L, None, Q)``. - ``-10`` (``'lfactor'``): return ``(L, None, None)``. - ``-11`` (``'rorthog'``): return ``(None, None, Q)``. :type absorb: {-1, -10, -11} or str, optional :param shift: Diagonal regularization shift. If True or negative, auto-compute from dtype machine epsilon. The shift is always applied as a relative shift scaled by the trace of each block. Default is True. :type shift: float, optional :param solve_triangular: Whether to use triangular solve (faster) or general solve to compute Q. Default is True. :type solve_triangular: bool, optional :returns: * **L** (*ArrayCommon or None*) -- The lower triangular factor. * **s** (*None*) -- Always None. * **Q** (*ArrayCommon or None*) -- The isometric factor. .. py:function:: qr_via_cholesky(x: symmray.array_common.ArrayCommon, *args, **kwargs) QR decomposition via Cholesky factorization of ``x^H @ x``. Computes ``x = Q @ R`` where ``Q`` is isometric and ``R`` is upper triangular. :param x: The 2D block-symmetric array to decompose. :type x: ArrayCommon :param absorb: How to return the factors: - ``1`` (``'right'``): return ``(Q, None, R)``. - ``11`` (``'rfactor'``): return ``(None, None, R)``. - ``10`` (``'lorthog'``): return ``(Q, None, None)``. :type absorb: {1, 11, 10} or str, optional :param shift: Diagonal regularization shift. If True or negative, auto-compute from dtype machine epsilon. The shift is always applied as a relative shift scaled by the trace of each block. Default is True. :type shift: float, optional :param solve_triangular: Whether to use triangular solve (faster) or general solve. Default is True. :type solve_triangular: bool, optional :returns: * **Q** (*ArrayCommon or None*) -- The isometric factor. * **s** (*None*) -- Always None. * **R** (*ArrayCommon or None*) -- The upper triangular factor.