How do you put a qubit, a resonator and a coupling in one model, sweep it, and hand it to a time-evolution solver?
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.
HilbertSpace is the tensor product of topic 24 with a memory of which factor is which, and the Hamiltonian it returns is an ordinary QuTiP object.Sweep the flux of a tunable transmon past a resonator at 7 GHz. The dashed lines are the bare qubit and the bare resonator, which cross; the solid lines are the dressed levels, which cannot. Move the coupling and watch the gap at the crossing grow in proportion; the marker reads off its size, and the strip under the plot compares it with the textbook 2g|⟨0|n|1⟩|.
Honest limits: the explorer uses the real transmon matrix elements, the charge coupling and the same level counts as the code below, so it agrees with scqubits’ ParameterSweep to the precision tabulated in the verifier. Two lines are drawn, the two lowest excited levels; the full sweep has more.
Topic 24 said the state space of two systems is the product of their spaces. Topic 31 built a qubit and a resonator that way by hand. scqubits does the same job with the circuit physics already in. You make each subsystem as an object, hand the list to HilbertSpace, and describe the coupling as a term: a strength, an operator on one subsystem and an operator on the other. truncated_dim on each subsystem says how many of its own levels to keep before the product is formed, so a transmon of 61 charge states and a resonator of 10 photon levels become a space of 4 × 10 = 40 and not 61 × 10: the first thing a subsystem does is diagonalise itself and throw away what the problem will not reach.
Dressed and bare. The coupling mixes the two subsystems, so the eigenstates of the whole are not products. scqubits labels each eigenstate by the bare state it overlaps most, and generate_lookup() builds that table. After it, energy_by_bare_index((1, 0)) means “the energy of the dressed state that is mostly transmon level 1 and resonator level 0”, which is how you read a frequency off a coupled system. The qubit frequency of the example, 5.67429 GHz, is 8.3 MHz below the bare transmon’s 5.68258 GHz: the resonator pushes it down, which is the same mechanism as the dispersive shift of topic 31.
A ParameterSweep is a labelled loop. It does what you would write by hand (change a parameter, rebuild, diagonalise, label, store) and caches the result, so that a later question such as “where is the avoided crossing” reads the stored table. You give it the Hilbert space, the values to sweep, and a function saying what to change; subsys_update_info says which subsystems need re-diagonalising when that parameter moves, which saves almost all the time. The scqubits Explorer puts sliders on top of one in a notebook. In the example, the crossing between the tunable transmon and the resonator comes at a flux of 0.310, where the closest approach of the two levels is 125.2 MHz.
The bridge to QuTiP. hs.hamiltonian() returns a QuTiP Qobj with the right dims, and the operators of each subsystem can be turned into Qobjs on the full space the same way. From there it is the solvers of topic 35: put the resonator’s loss and the qubit’s T1 in as collapse operators and drive the qubit with a list-format term. That hand-off is the reason the two packages are used together.
Why a real transmon is not an anharmonic oscillator. In topic 31 the transmon was a ladder with a fixed anharmonicity. Here it is the true circuit. Giving the anharmonic-oscillator model this transmon’s own frequency and anharmonicity and the same 0–1 coupling, the dispersive half-splitting is 1.545 MHz; scqubits, with the real matrix elements (the 1–2 element is 1.37 times the 0–1 element, not √2 = 1.41 times), gives 1.881 MHz. The model is 18% low. The dispersive shift is a small difference of large terms, so it is exactly the quantity that a simplified model gets wrong, and it is the reason a design is checked on the real circuit.
| Call | What it does |
|---|---|
HilbertSpace([sys1, sys2]) | the product space; each subsystem keeps its truncated_dim lowest levels |
hs.add_interaction(g_strength=g, op1=sys1.n_operator, op2=sys2.creation_operator, add_hc=False) | adds g · op1 ⊗ op2. The operator is handed over as a method, not called; add_hc=True would add the Hermitian conjugate instead of writing the second term |
hs.generate_lookup() | diagonalises the whole and labels each eigenstate by its largest bare overlap |
hs.energy_by_bare_index((i, j)) | the energy of the dressed state labelled by bare levels (i, j) |
hs.hamiltonian() | the full Hamiltonian as a QuTiP object |
ParameterSweep(hilbertspace=hs, paramvals_by_name={...}, update_hilbertspace=f, evals_count=n, subsys_update_info={...}) | runs the loop and stores the dressed spectrum, the lookup and, if asked, matrix elements |
The size of the avoided crossing. Close to resonance the two nearest levels repel by 2geff, where geff = g|⟨0|n|1⟩| is the coupling strength times the matrix element of the operator that couples (topic 31). With g = 50 MHz and |⟨0|n|1⟩| = 1.210 at the flux of the closest approach (0.310; the bare transmon crosses 7 GHz a little further on, at 0.312, where the element is 1.207), 2geff = 121.0 MHz; the sweep finds a closest approach of 125.2 MHz, 3.4% more, because the neighbouring levels and the counter-rotating terms push on it too. The numbers: a dressed qubit at 5.67429 GHz, the resonator at 7.00676 GHz with the qubit in |0⟩ and 7.003 GHz with it in |1⟩, a half-splitting of 1.881 MHz against 1.545 MHz for the anharmonic-oscillator model, and a swept closest approach of 125.2 MHz at flux 0.310.
Try it in Python. This is the most useful block of the part: it is the whole path from a circuit to a number to a QuTiP Hamiltonian. It needs the pinned environment of topic 33: with QuTiP 5.3.1, building the Hilbert space raises the csc_array error listed there. Run on scqubits 4.3.1 with QuTiP 5.2.0, 25 October 2026:
import numpy as np
import qutip as qt
import scqubits as scq
scq.settings.PROGRESSBAR_DISABLED = True
# a transmon and a resonator, joined by a charge coupling (GHz)
tm = scq.Transmon(EJ=15.0, EC=0.3, ng=0.0, ncut=30, truncated_dim=4) # keep its lowest four levels
res = scq.Oscillator(E_osc=7.0, truncated_dim=10) # ten photon levels
ev, vec = tm.eigensys(evals_count=4)
n01 = abs(vec[:, 0].conj() @ (tm.n_operator() @ vec[:, 1]))
g = 0.1 / n01 # chosen so the qubit-resonator coupling g*<0|n|1> is 0.1 GHz, as in topic 31
hs = scq.HilbertSpace([tm, res])
hs.add_interaction(g_strength=g, op1=tm.n_operator, op2=res.creation_operator, add_hc=False)
hs.add_interaction(g_strength=g, op1=tm.n_operator, op2=res.annihilation_operator, add_hc=False)
print("dimension:", hs.dimension)
hs.generate_lookup() # label the dressed levels by the bare ones they came from
f_q = hs.energy_by_bare_index((1, 0)) - hs.energy_by_bare_index((0, 0))
f_r0 = hs.energy_by_bare_index((0, 1)) - hs.energy_by_bare_index((0, 0)) # resonator, qubit in |0>
f_r1 = hs.energy_by_bare_index((1, 1)) - hs.energy_by_bare_index((1, 0)) # resonator, qubit in |1>
print(round(f_q, 5), round(f_r0, 5), round(f_r1, 5))
print("half the dispersive splitting:", round((f_r0 - f_r1) / 2 * 1e3, 3), "MHz")
# the bridge: the same Hamiltonian as a QuTiP object, ready for mesolve
H = hs.hamiltonian()
print(type(H).__name__, H.dims)
print(np.allclose(H.eigenenergies()[:3], hs.eigenvals(evals_count=3)))
# the anharmonic-oscillator model of topic 31, given this transmon's own frequency and anharmonicity and the same coupling
wq, alpha = ev[1] - ev[0], (ev[2] - ev[1]) - (ev[1] - ev[0])
N, L = 10, 4
a = qt.tensor(qt.destroy(N), qt.qeye(L))
b = qt.tensor(qt.qeye(N), qt.destroy(L))
Hd = 7.0 * a.dag() * a + wq * b.dag() * b + 0.5 * alpha * b.dag() * b.dag() * b * b + 0.1 * (a.dag() * b + a * b.dag())
E, S = Hd.eigenstates()
def lev(i, j):
k = qt.tensor(qt.basis(N, i), qt.basis(L, j))
return E[int(np.argmax([abs(k.overlap(s)) ** 2 for s in S]))]
chi_duffing = ((lev(1, 0) - lev(0, 0)) - (lev(1, 1) - lev(0, 1))) / 2
print("anharmonic-oscillator model, same coupling:", round(chi_duffing * 1e3, 3), "MHz")dimension: 40
5.67429 7.00676 7.003
half the dispersive splitting: 1.881 MHz
Qobj [[4, 10], [4, 10]]
True
anharmonic-oscillator model, same coupling: 1.545 MHzAnd the sweep. num_cpus=1 keeps it in one process, which is the safest choice on a laptop and for a first run. Run on scqubits 4.3.1 with QuTiP 5.2.0, 25 October 2026:
import numpy as np
import scqubits as scq
scq.settings.PROGRESSBAR_DISABLED = True
# a flux-tunable transmon swept past a resonator at 7 GHz
tt = scq.TunableTransmon(EJmax=40.0, EC=0.3, d=0.0, flux=0.0, ng=0.0, ncut=30, truncated_dim=4)
res = scq.Oscillator(E_osc=7.0, truncated_dim=8)
hs = scq.HilbertSpace([tt, res])
hs.add_interaction(g_strength=0.05, op1=tt.n_operator, op2=res.creation_operator, add_hc=False)
hs.add_interaction(g_strength=0.05, op1=tt.n_operator, op2=res.annihilation_operator, add_hc=False)
def update(flux): # what to change at each point of the sweep
tt.flux = flux
fluxes = np.linspace(0.0, 0.45, 91)
sweep = scq.ParameterSweep(hilbertspace=hs, paramvals_by_name={"flux": fluxes}, update_hilbertspace=update,
evals_count=8, subsys_update_info={"flux": [tt]}, num_cpus=1)
E = sweep.dressed_specdata.energy_table # dressed energies: one row per flux value
split = E[:, 2] - E[:, 1] # the two lowest excited levels: qubit-like and photon-like
i = int(np.argmin(split))
print("closest approach of the two levels: %.1f MHz at flux %.3f" % (split[i] * 1e3, fluxes[i]))closest approach of the two levels: 125.2 MHz at flux 0.310Every frequency plan, every dispersive-shift check, every hand-off to a pulse simulator
A processor design is a set of these Hilbert spaces, one per qubit with its resonator and one per pair of neighbours, and the questions a review asks (where is each qubit really, how much does it shift when its neighbour is excited, where are the avoided crossings the flux sweep has to avoid) are the table lookups of this topic. The scqubits-examples notebooks demo_hilbertspace.ipynb, demo_parametersweep.ipynb, demo_explorer.ipynb and demo_qutip_fluxoniumcz.ipynb work the same steps through; the last of them, by its name, simulates a fluxonium CZ gate with QuTiP. The newer Circuit class goes one step earlier: given a circuit description it builds the Hamiltonian itself (Chitta et al., New J. Phys. 2022).
Check yourself
In scqubits, what does hs.energy_by_bare_index((1, 0)) return for a coupled transmon and resonator?
The coupling mixes the subsystems, so eigenstates are not products. generate_lookup labels each eigenstate by the bare state it overlaps most, and energy_by_bare_index reads the table; the qubit frequency of the example, 5.67429 GHz, is 8.3 MHz below the bare one.
The tempting wrong answers.