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.

Introduction

Superconducting circuits are one of the most promising platforms for realizing large-scale quantum information processing. One of the most important features of superconducting circuits is the freedom they allow in their circuit design. Since the first demonstration of coherent control of a Cooper pair box[1], various types of superconducting circuits have been demonstrated.

The Hamiltonian of a superconducting circuit can be derived using the standard quantization procedure applied to the charge and flux variables in the circuit[2]. The Hamiltonian of an LC circuit is well known to be that of a harmonic oscillator. The Hamiltonians of various kinds of superconducting qubits have also been well studied[3-7] and these Hamiltonians can be numerically diagonalized to obtain eigenenergies and eigenstates that accurately reproduce experimental data. On the other hand, the Hamiltonian of circuits containing two or more components, e.g., qubit-qubit, qubit-oscillator, or oscillator-oscillator systems, are usually treated in such a way that the Hamiltonian of the individual components and the coupling among them are separately obtained[8-11]. This separate treatment of individual circuit components works reasonably well in most circuits. Even for flux qubit-oscillator circuits in the ultrastrong-coupling regime[12,13], where the coupling strength g is around 10% of the oscillator’s frequency and the qubit minimum frequency , or the deep–strong-coupling regime[14-16], where g is comparable to or larger than and , the experimental data can be well fitted by the quantum Rabi Hamiltonian[17-19], where a two-level atom and a harmonic oscillator are coupled by a dipole-dipole interaction. However, the more rigorous approach based on the standard quantization procedure has not been applied to such circuits except in a few specific studies[20-22], and the validity of describing a flux qubit-oscillator circuit by the quantum Rabi Hamiltonian has been demonstrated in only a few specific circuits[23,24].

In this paper, we apply the standard quantization procedure to a superconducting circuit in which a single-Josephson-junction flux qubit (an rf-SQUID qubit or a fluxonium-equivalent circuit) and an LC oscillator are inductively coupled to each other by a shared inductor (Fig. 1a), and derive the Hamiltonian of the circuit. Our model with a single Josephson junction should be sufficient for the purposes of the present study, and the results can be applied to circuits of commonly used multi-Josephson-junction flux qubits[25,26], including the fluxonium. Note that single-Josephson-junction flux qubits with different parameters have been experimentally demonstrated using superinductors, which have been realized by high-kinetic-inductance superconductors[27,28], granular aluminum[29], and Josephson-junction arrays[30]. The derived circuit Hamiltonian consists of terms associated with the LC oscillator, the flux qubit (and its higher energy levels), and the product of the two flux operators. Excluding the qubit’s energy levels higher than the first excited state, this circuit Hamiltonian takes the form of the quantum Rabi Hamiltonian, which describes a two-level system coupled to a harmonic oscillator. To investigate contributions from the qubit’s higher energy levels, we numerically calculate the transition frequencies of the circuit Hamiltonian. We find that the qubit’s higher energy levels mainly cause a negative shift of the entire spectrum, and that the calculated transition frequencies are well fitted by the quantum Rabi Hamiltonian even when the qubit-oscillator circuit is in the deep–strong-coupling regime. The circuit Hamiltonian can be transformed to one that has a capacitive coupling term by a unitary transformation. We show, however, that the spectrum of the circuit Hamiltonian cannot be fitted by the variant of the quantum Rabi Hamiltonian that has different Pauli operators in the qubit and coupling terms. This situation arises when we perform the mapping from circuit variables to Pauli operators for a circuit Hamiltonian that has a capacitive coupling term, i.e. a Hamiltonian expressed in the charge gauge. We explain the advantage of the flux gauge over the charge gauge in this regard based on the properties of the eigenstates of the flux qubit Hamiltonian.

Figure 1

Circuit diagrams. (a) A superconducting circuit in which a single-Josephson-junction flux qubit and an LC oscillator are inductively coupled to each other by a shared inductor. (b) Equivalent circuit obtained by applying the so-called Y- transformation to the inductor network in circuit (a). (c) The outer loop of circuit (a), which forms an LC oscillator. (d) The inner loop of circuit (a), which forms a single-Josephson-junction flux qubit.

Results

Circuit Hamiltonian

Following the standard quantization procedure, nodes are assigned to the circuit as shown in Fig. 1a. Before deriving the circuit Hamiltonian, the circuit in Fig. 1a is transformed to the one shown in Fig. 1b by applying the so-called Y- transformation, by which a -shaped network of electrical elements is converted to an equivalent Y-shaped network or vice versa, to the inductor network. Thus, node 3 surrounded by the inductors is eliminated. The inductances of the new set of inductors are given asandThe Lagrangian of the circuit can now be obtained relatively easily[2]:where and are the capacitance and the Josephson energy of the Josephson junction, is the critical current of the Josephson junction, is the superconducting flux quantum, and and are the node flux and its time derivative for node k. The node charges are defined as the conjugate momenta of the node fluxes as . After the Legendre transformation, the Hamiltonian is obtained aswhereandAs can be seen from Eq. (5), the Hamiltonian can be separated into three parts: the first part consisting of the charge and flux operators of node 1, the second part consisting of node 2 operators, and the third part containing the product of the two flux operators.

Separate treatment of the qubit-oscillator circuit

Let us consider an alternative treatment of the circuit shown in Fig. 1a, where the circuit is assumed to be naively divided into two well-defined components. The capacitor and the inductors in the outer loop of the circuit in Fig. 1a form an LC oscillator (Fig. 1c). The Josephson junction and the inductors in the inner loop form a single-Josephson-junction flux qubit (Fig. 1d). The LC oscillator and the single-Josephson-junction flux qubit share the inductor at the common part of the two loops. It is now instructive to investigate the relation of the following pairs of Hamiltonians: and the Hamiltonian of the LC oscillator shown in Fig. 1c, and the Hamiltonian of the flux qubit shown in Fig. 1d, and and the Hamiltonian of the inductive coupling between the LC oscillator and the flux qubit , where we have used the relations of the oscillator current , the qubit current , and the mutual inductance . Actually, only the inductances are different in the Hamiltonians of the two pictures: The inductances of the LC oscillator, the flux qubit, and the coupling Hamiltonians derived from the separate treatment are respectively , , and , while those in are , , and . Figure 2 shows the inductances in the separate treatment and in as functions of on condition that the inductance sums are kept constant at  pH and  pH. As approaches 0, or more specifically when , we obtain , , and . In this way, the Hamiltonian derived from the separate treatment has the same form as , and the inductances in the separate treatment approach those of .

Figure 2

(a) The inductances of the LC oscillator , (b) the inductance of the flux qubit , and (c) the inverse inductance of , obtained by Eqs. (9), (10), and (3) (solid lines) are plotted as functions of together with their counterparts in the separate treatment (dotted lines) on condition that the inductance sums are kept constant at  pH and  pH.

Comparison to the quantum Rabi Hamiltonian

We compare the Hamiltonian with the generalized quantum Rabi Hamiltonian:The first part represents the energy of the LC oscillator, where and are the creation and annihilation operators, respectively. The second part represents the energy of the flux qubit written in the energy eigenbasis at . The operators are the standard Pauli operators. The parameters and are the tunnel splitting and the energy bias between the two states with persistent currents flowing in opposite directions around the qubit loop. The third part represents the inductive coupling energy.

The relation between and is straightforward. The resonance frequency and the operators in can be analytically described by the variables and operators in as and , where is the zero-point-fluctuation current. To see the relation between and , we numerically calculated the eigenenergies of as functions of . In the calculation, we used the following parameters:  pH,  pH ( GHz), and  fF ( GHz). As shown in Fig. 3a the lowest two energy levels are well separated from the higher levels, which are more than 40 GHz higher in frequency. The lowest two energy levels of are well approximated by (Fig. 3b), which gives almost identical results obtained by the local basis reduction method[31], with and swapped. The eigen frequencies of the ground state and the first excited state are respectively fitted by and . Besides the offset , the fitting parameters are and the maximum persistent current , which is determined as the proportionality constant between the energy bias and the flux bias, . We also numerically calculated the expectation values of the flux operator, and , which are well approximated by and , respectively (Fig. 3c). The expectation values of the flux operator and are respectively fitted by and . Here, and are respectively the ground and excited states of the qubit Hamiltonian . The only fitting parameter is determined by the ratio (). The Pauli operator is therefore identified as being proportional to the flux operator . The relation between and can now be obtained by using the relations for the oscillator and qubit operators identified above. This way we find that the Hamiltonian can be expressed as , which is exactly the same form as , with the coupling strength . Note that directly derived from the Lagrangian is in the flux gauge, as the qubit-oscillator coupling term is of the form , which is optimal for our system with a single oscillator mode [32-34]. Excluding the qubit’s energy levels higher than the first excited states, takes the form of . In other words, once the circuit parameters, i.e. , , , , C, , and , are given, the corresponding parameters in (, , , and g) are set. The relation between and is summarized in Table 1.

Figure 3

(a) Numerically calculated energy levels of as functions of . (b) The lowest two energy levels of . (c) The expectation values of the flux operator and . In (b,c), the solid circles are obtained from numerical calculations of , while the lines are obtained from fitting by . The black and red colors respectively indicate states and . In the calculation, we used the following parameters:  pH,  pH ( GHz), and  fF ( GHz).

Table 1

Hamiltonians, operators, and variables in and their counterparts in .

H^_R	H^_circ
H^_LC	H^_1
H^_FQ	H^_2
H^_coup	H^_12
a^+a^^†	Φ^_1/(L_LCI_zpf)
(a^-a^^†)/i	q^_1/(CV_zpf)
σ^_x	-Φ^_2/Φ_2max
σ^_z	–
ω	1/L_LCC
ε	2I_p(Φ_x-0.5Φ_0)
Δ_q	Minimum qubit frequency (numerically evaluated)
g	(L_LC/L_12)I_zpfΦ_2max/ħ

and are the zero-point-fluctuation current and voltage, respectively.

The parameters , , and are obtained by numerically calculating the eigenenergies of as functions of and fitting the lowest two energy levels by . Note that there is no analytic expression for and .

As previously mentioned, considers only the lowest two energy levels of the flux qubit, while includes all energy levels. To investigate the effect of the qubit’s higher energy levels, we perform numerical calculations and compare the energy levels calculated by and . The details of the numerical diagonalization of are given in “Methods”. Since the contributions from the qubit’s higher energy levels are expected to become larger as the coupling strength increases, we consider parameters that cover a wide range of coupling strengths from the weak-coupling to the deep–strong-coupling regime. In the calculations, we fix  pH,  pH,  pF,  pH ( GHz), and  fF ( GHz), and sweep the flux bias around at various values of . These parameters are used in all the calculations for Figs. 4, 5, 6 and 7. Some of our calculations were performed using the QuTiP simulation package[35].

Figure 4

Transition frequencies of the qubit-oscillator circuit fitted with different approximate Hamiltonians. Numerically calculated transition frequencies from are plotted as circles while the fitting spectra are plotted as lines. Two mutual inductance values, corresponding to different qubit-oscillator coupling strengths, are used: = 20 pH (a,c,e,g) and = 350 pH (b,d,f,h). The spectra are fitted with (a–d) and (e–h). (a,b) The parameters of are obtained from the relations described in Table 1:  GHz,  GHz,  GHz, and  nA for  pH and  GHz,  GHz,  GHz, and  nA for  pH. (c,d) The parameters of are obtained by fitting the spectra of to :  GHz,  GHz,  GHz, and  nA for  pH and  GHz,  GHz,  GHz, and  nA for  pH. (e,f) The parameters of are obtained from the circuit parameters of in the similar way as in Table 1:  GHz,  GHz, 0.043 GHz, and  nA for  pH, and  GHz,  GHz, 0.492 GHz, and  nA for  pH. (g,h) The parameters of are obtained by fitting the spectra of to :  GHz,  GHz,  GHz, and  nA for  pH and  GHz,  GHz,  GHz, and  nA for  pH. Black, gray, orange, magenta, and red colors indicate transitions , , , , and , respectively.

Figure 5

Parameters of obtained from the relations described in Table 1, (red solid lines) and by fitting the spectra of up to the third excited state to (black dashed lines) as a function of . Panels (a–d) respectively correspond to the oscillator frequency , qubit frequency , coupling strength g, and persistent current in . (e) The average of the squares of the residuals in the least-squares method for obtaining the parameters of by fitting, (i = 1, 2, and 3).

Figure 6

Numerically calculated matrix elements for and (a) , (b) , where , , , and respectively represent the second, third, fourth, and fifth excited states of . Numerically calculated matrix elements for and (c) , (d) . Note that the x axis ranges are smaller than those in Fig. 4.

Figure 7

Expectation values of photon number, flux, and current operators for (a–d)  and for (e,f) . (a,e) Expectation numbers of photons in the oscillator in the ground state as a function of . The solid curve corresponds to , while the dashed line in panel (a) indicates the simple expression , which is valid in the limit . (b,f) [d,h] Expectation values of fluxes [Expectation values of currents ] as functions of the flux bias in the case  pH. Black and red colors respectively indicate and 2, and solid and dashed lines indicate the ground and the excited states, respectively. (c,g) Expectation value of flux in the ground state as a function of at .

Transition frequencies of the qubit-oscillator circuit corresponding to the transition numerically calculated from around the resonance frequency of the oscillator are plotted in Fig. 4 for (a)  pH and (b)  pH, where the indices i and j label the energy eigenstates according to their order in the energy-level ladder, with the index 0 denoting the ground state. The same circuit parameters as in Fig. 3 are used. In the case  pH, two characteristic features are observed: avoided level crossings between the qubit and oscillator transition signals approximately at and 0.503, and the dispersive shift that for example creates the separation between the frequencies of the transitions and , leading to the peak/dip feature at the symmetry point, i.e. . Note that transitions and around respectively correspond to transitions and , where “g” and “e” denote, respectively, the ground and excited states of the qubit, and “0” and “1” the number of photons in the oscillator’s Fock state. In the case  pH, the characteristic spectrum indicates that the qubit-oscillator circuit is in the deep–strong-coupling regime[15]. Transition frequencies of are also plotted in Fig. 4a,b. It should be mentioned that the parameters of are obtained in two different ways. Here, the parameters of are obtained from the relations described in Table 1. The overall shapes of the spectra of and look similar. On the other hand, the shift of the entire spectrum becomes as large as more than 200 MHz for  pH. To quantify the difference of the spectra between and , transition frequencies up to the third excited state numerically calculated from are fitted by . In the fitting, with an initial approximate set of parameters is numerically diagonalized and then the parameters, , , g, and , are varied to obtain the best fit. Figure 4c,d show that the same spectra of in Fig. 4a,b are well fitted by , but with a different parameter set.

Figure 5 shows parameters of obtained from the relations described in Table 1, and by fitting the spectra of up to the third excited state to as a function of . The parameters of obtained in these two different ways become quite different from each other at large values of due to contributions of the qubit’s higher energy levels. Interestingly, the average of the squares of the residuals in the least-squares method for obtaining the parameters of by fitting, (i = 1, 2, and 3), remains rather small, at most 25 MHz, which is consistent with the good fitting shown in Fig. 4c,d. In fact, the energy level structure of up to the seventh excited state can still be fitted well by the quantum Rabi Hamiltonian as shown in Supplementary Information S1[36].

Hamiltonian in the charge gauge

The circuit Hamiltonian is in the flux gauge, as the qubit-oscillator coupling term is of the form . The Hamiltonian can be transformed into the charge gauge:whereThe details of the calculations are given in Supplementary Information S2[36]. As can be seen in Eq. (12), the Hamiltonian can be separated into three parts: the first part consisting of the charge and flux operators of node 1, the second part consisting of the charge and flux operators of node 2, and the third part containing the product of the two charge operators. It is worth mentioning that the inductance in is equal to , which is the inductance of the flux qubit in the separate treatment:

The overall form of [Eq. (12)] is almost the same as that of [Eq. (5)], with the exception that the coupling term is of the form instead of . When we examined the mapping between and , we explained that the operator can be identified as . Similarly, if we calculate the matrix elements of the operator for the two lowest qubit states, we find that and (Fig. 6a). We can therefore identify the operator as the qubit operator .

It then seems natural to look for a mapping of to the variant of the quantum Rabi Hamiltonian:Compared with , only is replaced by . It should be noted, however, that there is a physical difference between the two Hamiltonians and that there is no simple mapping between them. The operator in is the same as one of the Pauli operators in . In contrast, the operator in is different from the two operators ( and ) in . As a result, and are physically different and for example produce different spectra.

In Fig. 4 we plot a few of the transition frequencies in the circuit’s spectrum along with the corresponding transition frequencies obtained from . When the parameters of are read off , similarly to what is shown in Table 1, Fig. 4e,f show that the fitting is poor in most parts of the spectrum, even for the relatively weak coupling case  pH. Here, , , , , and is numerically calculated as shown in Fig. 6a. Contrary to the open grey and red circles obtained from , the solid grey and red lines obtained from are larger and do not show the peaks and dips around the symmetry point. Even with a numerical optimization of the fitting parameters (Fig. 4g,h), only parts of the spectrum can be fitted well. In particular, the peaks and dips that occur in the spectrum at are not reproduced by .

Here it is useful to consider the two characteristic features in the spectrum, i.e. the avoided crossings and dispersive shifts, as being the result of energy shifts in the spectrum of an uncoupled circuit when the coupling term is added. The gap of an avoided level crossing, or in other words the Rabi splitting, is proportional to the matrix element of the coupling term between the relevant energy eigenstates of the uncoupled system. The dispersive shift of one energy level caused by another energy level is proportional to the square of the matrix element between the two energy eigenstates according to perturbation theory. The details of the dispersive shifts up to second order in perturbation theory are described in Supplementary information S3[36].

The difference between the spectra shown in the solid lines in Fig. 4a,e is attributed to the difference between the matrix elements of the qubit’s flux and charge operators, and . The flux bias dependences of the numerically calculated matrix elements of the qubit’s charge and flux operators are shown in Fig. 6. Regarding the matrix elements of the qubit’s flux operator (), those involving the higher qubit levels are smaller than those of . Here, , , , and respectively represent the second, third, fourth, and fifth excited states of . Regarding the matrix elements of the qubit’s charge operator (), on the other hand, some of those involving the higher qubit levels are significantly larger than those of in most of the flux bias range. This difference stems from the fact that the qubit Hamiltonian involves a double-well potential of the flux variables, as explained in Supplementary information S4[36].

The matrix elements and have a peak at , which directly leads to the peak/dip feature at the symmetry point as shown in the solid lines in Fig. 4a. On the other hand, the matrix elements and are almost constant in the flux bias range of Fig. 6, which is consistent with the absence of a peak/dip feature at the symmetry point in the solid lines in Fig. 4e. Instead, the matrix elements and , which have somewhat similar flux-bias dependence to those of and , have peaks at . As a result, for the flux gauge, the contribution from higher levels is just a small correction, while it cannot be neglected in the charge gauge. For this reason, the flux gauge turns out to be more convenient for purposes of mapping to the quantum Rabi Hamiltonian.

It is worth mentioning that we are dealing with an inductively coupled circuit, and one might think that the nature of the coupling, i.e. inductive or capacitive, will determine which gauge is more suitable, especially because the two gauges differ mainly by the form of the coupling term. We show that the deciding factor is the qubit rather than the coupling. As a result, if we have a flux qubit capacitively coupled to the oscillator, the flux gauge will be more suitable for the purpose of approximating the circuit Hamiltonian by the quantum Rabi Hamiltonian. It is also worth mentioning that our results should apply to fluxonium-resonator circuits, since the fluxonium Hamiltonian close to the degeneracy point also involves a double-well potential of the flux variables. It should be noted, however, that the fluxonium’s transition frequencies to the higher energy levels, and , are smaller than those of flux qubits, and, hence, the contribution from higher levels would be larger in fluxonium-resonator circuits.

We note here that the differences between the energy spectra obtained from , , and in the case , , and a quadratic-plus-quartic potential energy function were extensively studied in Ref.[24]. Contrary to the case of , the spectra obtained from and are almost identical for the lowest energy levels in the whole range of while that of deviates from the other two when . The Rabi splitting between the first two excited energy levels, which is symmetric and is proportional to g, is clearly visible for , while the dispersive shift due to the higher energy levels is proportional to , and can be observed only when . We also note that for  pH,  GHz, which is much smaller than  GHz. On the other hand, the gap of the avoided level crossings plotted as lines in Fig. 4e is not much smaller than that in Fig. 4a considering . The gap of the avoided level crossings for is while that for is . The factor partly explains why the avoided level crossings for and are not very different from each other.

Expectation values of the photon number and the field operator

One of the most paradoxical features of in the deep–strong-coupling regime is the non-negligible number of photons in the oscillator in the ground state[37]. In terms of the creation and annihilation operators, the photon number operator is . Considering the mapping between and , the photon number operator in is . As shown in Fig. 7a, a nonzero number of photons in the oscillator is obtained. Another paradoxical feature of in the deep–strong-coupling regime is a non-negligible expectation value of the field operator when the qubit flux bias is . The corresponding operator in is . Fig. 7b shows the expectation values of the fluxes, and , as functions of the flux bias for the ground and the first excited states in the case  pH. At , a nonzero expectation value is demonstrated. One may suspect that a nonzero is unphysical because it might result in a nonzero DC current through the inductor of an LC oscillator in an energy eigenstate. The relation between the flux and current operators is explicitly given in the circuit model (Fig. 1a):where () is defined as the operator of the current flowing from the ground node to node k. In the case , Eq. (19) reduces to and which means that the operator can indeed be understood as a current operator. As shown in Fig. 7c, at , the expectation value in the ground state is proportional to and is zero only when . Figure 7d shows the expectation values of the currents and as functions of in the case  pH. The expectation value is nonzero at , while is exactly zero at all values of . In this way, although is nonzero in the case , an unphysical DC current through the inductor of an LC oscillator is not predicted.

Next, we consider the case of . Since unitary transformations do not change the eigenenergies of Hamiltonians, and have exactly the same eigenenergies. On the other hand, the photon-number and flux operators do not commute with the gauge transformation:andThe corresponding photon number operator in is , where is the resonance frequency of . The corresponding field operator in is . As shown in Fig. 7e, the expectation number of photons in the charge gauge is much smaller than that in the flux gauge, and the non-negligible expectation number of photons in the oscillator in the ground state arises only in the flux-gauge. We note here that the resonance frequency of is larger than GHz:  GHz in the case  pH. We can see that the expectation value of the operator in the case  pH is zero at all values of (Fig. 7f), which are also different from the case of the flux-gauge Hamiltonian. The expectation value of the current operator is zero at all values of , which is true in the flux gauge as well.

Discussion

We have derived the Hamiltonian of a superconducting circuit that comprises a single-Josephson-junction flux qubit and an LC oscillator using the standard quantization procedure. Excluding the qubit’s higher energy levels, the derived circuit Hamiltonian takes the form of the quantum Rabi Hamiltonian. We show that the Hamiltonian derived from the separate treatment, where the circuit is assumed to be naively divided into the two well-defined components, has the same form as the circuit Hamiltonian, and the inductances in the separate treatment approach those of the circuit Hamiltonian as approaches 0. The qubit’s higher energy levels mainly cause a negative shift of the entire spectrum, but the energy level structure can still be fitted well by the quantum Rabi Hamiltonian even when the qubit-oscillator circuit is in the deep–strong-coupling regime. We also show that although the circuit Hamiltonian can be transformed to a Hamiltonian containing a capacitive coupling term, the resulting circuit Hamiltonian cannot be approximated by the capacitive-coupling variant of the quantum Rabi Hamiltonian.

Methods

As a simple example, let us consider in the main text:where, , , and and are operators of dimensionless magnetic flux and charge, respectively, and satisfy . Using the relation , the Hamiltonian can be rewritten asLet us calculate wavefunctions and their eigenenergies E of this Hamiltonian. We expand the wavefunction with plane waves asHere, is the length of the n-space. Considering the periodic boundary condition , the wave number k is given byIn numerical calculations in this work, we have used 32 waves for the qubit and 64 waves for the oscillator. Then, the equation for determining and E is obtained aswhereThe set of wavefunctions and eigenenergies are obtained by solving Eq. (26). Numerically, we can get by the fast Fourier transform (FFT) for discretized n-space. After solving the eigenvalue equation in Eq. (26), the wavefunctions are also obtained from by the FFT. For the calculation of in the main text, we use the following equation instead of Eq. (26),where , , and . Note that a fluxonium-resonator circuit that has the identical circuit diagram but different circuit parameters was numerically diagonalized using the harmonic oscillator basis[22].

Supplementary Information.