quimb.calc ========== .. py:module:: quimb.calc .. autoapi-nested-parse:: Functions for more advanced calculations of quantities and properties of quantum objects. Attributes ---------- .. autoapisummary:: quimb.calc.entropy_subsys quimb.calc.mutual_information quimb.calc.tr_sqrt_subsys quimb.calc.logarithmic_negativity quimb.calc.pauli_decomp quimb.calc.bell_decomp Functions --------- .. autoapisummary:: quimb.calc.fidelity quimb.calc.purify quimb.calc.kraus_op quimb.calc.projector quimb.calc.measure quimb.calc.simulate_counts quimb.calc.dephase quimb.calc.entropy quimb.calc.mutinf quimb.calc.check_dims_and_indices quimb.calc.mutinf_subsys quimb.calc.schmidt_gap quimb.calc.tr_sqrt quimb.calc.partial_transpose quimb.calc.partial_transpose_norm quimb.calc.logneg quimb.calc.logneg_subsys quimb.calc.negativity quimb.calc.concurrence quimb.calc.one_way_classical_information quimb.calc.quantum_discord quimb.calc.trace_distance quimb.calc.cprint quimb.calc.decomp quimb.calc.correlation quimb.calc.pauli_correlations quimb.calc.ent_cross_matrix quimb.calc.qid quimb.calc.is_degenerate quimb.calc.is_eigenvector quimb.calc.page_entropy quimb.calc.heisenberg_energy Module Contents --------------- .. py:function:: fidelity(p1, p2, squared=False) Fidelity between two quantum states. By default, the unsquared fidelity is used, equivalent in the pure state case to ``||``. :param p1: First state. :type p1: vector or operator :param p2: Second state. :type p2: vector or operator :param squared: Whether to use the squared convention or not, by default False. :type squared: bool, optional :rtype: float .. py:function:: purify(rho) Take state rho and purify it into a wavefunction of squared dimension. :param rho: Density operator to purify. :type rho: operator :returns: The purified ket. :rtype: vector .. py:function:: kraus_op(rho, Ek, dims=None, where=None, check=False) Operate on a mixed state with a set of kraus operators: .. math:: \sigma = \sum_k E_k \rho E_k^{\dagger} :param rho: The density operator to perform operation on. :type rho: (d, d) density matrix :param Ek: The Kraus operator(s). :type Ek: (K, d, d) array or sequence of K (d, d) arrays :param dims: The subdimensions of ``rho``. :type dims: sequence of int, optional :param where: Which of the subsystems to apply the operation on. :type where: int or sequence of int, optional :param check: Where to check ``sum_k(Ek.H @ Ek) == 1``. :type check: bool, optional :returns: **sigma** -- The state after the kraus operation. :rtype: density matrix .. rubric:: Examples The depolarising channel: .. code-block:: python3 In [1]: import quimb as qu In [2]: rho = qu.rand_rho(2) In [3]: I, X, Y, Z = (qu.pauli(s) for s in 'IXYZ') In [4]: [qu.expec(rho, A) for A in (X, Y, Z)] Out[4]: [-0.652176185230622, -0.1301762132792484, 0.22362918368272583] In [5]: p = 0.1 In [6]: Ek = [ ...: (1 - p)**0.5 * I, ...: (p / 3)**0.5 * X, ...: (p / 3)**0.5 * Y, ...: (p / 3)**0.5 * Z ...: ] In [7]: sigma = qu.kraus_op(rho, Ek) In [8]: [qu.expec(sigma, A) for A in (X, Y, Z)] Out[8]: [-0.5652193605332058, -0.11281938484201527, 0.1938119591916957] .. py:function:: projector(A, eigenvalue=1.0, tol=1e-12, autoblock=False) Compute the projector for the target ``eigenvalue`` of operator ``A``. :param A: The hermitian observable. If a tuple is supplied, assume this is the eigendecomposition of the observable. :type A: operator or tuple[array, operator]. :param eigenvalue: The target eigenvalue to construct the projector for, default: 1. :type eigenvalue: float, optional :param tol: The tolerance with which to select all eigenvalues. :type tol: float, optional :param autoblock: Whether to use automatic symmetry exploitation. :type autoblock: bool, optional :returns: **P** -- The projector onto the target eigenspace. :rtype: dense matrix .. py:function:: measure(p, A, eigenvalue=None, tol=1e-12) Measure state ``p`` with observable ``A``, collapsing the state in the process and returning the relevant eigenvalue. .. math:: \left| \psi \right\rangle \rightarrow \frac{ P_i \left| \psi \right\rangle }{ \sqrt{\left\langle \psi \right| P_i \left| \psi \right\rangle} } and .. math:: \rho \rightarrow \frac{ P_i \rho P_i^{\dagger} }{ \text{Tr} \left[ P_i \rho \right] } along with the corresponding eigenvalue :math:`\lambda_i`. :param p: The quantum state to measure. :type p: vector or matrix :param A: The hermitian operator corresponding to an observable. You can supply also a pre-diagonalized operator here as a tuple of eigenvalues and eigenvectors. :type A: matrix or tuple[array, matrix] :param tol: The tolerance within which to group eigenspaces. :type tol: float, optional :param eigenvalue: If specified, deterministically collapse to this result. Otherwise randomly choose a result as in 'real-life'. :type eigenvalue: float, optional :returns: * **result** (*float*) -- The result of the measurement. * **p_after** (*vector or matrix*) -- The quantum state post measurement. .. py:function:: simulate_counts(p, C, phys_dim=2, seed=None) Simulate measuring each qubit of ``p`` in the computational basis, producing output like that of ``qiskit``. :param p: The quantum state, assumed to be normalized, as either a ket or density operator. :type p: vector or operator :param C: The number of counts to perform. :type C: int :param phys_dim: The assumed size of the subsystems of ``p``, defaults to 2 for qubits. :type phys_dim: int, optional :returns: **results** -- The counts for each bit string measured. :rtype: dict[str, int] .. rubric:: Examples Simulate measuring the state of each qubit in a GHZ-state: .. code:: python3 >>> import quimb as qu >>> psi = qu.ghz_state(3) >>> qu.simulate_counts(psi, 1024) {'000': 514, '111': 510} .. py:function:: dephase(rho, p, rand_rank=None) Dephase ``rho`` by amount ``p``, that is, mix it with the maximally mixed state: rho -> (1 - p) * rho + p * I / d :param rho: The state. :type rho: operator :param p: The final proportion of identity. :type p: float :param rand_rank: If given, dephase with a random diagonal operator with this many non-zero entries. If float, proportion of full size. :type rand_rank: int or float, optional :returns: **rho_dephase** -- The dephased density operator. :rtype: operator .. py:function:: entropy(a, rank=None) Compute the (von Neumann) entropy. :param a: Positive operator or list of positive eigenvalues. :type a: operator or 1d array :param rank: If operator has known rank, then a partial decomposition can be used to accelerate the calculation. :type rank: int (optional) :returns: The von Neumann entropy. :rtype: float .. seealso:: :py:obj:`mutinf`, :py:obj:`entropy_subsys`, :py:obj:`entropy_subsys_approx` .. py:data:: entropy_subsys Calculate the entropy of a pure states' subsystem, optionally switching to an approximate lanczos method when the subsystem is very large. :param psi_ab: Bipartite state. :type psi_ab: vector :param dims: The sub-dimensions of the state. :type dims: sequence of int :param sysa: The indices of which dimensions to calculate the entropy for. :type sysa: sequence of int :param approx_thresh: The size of sysa at which to switch to the approx method. Set to ``None`` to never use the approximation. :type approx_thresh: int, optional :param \*\*approx_opts: Supplied to :func:`entropy_subsys_approx`, if used. :returns: The subsytem entropy. :rtype: float .. seealso:: :py:obj:`entropy`, :py:obj:`entropy_subsys_approx`, :py:obj:`mutinf_subsys` .. py:function:: mutinf(p, dims=(2, 2), sysa=0, rank=None) Find the mutual information for a bipartition of a state. That is, ``H(A) + H(B) - H(AB)``, for von Neumann entropy ``H``, and two subsystems A and B. :param p: State, can be vector or operator. :type p: vector or operator :param dims: Internal dimensions of state. :type dims: tuple(int), optional :param sysa: Index of first subsystem, A. :type sysa: int, optional :param sysb: Index of second subsystem, B. :type sysb: int, optional :param rank: If known, the rank of rho_ab, to speed calculation of ``H(AB)`` up. For example, if ``p`` comes from tracing out three qubits from a system, then its rank is 2^3 = 8 etc. :type rank: int, optional :rtype: float .. seealso:: :py:obj:`entropy`, :py:obj:`mutinf_subsys`, :py:obj:`entropy_subsys_approx` .. py:data:: mutual_information .. py:function:: check_dims_and_indices(dims, *syss) Make sure all indices found in the tuples ``syss`` are in ``range(len(dims))``. .. py:function:: mutinf_subsys(psi_abc, dims, sysa, sysb, approx_thresh=2**13, **approx_opts) Calculate the mutual information of two subsystems of a pure state, possibly using an approximate lanczos method for large subsytems. :param psi_abc: Tri-partite pure state. :type psi_abc: vector :param dims: The sub dimensions of the state. :type dims: sequence of int :param sysa: The index(es) of the subsystem(s) to consider part of 'A'. :type sysa: sequence of int :param sysb: The index(es) of the subsystem(s) to consider part of 'B'. :type sysb: sequence of int :param approx_thresh: The size of subsystem at which to switch to the approximate lanczos method. Set to ``None`` to never use the approximation. :type approx_thresh: int, optional :param approx_opts: Supplied to :func:`entropy_subsys_approx`, if used. :returns: The mutual information. :rtype: float .. seealso:: :py:obj:`mutinf`, :py:obj:`entropy_subsys`, :py:obj:`entropy_subsys_approx`, :py:obj:`logneg_subsys` .. py:function:: schmidt_gap(psi_ab, dims, sysa) Find the schmidt gap of the bipartition of ``psi_ab``. That is, the difference between the two largest eigenvalues of the reduced density operator. :param psi_ab: Bipartite state. :type psi_ab: vector :param dims: The sub-dimensions of the state. :type dims: sequence of int :param sysa: The indices of which dimensions to calculate the entropy for. :type sysa: sequence of int :rtype: float .. py:function:: tr_sqrt(A, rank=None) Return the trace of the sqrt of a positive semidefinite operator. .. py:data:: tr_sqrt_subsys Compute the trace sqrt of a subsystem, possibly using an approximate lanczos method when the subsytem is big. :param psi_ab: Bipartite state. :type psi_ab: vector :param dims: The sub-dimensions of the state. :type dims: sequence of int :param sysa: The indices of which dimensions to calculate the trace sqrt for. :type sysa: sequence of int :param approx_thresh: The size of sysa at which to switch to the approx method. Set to ``None`` to never use the approximation. :type approx_thresh: int, optional :param \*\*approx_opts: Supplied to :func:`tr_sqrt_subsys_approx`, if used. :returns: The subsytem entropy. :rtype: float .. seealso:: :py:obj:`tr_sqrt`, :py:obj:`tr_sqrt_subsys_approx`, :py:obj:`partial_transpose_norm` .. py:function:: partial_transpose(p, dims=(2, 2), sysa=0) Partial transpose of a density operator. :param p: The state to partially transpose. :type p: operator or vector :param dims: The internal dimensions of the state. :type dims: tuple(int), optional :param sysa: The indices of 'system A', everything else assumed to be 'system B'. :type sysa: sequence of int :rtype: operator .. seealso:: :py:obj:`logneg`, :py:obj:`negativity` .. py:function:: partial_transpose_norm(p, dims, sysa) Compute the norm of the partial transpose for (log)-negativity, taking a shortcut (trace sqrt of reduced subsytem), when system is a vector. .. py:function:: logneg(p, dims=(2, 2), sysa=0) Compute logarithmic negativity between two subsytems. This is defined as log_2( | rho_{AB}^{T_B} | ). This only handles bipartitions (including pure states efficiently), and will not trace anything out. :param p: State to compute logarithmic negativity for. :type p: ket vector or density operator :param dims: The internal dimensions of ``p``. :type dims: tuple(int), optional :param sysa: Index of the first subsystem, A, relative to ``dims``. :type sysa: int, optional :rtype: float .. seealso:: :py:obj:`negativity`, :py:obj:`partial_transpose`, :py:obj:`logneg_subsys_approx` .. py:data:: logarithmic_negativity .. py:function:: logneg_subsys(psi_abc, dims, sysa, sysb, approx_thresh=2**13, **approx_opts) Compute the logarithmic negativity between two subsystems of a pure state, possibly using an approximate lanczos for large subsystems. Uses a special method if the two subsystems form a bipartition of the state. :param psi_abc: Tri-partite pure state. :type psi_abc: vector :param dims: The sub dimensions of the state. :type dims: sequence of int :param sysa: The index(es) of the subsystem(s) to consider part of 'A'. :type sysa: sequence of int :param sysb: The index(es) of the subsystem(s) to consider part of 'B'. :type sysb: sequence of int :param approx_thresh: The size of subsystem at which to switch to the approximate lanczos method. Set to ``None`` to never use the approximation. :type approx_thresh: int, optional :param approx_opts: Supplied to :func:`~quimb.logneg_subsys_approx`, if used. :returns: The logarithmic negativity. :rtype: float .. seealso:: :py:obj:`logneg`, :py:obj:`mutinf_subsys`, :py:obj:`logneg_subsys_approx` .. py:function:: negativity(p, dims=(2, 2), sysa=0) Compute negativity between two subsytems. This is defined as (| rho_{AB}^{T_B} | - 1) / 2. If ``len(dims) > 2``, then the non-target dimensions will be traced out first. :param p: State to compute logarithmic negativity for. :type p: ket vector or density operator :param dims: The internal dimensions of ``p``. :type dims: tuple(int), optional :param sysa: Index of the first subsystem, A, relative to ``dims``. :type sysa: int, optional :rtype: float .. seealso:: :py:obj:`logneg`, :py:obj:`partial_transpose`, :py:obj:`negativity_subsys_approx` .. py:function:: concurrence(p, dims=(2, 2), sysa=0, sysb=1) Concurrence of two-qubit state. If ``len(dims) > 2``, then the non-target dimensions will be traced out first. :param p: State to compute concurrence for. :type p: ket vector or density operator :param dims: The internal dimensions of ``p``. :type dims: tuple(int), optional :param sysa: Index of the first subsystem, A, relative to ``dims``. :type sysa: int, optional :param sysb: Index of the first subsystem, B, relative to ``dims``. :type sysb: int, optional :rtype: float .. py:function:: one_way_classical_information(p_ab, prjs, precomp_func=False) One way classical information for two qubit density operator. :param p_ab: State of two qubits :type p_ab: operator :param prjs: The POVMs. :type prjs: sequence of matrices :param precomp_func: Whether to return a pre-computed function, closed over the actual state. :type precomp_func: bool, optional :returns: The one-way classical information or the function to compute it for the given state which takes a set of POVMs as its single argument. :rtype: float or callable .. py:function:: quantum_discord(p, dims=(2, 2), sysa=0, sysb=1, method='COBYLA', tol=1e-12, maxiter=2**14) Quantum Discord for two qubit density operator. If ``len(dims) > 2``, then the non-target dimensions will be traced out first. :param p: State to compute quantum discord for. :type p: ket vector or density operator :param dims: The internal dimensions of ``p``. :type dims: tuple(int), optional :param sysa: Index of the first subsystem, A, relative to ``dims``. :type sysa: int, optional :param sysb: Index of the first subsystem, B, relative to ``dims``. :type sysb: int, optional :rtype: float .. py:function:: trace_distance(p1, p2) Trace distance between two states: .. math:: \delta(\rho, \sigma) = \frac{1}{2} \left| \rho - \sigma \right|_\mathrm{tr} If two wavefunctions are supplied the trace distance will be computed via the more efficient expression: .. math:: \delta(|\psi\rangle\langle\psi|, |\phi\rangle\langle\phi|) = \sqrt{1 - \langle \psi | \phi \rangle^2} :param p1: The first state. :type p1: ket or density operator :param p2: The second state. :type p2: ket or density operator :rtype: float .. py:function:: cprint(psi, prec=6) Print a state in the computational basis. :param psi: The pure state. :type psi: vector :param prec: How many significant figures to print. :type prec: int, optional .. rubric:: Examples >>> cprint(rand_ket(2**2)) (-0.0751069-0.192635j) |00> (0.156837-0.562213j) |01> (-0.307381+0.0291168j) |10> (-0.14302+0.707661j) |11> >>> cprint(w_state(4)) (0.5+0j) |0001> (0.5+0j) |0010> (0.5+0j) |0100> (0.5+0j) |1000> .. py:function:: decomp(a, fn, fn_args, fn_d, nmlz_func, mode='p', tol=0.001) Decomposes an operator via the Hilbert-schmidt inner product. Can both print the decomposition or return it. :param a: Operator to decompose. :type a: ket or density operator :param fn: Function to generate operator/state to decompose with. :type fn: callable :param fn_args: Sequence of args whose permutations will be supplied to ``fn``. :param fn_d: The dimension of the operators that `fn` produces. :type fn_d: int :param nmlz_func: Function to produce a normlization coefficient given the ``n`` permutations of operators that will be produced. :type nmlz_func: callable :param mode: String, include ``'p'`` to print the decomp and/or ``'c'`` to return OrderedDict, sorted by size of contribution. :param tol: Print operators with contirbution above ``tol`` only. :returns: Pauli operator name and expec with ``a``. :rtype: None or OrderedDict .. seealso:: :py:obj:`pauli_decomp`, :py:obj:`bell_decomp` .. py:data:: pauli_decomp Decompose an operator into Paulis. .. py:data:: bell_decomp Decompose an operator into bell-states. .. py:function:: correlation(p, A, B, sysa, sysb, dims=None, sparse=None, precomp_func=False) Calculate the correlation between two sites given two operators. :param p: State to compute correlations for, ignored if ``precomp_func=True``. :type p: ket or density operator :param A: Operator to act on first subsystem. :type A: operator :param B: Operator to act on second subsystem. :type B: operator :param sysa: Index of first subsystem. :type sysa: int :param sysb: Index of second subsystem. :type sysb: int :param dims: Internal dimensions of ``p``, will be assumed to be qubits if not given. :type dims: tuple of int, optional :param sparse: Whether to compute with sparse operators. :type sparse: bool, optional :param precomp_func: Whether to return result or single arg function closed over precomputed operator. :type precomp_func: bool, optional :returns: The correlation, - , or a function to compute for an arbitrary state. :rtype: float or callable .. py:function:: pauli_correlations(p, ss=('xx', 'yy', 'zz'), sysa=0, sysb=1, sum_abs=False, precomp_func=False) Calculate the correlation between sites for a list of operator pairs choisen from the pauli matrices. :param p: State to compute correlations for. Ignored if ``precomp_func=True``. :type p: ket or density operator :param ss: List of pairs specifiying pauli matrices. :type ss: tuple or str :param sysa: Index of first site. :type sysa: int, optional :param sysb: Index of second site. :type sysb: int, optional :param sum_abs: Whether to sum over the absolute values of each correlation :type sum_abs: bool, optional :param precomp_func: whether to return the values or a single argument function closed over precomputed operators etc. :type precomp_func: bool, optional :returns: Either the value(s) of each correlation or the function(s) to compute the correlations for an arbitrary state, depending on ``sum_abs`` and ``precomp_func``. :rtype: list of float, list of callable, float or callable .. py:function:: ent_cross_matrix(p, sz_blc=1, ent_fn=logneg, calc_self_ent=True, upscale=False) Calculate the pair-wise function ent_fn between all sites or blocks of a state. :param p: State. :type p: ket or density operator :param sz_blc: Size of the blocks to partition the state into. If the number of individual sites is not a multiple of this then the final (smaller) block will be ignored. :type sz_blc: int :param ent_fn: Bipartite function, notionally entanglement :type ent_fn: callable :param calc_self_ent: Whether to calculate the function for each site alone, purified. If the whole state is pure then this is the entanglement with the whole remaining system. :type calc_self_ent: bool :param upscale: Whether, if sz_blc != 1, to upscale the results so that the output array is the same size as if it was. :type upscale: bool, optional :returns: array of pairwise ent_fn results. :rtype: 2D-array .. py:function:: qid(p, dims, inds, precomp_func=False, sparse_comp=True, norm_func=norm, power=2, coeff=1) .. py:function:: is_degenerate(op, tol=1e-12) Check if operator has any degenerate eigenvalues, determined relative to mean spacing of all eigenvalues. :param op: Operator or assumed eigenvalues to check degeneracy for. :type op: operator or 1d-array :param tol: How much closer than evenly spaced the eigenvalue gap has to be to count as degenerate. :type tol: float :returns: **n_dgen** -- Number of degenerate eigenvalues. :rtype: int .. py:function:: is_eigenvector(x, A, tol=1e-14) Determines whether a vector is an eigenvector of an operator. :param x: Vector to check. :type x: vector :param A: Matrix to check. :type A: operator :param tol: The variance must be smaller than this value. :type tol: float, optional :returns: Whether ``A @ x = l * x`` for some scalar ``l``. :rtype: bool .. py:function:: page_entropy(sz_subsys, sz_total) Calculate the page entropy, i.e. expected entropy for a subsytem of a random state in Hilbert space. :param sz_subsys: Dimension of subsystem. :type sz_subsys: int :param sz_total: Dimension of total system. :type sz_total: int :returns: **s** -- Entropy in bits. :rtype: float .. py:function:: heisenberg_energy(L) Get the analytic isotropic heisenberg chain ground energy for length L. Useful for testing. Assumes the heisenberg model is defined with spin operators not pauli matrices (overall factor of 2 smaller). Taken from [1]. [1] Nickel, Bernie. "Scaling corrections to the ground state energy of the spin-½ isotropic anti-ferromagnetic Heisenberg chain." Journal of Physics Communications 1.5 (2017): 055021 :param L: The length of the chain. :type L: int :returns: **energy** -- The ground state energy. :rtype: float