The Machinery · Topic 34 of 41 · The Workbench: the two programs

QuTiP, part one: quantum objects

What is a Qobj, and which dozen calls do you use every day?

  1. 01
  2. 02
  3. 03
  4. 04
  5. 05
  6. 06
  7. 07
  8. 08
  9. 09
  10. 10
  11. 11
  12. 12
  13. 13
  14. 14
  15. 15
  16. 16
  17. 17
  18. 18
  19. 19
  20. 20
  21. 21
  22. 22
  23. 23
  24. 24
  25. 25
  26. 26
  27. 27
  28. 28
  29. 29
  30. 30
  31. 31
  32. 32
  33. 33
  34. 34
  35. 35
  36. 36
  37. 37
  38. 38
  39. 39
  40. 40
  41. 41
About this section: The Workbench: the two programs

Installation first, with the exact error messages you are likely to meet, and a map of where each program’s own documentation lives. Then QuTiP in two topics (objects, then solvers) and scqubits in two (a qubit, then coupled systems and sweeps). The last of them joins the two programs, which is how they are used together. The code is meant to be copied and changed.

See it

Qobj the numbers (a matrix) + dims + type + isherm remembers what it is A * B product A.dag() dagger tensor(A, B) psi.ptrace(1) expect(A, psi) H.eigenstates()
A Qobj is a matrix that knows what it is. It carries its subsystem structure, its kind and whether it is Hermitian, which is why a mistake (a 3-level operator on a 2-level state) is an error message and not a silent wrong answer.

Try it: look inside some Qobjs

Choose a line of QuTiP from the list. The panel shows what that line builds: the matrix, its dims, shape, type and whether it is Hermitian, and its eigenvalues. Every entry shown was checked against the real library’s own printout by the verifier behind this page.

The intuition

You met the idea of the arrays in topic 24: a state is a column, an operator is a table, a system of two parts is a Kronecker product. QuTiP’s Qobj is those arrays with a label attached. The label says what kind of thing it is (ket, bra, oper for an operator, or super for a superoperator, an operator that acts on density matrices), how it divides into subsystems (dims), and whether it is Hermitian (isherm). That is all a Qobj adds to a NumPy array, and it is enough to make the usual mistakes into errors.

The dozen calls. Build a state with basis(N, k) (the k-th level of an N-level system, counting from 0), add and .unit() to normalise a superposition, and turn a ket into a density matrix with ket2dm. Build operators with sigmax(), sigmay(), sigmaz(), destroy(N), create(N), num(N) and qeye(N). Combine with * (matrix product), +, .dag(), and tensor. Ask questions with expect(operator, state), .tr(), .eigenstates() and .ptrace(which). That set covers almost every example in this part.

The traps, each one of which costs a newcomer an afternoon once. * between two Qobjs is the matrix product, never an elementwise one. ptrace takes the subsystems to keep, not the ones to remove. The first factor of a tensor is subsystem 0, the leftmost. basis(2, 0) is the first basis vector, which in this part is the ground state, and sigmaz() is +1 on it, so a qubit whose ground state is |0⟩ has -0.5 * w * sigmaz() or w * num(2) as its Hamiltonian. Equality with == compares the numbers within a small tolerance, so it is safe to use for checks like the commutator below. And .full() hands you a plain NumPy array when you want to leave QuTiP.

Mixed states are first-class. A pure state is a ket; a state with classical uncertainty is a density matrix, an operator with unit trace. expect accepts either, and ptrace turns the first into the second: tracing out one half of an entangled pair leaves the other half mixed, which is the content of topic 2. The example below shows it on a dressed state of the qubit and the resonator of topic 31.

The mathematics, as code

Each line of the mathematics of the earlier topics has a line of QuTiP. The right-hand column is what you type; the left is what it means.

MathematicsQuTiP
|k⟩ in an N-level systembasis(N, k) (the projector |k⟩⟨k| is fock_dm(N, k))
(|0⟩ + |1⟩)/√2(basis(2, 0) + basis(2, 1)).unit()
ρ = |ψ⟩⟨ψ|ket2dm(psi)
⟨ψ|A|ψ⟩, or Tr(ρA)expect(A, psi), or expect(A, rho)
A ⊗ Btensor(A, B)
TrB ρrho.ptrace(0) (keep subsystem 0)
[A, B]commutator(A, B)
e−iHt(-1j * H * t).expm()
H|n⟩ = En|n⟩E, states = H.eigenstates()
S = −Tr(ρ log₂ ρ)entropy_vn(rho, base=2)

A worked check. Take the qubit at 6 GHz beside the resonator at 7 GHz with g = 0.1 GHz, as in topic 31. In the dressed state with energy 7.0099 GHz, the qubit is excited with probability 0.0097. Perturbation theory predicts (g/Δ)² = 0.01, with Δ = 1 GHz, and the exact figure is 3% below it. The reduced state of the qubit is mixed, with an entropy of 0.0789 bits, which is the amount by which that dressed state is entangled. The numbers are the dressed state at 7.0099 GHz has qubit population 0.0097 and entropy 0.0789 bits, with E0..3 = 0, 5.9901, 7.0099 and 12.98038 GHz.

Try it in Python. The first lines are the table above in use; the last is the mistake the figure promised an error for. Run on QuTiP 5.2.0, 25 October 2026:

import numpy as np
import qutip as qt

# states and expectation values
plus = (qt.basis(2, 0) + qt.basis(2, 1)).unit()
print(round(qt.expect(qt.sigmax(), plus), 6))          # <+|sigma_x|+> = 1
rho = qt.ket2dm(plus)                                    # |+><+| as a density matrix
print(round(rho.tr(), 6), round((rho * rho).tr(), 6))    # trace 1; purity 1 because it is a pure state

# operator algebra is exact, and == compares with a tolerance
sx, sy, sz = qt.sigmax(), qt.sigmay(), qt.sigmaz()
print(qt.commutator(sx, sy) == 2j * sz)

# a resonator (6 levels) and a qubit: dims remember which factor is which
N = 6
a = qt.tensor(qt.destroy(N), qt.qeye(2))                 # acts on the resonator, the first factor
sm = qt.tensor(qt.qeye(N), qt.destroy(2))                # acts on the qubit, the second factor
print(a.dims, a.shape)
H = 7.0 * a.dag() * a + 6.0 * sm.dag() * sm + 0.1 * (a.dag() * sm + a * sm.dag())
E, states = H.eigenstates()
print(np.round(E[:4], 5))
rho_q = states[2].ptrace(1)                              # keep subsystem 1 (the qubit), trace the resonator out
print(round(qt.expect(qt.num(2), rho_q), 4))             # the qubit's excited-state population in that dressed state
print(round(qt.entropy_vn(rho_q, base=2), 4))            # above 0: the dressed state is entangled

# what a mistake looks like: a 3-level operator on a 2-level state
try:
    qt.destroy(3) * qt.basis(2, 0)
except Exception as e:
    print(type(e).__name__ + ":", e)
1.0
1.0 1.0
True
[[6, 2], [6, 2]] (12, 12)
[ 0.       5.9901   7.0099  12.98038]
0.0097
0.0789
TypeError: incompatible dimensions [[3], [3]] and [[2], [1]]

Where it actually matters

The habit that prevents the most wasted time The value of the Qobj is the errors it raises. A dimension mismatch, a missing dag and an operator on the wrong subsystem are all the kind of mistake that, in a bare NumPy script, gives a plausible wrong number. In QuTiP they give a message like the last line above. The user guide chapters that go with this topic are “Basic Operations on Quantum Objects”, “Manipulating States and Operators” and “Using Tensor Products and Partial Traces”. Read the chapter on tensor products before you build your first two-system model, because every model in the next topics is one.