Hamiltonian of a flux qubit-LC oscillator circuit in the deep–strong-coupling regime

We derive the Hamiltonian of a superconducting circuit that comprises a single-Josephson-junction flux qubit inductively coupled to an LC oscillator, and we compare the derived circuit Hamiltonian with the quantum Rabi Hamiltonian, which describes a two-level system coupled to a harmonic oscillator. We show that there is a simple, intuitive correspondence between the circuit Hamiltonian and the quantum Rabi Hamiltonian. While there is an overall shift of the entire spectrum, the energy level structure of the circuit Hamiltonian up to the seventh excited states can still be fitted well by the quantum Rabi Hamiltonian even in the case where the coupling strength is larger than the frequencies of the qubit and the oscillator, i.e., when the qubit-oscillator circuit is in the deep–strong-coupling regime. We also show that although the circuit Hamiltonian can be transformed via a unitary transformation to a Hamiltonian containing a capacitive coupling term, the resulting circuit Hamiltonian cannot be approximated by the variant of the quantum Rabi Hamiltonian that is obtained using an analogous procedure for mapping the circuit variables onto Pauli and harmonic oscillator operators, even for relatively weak coupling. This difference between the flux and charge gauges follows from the properties of the qubit Hamiltonian eigenstates.


S1. FITTING OF H circ UP TO THE SEVENTH EXCITED STATE
Transition frequencies of the qubit-oscillator circuit numerically calculated fromĤ circ up to the seventh excited state are plotted in Figure S1 for L c = 350 pH using the same circuit parameters as in Fig. 3. The calculated spectra can still be fitted well byĤ R . In the fitting, H R are numerically diagonalized and then the parameters ω, ∆ q , g, and I p , are varied to obtain the best fit. The average of the squares of the residuals in the least-squares method of obtaining the parameters ofĤ R by fitting, [δω 0i /(2π)] 2 (i = 1,2,3,4,5,6, and 7), is 152 MHz 2 .
This value is almost 10 times larger than the case of the fitting up to the third excited state, but still within a moderate level considering that the largest transition frequency is more than 20 GHz.

S2. GAUGE TRANSFORMATION
The circuit Hamiltonian in the flux gaugeĤ circ can be written aŝ The Hamiltonian can now be transformed into the charge gauge rather easily. First we note that the unitary operatorÛ transforms the flux and charge operators as follows: Then, if we set α = L LC /L 12 , the circuit Hamiltonian is transformed intô

S3. ENERGY SHIFTS UP TO SECOND ORDER IN PERTURBATION THEORY
To investigate the reason for the poor reproducibility of the characteristic spectra ofĤ ′ circ byĤ ′ R , we calculate the energy shifts of the four lowest energy levels using perturbation theory. Let us consider the HamiltonianĤ circ,λ ≡Ĥ 1 +Ĥ 2 + λĤ 12 . Here, we assume that the eigenenergies and the eigenstates ofĤ circ,λ can be described by a power series in λ as ni and ni (0) are the eigenenergies and the eigenstates of the non-interacting HamiltonianĤ 1 +Ĥ 2 , n is the number of photons in the oscillator, and i = g, e represents the eigenstate ofĤ 2 . By taking the eigenenergies and the eigenstates ofĤ circ,λ and multiplying by mj (k) from the left, the correction terms can be written as Here, χ ni,mj is the energy shift of the state |ni⟩ due to the interaction with the state |mj⟩, where m = 0, 1, 2, . . . , j = g, e, f, h, . . . , and |f ⟩ and |h⟩ respectively represent the second and third excited states ofĤ 2 . We find that E (1) ni = 0 for all the combinations of n = 0, 1 and i = g, e, and the total energy shifts up to the second-order perturbation χ ni are simply determined by the second-order correction terms: χ ni = m,j χ ni,mj . We also find that the third-or higher-order correction terms are much smaller than the second-order correction term: the numerator of the nth order correction term is the product of n matrix elements ofĤ 12 , which are at most several hundred MHz, while the denominator of the nth order correction term is the product of n − 1 energy-level differences, which are in the order of GHz or more.
Figures S2(a)-(d) show nonzero energy level shifts χ ni,mj of the non-interacting Hamilto-nianĤ 1 +Ĥ 2 due to the interaction termĤ 12 for the lowest four eigenstates |ni⟩ = |0g⟩, |0e⟩, |1g⟩, and |1e⟩ of the non-interacting HamiltonianĤ 1 +Ĥ 2 for L c = 20 pH. Although some energy shifts χ ni,mj appear to be equal to zero everywhere, they are in fact finite but extremely small. The level shifts χ ni,mj that involve the higher qubit levels using dashed lines are all close to zero, and the net level shifts are almost completely determined by the eigenstates involving the two lowest energy levels ofĤ 2 . Figure S2(e) shows the net change of the transition frequencies δω g 01 = χ 1g − χ 0g and δω e 01 = χ 1e − χ 0e from the resonance frequency of the oscillator described byĤ 1 , which reproduce the peaks and dips that occur in the spectrum at Φ x /Φ 0 = 0.5 shown in Fig. 4a in the main text.
We now consider the case of the circuit HamiltonianĤ ′ circ . Figures S3(a)-(d) show nonzero energy level shifts χ ni,mj of the non-interacting HamiltonianĤ ′ 1 +Ĥ ′ 2 due to the interaction with the state |mj⟩ for eigenstates |ni⟩ (n = 0, 1, i = g, e) for L c = 20 pH. Figure S3(e) shows the net change of the transition frequencies δω g 01 = χ 1g − χ 0g and δω e 01 = χ 1e − χ 0e from the resonance frequency of the oscillator described byĤ ′ 1 . Although the qualitative agreement is still not excellent, the peaks and dips that occur in the spectrum at Φ x /Φ 0 are qualitatively reproduced. The agreement of δω g 01 is better than that of δω e 01 . Contrary to and δω e 01 = χ 1e − χ 0e , obtained by combining the shifts in Panels (a)-(d). Numerically calculated transition frequencies ofĤ circ , ω 02 −ω ′ and ω 13 −ω ′ , are plotted as circles, where ω ′ is the resonance frequency ofĤ ′ circ .
the case of the circuit HamiltonianĤ circ , the level shifts χ ni,mj that involve the higher qubit levels j = f, h plotted using dashed lines are larger than those including j = g, e in most of the flux bias conditions, and the level shifts are mainly determined by the third and the fourth lowest energy levels ofĤ ′ 2 . These results explain the difference between the spectra ofĤ ′ circ and those ofĤ ′ R , in whichĤ ′ F Q contains only the lowest two energy levels, as shown in Fig. 4e in the main text: The spectrum represented by solid grey and red lines obtained fromĤ ′ R are similar to those of an uncoupled qubit-oscillator circuit, whereas the spectrum represented by the open grey and red circles obtained fromĤ ′ circ show the peaks and dips around the symmetry point and the overall frequency is lower.

S4. MATRIX ELEMENTS OF THE CHARGE AND FLUX OPERATORS
In Section S3, the energy level shifts of the non-interacting Hamiltonian caused by the interaction term were studied using perturbation theory. The energy shifts of the states |ni⟩ are given in Eq. (S9). Considering that the states ni (0) and mj (0) are separable, the following matrix elements can be written as products of matrix elements of the oscillator and the qubit: and The matrix elements of the oscillator operatorsΦ 1 andq 1 are analytically given as and Only matrix elements that satisfy m = n ± 1 are nonzero. The matrix elements of the qubit operatorsΦ 2 andq 2 for i = g, e and j = g, e, f, h, k, l, where |k⟩ and |l⟩ respectively represent the fourth and fifth excited states ofĤ 2 , are numerically calculated as shown in Fig. 6 in the main text. Note that the matrix elements of the qubit operators at ε = 0 and with a quadratic-plus-quartic potential energy function were studied in Ref. [1].
Regarding the matrix elements j Φ 2 i (i = g, e), those involving the higher qubit levels j = f, h, k, l are smaller than those of j = g, e. Together with the fact that the energy mj of j = f, h, k, l is larger than those of j = g, e, which appears in the denominator of Eq. (S9), this result explains why the level shifts are almost completely determined by the eigenstates involving the two lowest energy levels ofĤ 2 . Regarding the matrix elements |⟨j |q 2 | i⟩| (i = g, e) on the other hand, some of those involving the higher qubit levels j = f, h, k, l are significantly larger than those of j = g, e in most of the flux mj is larger than those of j = g, e, the large matrix elements result in a situation where the level shifts are mainly determined by the third and fourth lowest energy levels ofĤ ′ 2 . We note here that the matrix elements involving the levels j = k, l are smaller than those involving the levels j = f, h and the energy difference is larger than those of j = f, h. Hence, the level shifts caused by higher levels with j = k, l are smaller than those with j = f, h.
The relation between the different matrix elements can be intuitively understood using an analogy to a basic quantum physics problem. The qubit HamiltoniansĤ 2 andĤ ′ 2 have the same form as the Hamiltonian of a single particle in a trapping potential, withΦ 2 and q 2 playing the roles of the position and momentum variables, respectively. Around the symmetry point, the trapping potential is a double-well potential. The energy eigenstates in the double-well potential are superpositions of the energy eigenstates in the two separate wells, which are approximately harmonic oscillator potentials. We therefore have superpositions of harmonic oscillator states whose centers are separated by a distance that is larger than their widths. Figures S4(a)-(d) show wave functions of the eigenstates ofĤ 2 , Ψ i (i = g, e, f, h), as functions of Φ 2 at the flux bias point Φ x /Φ 0 = 0.5. Matrix elements of the position operator can be written using the wave-function representation: Matrix elements that combine a pair of functions Ψ j and Φ 2 Ψ i [Figs. S4(e)-(h)] with similar shapes are large. At the symmetry point, the matrix elements e Φ 2 g = g Φ 2 e ∼ h Φ 2 f = f Φ 2 h are large and the others are small or zero. Away from the symmetry point, the matrix elements g Φ 2 g , e Φ 2 e , f Φ 2 f , and h Φ 2 h also become large. The matrix elements between states with i = g, e and states with j = f, h are therefore small compared to matrix elements involving only g and e. Similarly, the matrix elements of the momentum operator can be written using the wave-function representation: Matrix elements that have a pair of functions Ψ j and ∂Ψ i /∂Φ 2 [Figs. S4(i)-(l)] with similar shapes are large. At the symmetry point, the matrix elements |⟨h |q 2 | g⟩| = |⟨g |q 2 | h⟩| ∼ |⟨f |q 2 | e⟩| = |⟨e |q 2 | f ⟩| are large and the others are small or zero. Away from the symmetry point, the matrix elements |⟨f |q 2 | g⟩|, |⟨g |q 2 | f ⟩|, |⟨h |q 2 | e⟩|, and |⟨e |q 2 | h⟩| also become large.
The matrix elements between states with i = g, e and states with j = f, h are therefore large compared to matrix elements involving only g and e.