I am trying to calculate mutual entropies using QuTiP, but I am being unsuccessful so far. More specifically, I consider a 2^n x 2^n matrix representing the density operator of a n-qubit bipartite system AB made of system A (first m < n qubits) and B (remaining n-m qubits). No tutorial nor material on the internet addressed this specific task.

For simplicity, let us consider a 1-qubit system A and a 2-qubit system B and a density operator of dimension 8x8 representing AB in computational basis.

More practically in python, let rhoAB = Qobj=(np.random.rand(8,8)), and assume that this is a valid density operator.

How should I call entropy_mutual so that I can get this measure between A and B, in particular, regarding the arguments selA and selB? Ideally, I would call something like entopy_mutual(rhoAB, selA=[1], selB=[2,3]) but this not the approach how the function interprets the subsystems and their respective dimensions.


Best to look at the source code when the documentation isn't helpful enough. The definition of entropy_mutual is

def entropy_mutual(rho, selA, selB, base=e, sparse=False):
    Calculates the mutual information S(A:B) between selection
    components of a system density matrix.

    rho : qobj
        Density matrix for composite quantum systems
    selA : int/list
        `int` or `list` of first selected density matrix components.
    selB : int/list
        `int` or `list` of second selected density matrix components.
    base : {e,2}
        Base of logarithm.
    sparse : {False,True}
        Use sparse eigensolver.

    ent_mut : float
       Mutual information between selected components.

    if isinstance(selA, int):
        selA = [selA]
    if isinstance(selB, int):
        selB = [selB]
    if rho.type != 'oper':
        raise TypeError("Input must be a density matrix.")
    if (len(selA) + len(selB)) != len(rho.dims[0]):
        raise TypeError("Number of selected components must match " +
                        "total number.")

    rhoA = ptrace(rho, selA)
    rhoB = ptrace(rho, selB)
    out = (entropy_vn(rhoA, base, sparse=sparse) +
           entropy_vn(rhoB, base, sparse=sparse) -
           entropy_vn(rho, base, sparse=sparse))
    return out

So we see selA and selB are passed as arguments to compute the partial trace. I am not too familiar with qutip but here is an example computing $S(A:B)$ for $\rho_{AB}$ where $A$ is a qubit system and $B$ is a two-qubit system.

import qutip as qtp
# note there is a rand_dm function
# We should also let qutip know how are systems are partitioned
# This is so it knows how to correctly compute the partial trace
rho = qtp.rand_dm(8, dims=[[2,4],[2,4]])

With the above example we could also specify the second system as two-qubits instead of a four dimensional system i.e.

rho = qtp.rand_dm(8, dims=[[2,2,2],[2,2,2]])
  • $\begingroup$ Thank you! I have figured it out myself a couple of hours before your answer and indeed the catch is to initialize rho calling the arg "dims" explicitly in order to partition it in two (or more) subsystems. $\endgroup$
    – EmFed
    May 1 '20 at 21:08

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.