quimb.linalg.approx_spectral ============================ .. py:module:: quimb.linalg.approx_spectral .. autoapi-nested-parse:: Use stochastic Lanczos quadrature to approximate spectral function sums of any operator which has an efficient representation of action on a vector. Functions --------- .. autoapisummary:: quimb.linalg.approx_spectral.lazy_ptr_linop quimb.linalg.approx_spectral.lazy_ptr_ppt_linop quimb.linalg.approx_spectral.inner quimb.linalg.approx_spectral.norm_fro quimb.linalg.approx_spectral.norm_fro_approx quimb.linalg.approx_spectral.random_rect quimb.linalg.approx_spectral.construct_lanczos_tridiag quimb.linalg.approx_spectral.lanczos_tridiag_eig quimb.linalg.approx_spectral.calc_trace_fn_tridiag quimb.linalg.approx_spectral.ext_per_trim quimb.linalg.approx_spectral.nbsum quimb.linalg.approx_spectral.std quimb.linalg.approx_spectral.calc_est_fit quimb.linalg.approx_spectral.calc_est_window quimb.linalg.approx_spectral.single_random_estimate quimb.linalg.approx_spectral.calc_stats quimb.linalg.approx_spectral.get_single_precision_dtype quimb.linalg.approx_spectral.get_equivalent_real_dtype quimb.linalg.approx_spectral.plot_approx_spectral_info quimb.linalg.approx_spectral.approx_spectral_function quimb.linalg.approx_spectral.tr_abs_approx quimb.linalg.approx_spectral.tr_exp_approx quimb.linalg.approx_spectral.tr_sqrt_approx quimb.linalg.approx_spectral.xlogx quimb.linalg.approx_spectral.tr_xlogx_approx quimb.linalg.approx_spectral.entropy_subsys_approx quimb.linalg.approx_spectral.tr_sqrt_subsys_approx quimb.linalg.approx_spectral.norm_ppt_subsys_approx quimb.linalg.approx_spectral.logneg_subsys_approx quimb.linalg.approx_spectral.negativity_subsys_approx quimb.linalg.approx_spectral.gen_bipartite_spectral_fn Module Contents --------------- .. py:function:: lazy_ptr_linop(psi_ab, dims, sysa, **linop_opts) A linear operator representing action of partially tracing a bipartite state, then multiplying another 'unipartite' state:: ( | ) +-------+ | psi_a | ______ +_______+ / \ a| |b | +-------------+ | | psi_ab.H | | +_____________+ | | +-------------+ | | psi_ab | | +_____________+ | a| |b | | \______/ :param psi_ab: State to partially trace and dot with another ket, with size ``prod(dims)``. :type psi_ab: ket :param dims: The sub dimensions of ``psi_ab``. :type dims: sequence of int, optional :param sysa: Index(es) of the 'a' subsystem(s) to keep. :type sysa: int or sequence of int, optional .. py:function:: lazy_ptr_ppt_linop(psi_abc, dims, sysa, sysb, **linop_opts) A linear operator representing action of partially tracing a tripartite state, partially transposing the remaining bipartite state, then multiplying another bipartite state:: ( | ) +--------------+ | psi_ab | +______________+ _____ a| ____ b| / \ | / a\ | |c | | | +-------------+ | | | | psi_abc.H | | \ / +-------------+ | X | / \ +-------------+ | | | | psi_abc | | | | +-------------+ | | \____/a |b |c | a| | \_____/ :param psi_abc: State to partially trace, partially transpose, then dot with another ket, with size ``prod(dims)``. ``prod(dims[sysa] + dims[sysb])``. :type psi_abc: ket :param dims: The sub dimensions of ``psi_abc``. :type dims: sequence of int :param sysa: Index(es) of the 'a' subsystem(s) to keep, with respect to all the dimensions, ``dims``, (i.e. pre-partial trace). :type sysa: int or sequence of int, optional :param sysa: Index(es) of the 'b' subsystem(s) to keep, with respect to all the dimensions, ``dims``, (i.e. pre-partial trace). :type sysa: int or sequence of int, optional .. py:function:: inner(a, b) Inner product between two vectors .. py:function:: norm_fro(a) 'Frobenius' norm of a vector. .. py:function:: norm_fro_approx(A, **kwargs) Calculate the approximate frobenius norm of any hermitian linear operator: .. math:: \mathrm{Tr} \left[ A^{\dagger} A \right] :param A: Operator with a dot method, assumed to be hermitian, to estimate the frobenius norm of. :type A: linear operator like :param kwargs: Supplied to :func:`approx_spectral_function`. :rtype: float .. py:function:: random_rect(shape, dist='rademacher', orthog=False, norm=True, seed=False, dtype=complex) Generate a random array optionally orthogonal. :param shape: The shape of array. :type shape: tuple of int :param dist: Distribution of the random variables. :type dist: {'guassian', 'rademacher'} :param orthog: Orthogonalize the columns if more than one. :type orthog: bool or operator. :param norm: Explicitly normalize the frobenius norm to 1. :type norm: bool .. py:function:: construct_lanczos_tridiag(A, K, v0=None, bsz=1, k_min=10, orthog=False, beta_tol=1e-06, seed=False, v0_opts=None) Construct the tridiagonal lanczos matrix using only matvec operators. This is a generator that iteratively yields the alpha and beta digaonals at each step. :param A: The operator to approximate, must implement ``.dot`` method to compute its action on a vector. :type A: dense array, sparse matrix or linear operator :param K: The maximum number of iterations and thus rank of the matrix to find. :type K: int, optional :param v0: The starting vector to iterate with, default to random. :type v0: vector, optional :param bsz: The block size (number of columns) of random vectors to iterate with. :type bsz: int, optional :param k_min: The minimum size of the krylov subspace for form. :type k_min: int, optional :param orthog: If True, perform full re-orthogonalization for each new vector. :type orthog: bool, optional :param beta_tol: The 'breakdown' tolerance. If the next beta ceofficient in the lanczos matrix is less that this, implying that the full non-null space has been found, terminate early. :type beta_tol: float, optional :param seed: If True, seed the numpy random generator with a system random int. :type seed: bool, optional :Yields: * **alpha** (*sequence of float of length k*) -- The diagonal entries of the lanczos matrix. * **beta** (*sequence of float of length k*) -- The off-diagonal entries of the lanczos matrix, with the last entry the 'look' forward value. * **scaling** (*float*) -- How to scale the overall weights. .. py:function:: lanczos_tridiag_eig(alpha, beta, check_finite=True) Find the eigen-values and -vectors of the Lanczos triadiagonal matrix. :param alpha: The diagonal. :type alpha: array of float :param beta: The k={-1, 1} off-diagonal. Only first ``len(alpha) - 1`` entries used. :type beta: array of float .. py:function:: calc_trace_fn_tridiag(tl, tv, f, pos=True) Spectral ritz function sum, weighted by ritz vectors. .. py:function:: ext_per_trim(x, p=0.6, s=1.0) Extended percentile trimmed-mean. Makes the mean robust to asymmetric outliers, while using all data when it is nicely clustered. This can be visualized roughly as:: |--------|=========|--------| x x xx xx xxxxx xxx xx x x x Where the inner range contains the central ``p`` proportion of the data, and the outer ranges entends this by a factor of ``s`` either side. :param x: Data to trim. :type x: array :param p: For example, p=0.5 gives the inter-quartile range. :type p: Proportion of data used to define the 'central' percentile. :param s: away from the central percentile itself. :type s: Include data up to this factor times the central 'percentile' range :param Returns: :param xt: Trimmed data. :type xt: array .. py:function:: nbsum(xs) .. py:function:: std(xs) Simple standard deviation - don't invoke numpy for small lists. .. py:function:: calc_est_fit(estimates, conv_n, tau) Make estimate by fitting exponential convergence to estimates. .. py:function:: calc_est_window(estimates, conv_n) Make estimate from mean of last ``m`` samples, following: 1. Take between ``conv_n`` and 12 estimates. 2. Pair the estimates as they are alternate upper/lower bounds 3. Compute the standard error on the paired estimates. .. py:function:: single_random_estimate(A, K, bsz, beta_tol, v0, f, pos, tau, tol_scale, k_min=10, verbosity=0, *, seed=None, v0_opts=None, info=None, **lanczos_opts) .. py:function:: calc_stats(samples, mean_p, mean_s, tol, tol_scale) Get an estimate from samples. .. py:function:: get_single_precision_dtype(dtype) .. py:function:: get_equivalent_real_dtype(dtype) .. py:function:: plot_approx_spectral_info(info) .. py:function:: approx_spectral_function(A, f, tol=0.01, *, bsz=1, R=1024, R_min=3, tol_scale=1, tau=0.0001, k_min=10, k_max=512, beta_tol=1e-06, mpi=False, mean_p=0.7, mean_s=1.0, pos=False, v0=None, verbosity=0, single_precision='AUTO', info=None, progbar=False, plot=False, **lanczos_opts) Approximate a spectral function, that is, the quantity ``Tr(f(A))``. :param A: Operator to approximate spectral function for. Should implement ``A.dot(vec)``. :type A: dense array, sparse matrix or LinearOperator :param f: Scalar function with which to act on approximate eigenvalues. :type f: callable :param tol: Relative convergence tolerance threshold for error on mean of repeats. This can pretty much be relied on as the overall accuracy. See also ``tol_scale`` and ``tau``. Default: 1%. :type tol: float, optional :param bsz: Number of simultenous vector columns to use at once, 1 equating to the standard lanczos method. If ``bsz > 1`` then ``A`` must implement matrix-matrix multiplication. This is a more performant way of essentially increasing ``R``, at the cost of more memory. Default: 1. :type bsz: int, optional :param R: The number of repeats with different initial random vectors to perform. Increasing this should increase accuracy as ``sqrt(R)``. Cost of algorithm thus scales linearly with ``R``. If ``tol`` is non-zero, this is the maximum number of repeats. :type R: int, optional :param R_min: The minimum number of repeats to perform. Default: 3. :type R_min: int, optional :param tau: The relative tolerance required for a single lanczos run to converge. This needs to be small enough that each estimate with a single random vector produces an unbiased sample of the operators spectrum.. :type tau: float, optional :param k_min: The minimum size of the krylov subspace to form for each sample. :type k_min: int, optional :param k_max: The maximum size of the kyrlov space to form. Cost of algorithm scales linearly with ``K``. If ``tau`` is non-zero, this is the maximum size matrix to form. :type k_max: int, optional :param tol_scale: This sets the overall expected scale of each estimate, so that an absolute tolerance can be used for values near zero. Default: 1. :type tol_scale: float, optional :param beta_tol: The 'breakdown' tolerance. If the next beta ceofficient in the lanczos matrix is less that this, implying that the full non-null space has been found, terminate early. Default: 1e-6. :type beta_tol: float, optional :param mpi: Whether to parallelize repeat runs over MPI processes. :type mpi: bool, optional :param mean_p: Factor for robustly finding mean and err of repeat estimates, see :func:`ext_per_trim`. :type mean_p: float, optional :param mean_s: Factor for robustly finding mean and err of repeat estimates, see :func:`ext_per_trim`. :type mean_s: float, optional :param v0: Initial vector to iterate with, sets ``R=1`` if given. If callable, the function to produce a random intial vector (sequence). :type v0: vector, or callable :param pos: If True, make sure any approximate eigenvalues are positive by clipping below 0. :type pos: bool, optional :param verbosity: How much information to print while computing. :type verbosity: {0, 1, 2}, optional :param single_precision: Try and convert the operator to single precision. This can lead to much faster operation, especially if a GPU is available. Additionally, double precision is not really needed given the stochastic nature of the algorithm. :type single_precision: {'AUTO', False, True}, optional :param lanczos_opts: Supplied to :func:`~quimb.linalg.approx_spectral.single_random_estimate` or :func:`~quimb.linalg.approx_spectral.construct_lanczos_tridiag`. :returns: The approximate value ``Tr(f(a))``. :rtype: scalar .. seealso:: :py:obj:`construct_lanczos_tridiag` .. py:function:: tr_abs_approx(*args, **kwargs) .. py:function:: tr_exp_approx(*args, **kwargs) .. py:function:: tr_sqrt_approx(*args, **kwargs) .. py:function:: xlogx(x) .. py:function:: tr_xlogx_approx(*args, **kwargs) .. py:function:: entropy_subsys_approx(psi_ab, dims, sysa, backend=None, **kwargs) Approximate the (Von Neumann) entropy of a pure state's subsystem. :param psi_ab: Bipartite state to partially trace and find entopy of. :type psi_ab: ket :param dims: The sub dimensions of ``psi_ab``. :type dims: sequence of int, optional :param sysa: Index(es) of the 'a' subsystem(s) to keep. :type sysa: int or sequence of int, optional :param kwargs: Supplied to :func:`approx_spectral_function`. .. py:function:: tr_sqrt_subsys_approx(psi_ab, dims, sysa, backend=None, **kwargs) Approximate the trace sqrt of a pure state's subsystem. :param psi_ab: Bipartite state to partially trace and find trace sqrt of. :type psi_ab: ket :param dims: The sub dimensions of ``psi_ab``. :type dims: sequence of int, optional :param sysa: Index(es) of the 'a' subsystem(s) to keep. :type sysa: int or sequence of int, optional :param kwargs: Supplied to :func:`approx_spectral_function`. .. py:function:: norm_ppt_subsys_approx(psi_abc, dims, sysa, sysb, backend=None, **kwargs) Estimate the norm of the partial transpose of a pure state's subsystem. .. py:function:: logneg_subsys_approx(psi_abc, dims, sysa, sysb, **kwargs) Estimate the logarithmic negativity of a pure state's subsystem. :param psi_abc: Pure tripartite state, for which estimate the entanglement between 'a' and 'b'. :type psi_abc: ket :param dims: The sub dimensions of ``psi_abc``. :type dims: sequence of int :param sysa: Index(es) of the 'a' subsystem(s) to keep, with respect to all the dimensions, ``dims``, (i.e. pre-partial trace). :type sysa: int or sequence of int, optional :param sysa: Index(es) of the 'b' subsystem(s) to keep, with respect to all the dimensions, ``dims``, (i.e. pre-partial trace). :type sysa: int or sequence of int, optional :param kwargs: Supplied to :func:`approx_spectral_function`. .. py:function:: negativity_subsys_approx(psi_abc, dims, sysa, sysb, **kwargs) Estimate the negativity of a pure state's subsystem. :param psi_abc: Pure tripartite state, for which estimate the entanglement between 'a' and 'b'. :type psi_abc: ket :param dims: The sub dimensions of ``psi_abc``. :type dims: sequence of int :param sysa: Index(es) of the 'a' subsystem(s) to keep, with respect to all the dimensions, ``dims``, (i.e. pre-partial trace). :type sysa: int or sequence of int, optional :param sysa: Index(es) of the 'b' subsystem(s) to keep, with respect to all the dimensions, ``dims``, (i.e. pre-partial trace). :type sysa: int or sequence of int, optional :param kwargs: Supplied to :func:`approx_spectral_function`. .. py:function:: gen_bipartite_spectral_fn(exact_fn, approx_fn, pure_default) Generate a function that computes a spectral quantity of the subsystem of a pure state. Automatically computes for the smaller subsystem, or switches to the approximate method for large subsystems. :param exact_fn: The function that computes the quantity on a density matrix, with signature: ``exact_fn(rho_a, rank=...)``. :type exact_fn: callable :param approx_fn: The function that approximately computes the quantity using a lazy representation of the whole system. With signature ``approx_fn(psi_ab, dims, sysa, **approx_opts)``. :type approx_fn: callable :param pure_default: The default value when the whole state is the subsystem. :type pure_default: float :returns: **bipartite_spectral_fn** -- The function, with signature: ``(psi_ab, dims, sysa, approx_thresh=2**13, **approx_opts)`` :rtype: callable