Spectroscopy of flux-driven Kerr parametric oscillators by reflection coefficient measurement

Dressed states formed in a driven quantum system are typically characterized by spectroscopic measurements. For example, if a two-level system is driven on resonance, the Mollow triplet [1] is observed due to the Rabi splitting of the two energy levels hybridized with the photon number states. If, on the other hand, the drive is far detuned, the shift of the energy level due to the Autler–Townes effect, namely, the ac Stark effect [2] is observed.

The experiments originally done in atomic system have not only been reproduced, but also extended to an unexplored parameter regime using superconducting artificial atoms, referred to as a circuit quantum electrodynamics (c-QED) system [3, 4]. For example, the ac Stark shift in a c-QED system was first reported in [5], and extended to the strong dispersive regime, where it is possible to determine the population of each photon number state in the resonator coupled to the qubit [6].

Compared to the experiments using one-photon drive, spectroscopy of a two-photon driven system has been relatively unexplored in c-QED systems. The two-photon drive plays an essential role in a device called a Kerr parametric oscillator (KPO) [7, 8], which has been attracting attention for applications in quantum information processing, such as quantum annealing [9, 10] and gate-model quantum computing [11–13]. KPOs are parametric oscillators, in which the Kerr nonlinearity is larger than the photon loss rate, particularly in the single-photon Kerr regime [14]. The KPOs can be implemented using Josephson parametric oscillators [15–19] or charge-driven transmons with a superconducting nonlinear asymmetric inductive element (SNAIL) [20–22]. In both cases, an external field nearly twice the resonance frequency of the KPOs is used for parametric pumping and works as the two-photon drive.

Spectroscopic measurement for studying the energy-level structure up to the tenth excited states of the KPO of charge-driven transmons has been reported by Frattini et al [21]. By measuring the excitation energy from the ground states as a function of the two-photon drive amplitude, they revealed the crossover from a nondegenerate spectrum to a pairwise kissing spectrum formed by two excited states in the double-well meta potential of the KPO, which leads to a staircase pattern in the coherent state lifetime. This type of spectroscopic measurement for revealing the energy spectra of KPOs is important not only for understanding its fundamental physics, but also for its application to quantum information processing. For example, in [23], a gate operation of the cat qubit using excited energy levels was proposed. Moreover, a comparison of the spectroscopic measurement with the theory generally gives the calibration of the drive amplitude at the quantum chip [5, 24, 25]. This is particularly important in KPOs for precisely determining the two-photon drive amplitude to create an optimal schedule to generate Schrödinger's cat state in a short time [13, 26, 27].

Some of the authors of the present study theoretically proposed an alternative method of spectroscopy for revealing the energy level structure of a KPO under two-photon drive [28]. The method is based on the measurement of reflection coefficient of the KPO and is simpler than that used in previous studies [21, 22] because it is only performed by continuous wave measurement without using microwave pulses. In the present paper, we experimentally demonstrate the proposed method by performing reflectometry measurements of the KPO under the two-photon drive induced by the flux modulation. We compare the results with numerical simulations on the basis of the theory to interpret the experimental data. We also show that the experimental data can be used to determine the calibration of the two-photon drive amplitude, which is important for the application of KPOs in quantum information processing.

The KPO used in this experiment is a lumped-element-type device. An optical microscope image and an equivalent circuit of the KPO are shown in figures 1(a) and (b), respectively. The left end face of the KPO is connected to a signal line via a coupling capacitor ${C_{{\text{in}}}}$, from which a probe signal is applied. The bottom end face is connected to a coupler for interaction with other KPOs, which are not used in this experiment. The top end face is connected to the ground plane via two Josephson junctions (JJs) and a superconducting quantum interference device (SQUID) in series, the latter of which is inductively coupled to the pump line. The JJ is placed in series with the SQUID to reduce the magnitude of the Kerr nonlinearity $\chi$ of the KPO [29]. The SQUID contains asymmetric JJs, whose critical currents are designed to be 375 nA and 625 nA, and those of the series junctions are designed to both be 1000 nA. The actual critical current of each JJ was verified by measuring the room-temperature resistance of a test structure. The KPO is placed in a dilution refrigerator and cooled to below 10 mK.

Zoom In Zoom Out Reset image size

We measure the reflection coefficient of the KPO without the two-photon drive to extract the parameters of the KPO. The reflection coefficient $\Gamma$ is expressed by the following equation,

$$\begin{align}\Gamma = 1 - \frac{{{\kappa _{\text{e}}}}}{{i\left( {{\omega _{{\text{ref}}}} - {\omega _{\text{r}}}} \right) + \frac{{{\kappa _{\text{e}}} + {\kappa _{\text{i}}}}}{2}}},\end{align} \tag{ 1 }$$

where ${\omega _{\text{r}}}$ is the resonance frequency of the KPO, ${\omega _{{\text{ref}}}}$ is the probe frequency, and ${\kappa _{\text{e}}}$ (${\kappa _{\text{i}}}$) is an external (internal) loss rate. Figure 1(c) shows the flux-bias dependence of the resonance frequency. The resonance frequency of the KPO changes periodically with the flux bias induced by the DC current applied to the pump line. In the following measurement, we fix the DC flux bias to $0.385{{{\Phi }}_0}$, where ${{{\Phi }}_0}$ is the magnetic flux quantum. At this flux bias, the resonance frequency ${\omega _{\text{r}}}/2\pi = 8.9653{\text{ GHz}}$, the external loss rate ${\kappa _{\text{e}}}/2\pi = 0.27{\text{ MHz}}$ and the internal loss rate ${\kappa _{\text{i}}}/2\pi = 0.45{\text{ MHz}}$ are estimated by fitting the measured $\Gamma$ to equation (1). Note that the internal loss rate obtained from this reflection measurement may include the contribution from pure dephasing; the actual internal loss rate $\kappa _{\text{i}}^*$ is related to ${\kappa _{\text{i}}}$ as $\kappa _{\text{i}}^* = {\kappa _{\text{i}}} - 2\gamma$, where $\gamma$ is the pure dephasing rate of the KPO.

Next, we measure $\Gamma$ under the two-photon drive using a vector network analyzer (VNA). The two-photon drive generated by an additional microwave source at room temperature is applied to the pump line on the chip (figure 1(b)). The frequency of the two-photon drive ${\omega _{\text{p}}}$ is approximately twice the resonance frequency of the KPO, and the detuning is defined as the difference between the resonance frequency and half of the two-photon-drive frequency, that is, ${{\Delta }} = {\omega _{\text{r}}} - {\omega _{\text{p}}}/2$.

The Hamiltonian of the KPO in a frame rotating at ${\omega _{\text{p}}}/2$ is expressed as [17]

$$\begin{align}\frac{\mathcal{H}}{\hbar } = {{\Delta }}{a^\dagger }a - \frac{\chi }{2}{a^\dagger }{a^\dagger }aa + \beta \left( {{a^\dagger }{a^\dagger } + aa} \right),\end{align} \tag{ 2 }$$

where $a$ is the photon annihilation operator, and $\beta$ is the amplitude of the two-photon drive. The Kerr nonlinearity $\chi /2\pi = 17.0{\text{ MHz}}$ is determined by measuring the probe-power dependence of $\Gamma$ [16]. The Fock state $|n\rangle$, where $n$ is the photon number, is the eigenstate with eigenenergy of $n{{\Delta }} - n\left( {n - 1} \right)\chi /2$ when the two-photon drive is absent. To examine situations with varying degeneracy of eigenstates, we use three different ${{\Delta ^{\prime} s}}$; ${{\Delta }} = 2\pi \times 8.20{\text{ MHz }}\left( {\sim}{\chi /2} \right),{\text{ }}2{{\pi }} \times 0.05{\text{ MHz }}( {<}{\kappa _{{\text{total}}}}{\text{ }}),{\text{ and }}2\pi \times \left( { - 8.10} \right)$ ${\text{ MHz }}\left( {\sim}{-{{\chi }}/2} \right)$, where ${\kappa _{{\text{total}}}}\left( { = {\kappa _{\text{e}}} + {\kappa _{\text{i}}}} \right)$ is the total photon loss rate. If ${{\Delta }} = {{\chi }}/2{\text{ }}\left( {{{\Delta }} = 0} \right)$, the Fock states $\left| 0 \right\rangle$ and $\left| 2 \right\rangle$ ($\left| 1 \right\rangle$) are degenerate at $\beta = 0$. In the case of ${{\Delta }} = - {{\chi }}/2$, there are no degenerate states at $\beta = 0$.

Figure 2 shows the amplitude of $\Gamma$ as a function of ${\omega _{{\text{ref}}}}$ and the two-photon drive power at the generator output ${\bar p_{RT}}$ for different detunings [(a)${{\Delta }}/2\pi = + 8.20{\text{ MHz}}$, (b)$+ 0.05{\text{ MHz}}$ , (c)$- 8.10{\text{ MHz }}$]. Here, ${\bar p_{RT}}$ was changed for each scan of ${\omega _{{\text{ref}}}}$ using the VNA. In all cases, we observed single dips in $\left| \Gamma \right|$ corresponding to the transition from the state $\left| 0 \right\rangle$ to the state $\left| 1 \right\rangle$ for sufficiently small ${\bar p_{RT}}$. As ${\bar p_{RT}}$ is increased, the spectra exhibit other transitions, and the patterns and the number of transitions [6, 2, and 3 for (a), (b), and (c), respectively] vary with the detuning. In addition, amplification peaks ($\left| \Gamma \right| > 1$) are observed in figures 2(a) and (c), while they are absent in figure 2(b). For ${{\Delta }}/2\pi = + 8.20{\text{ MHz}}$, the transition observed at $\left( {{\omega _{{\text{ref}}}} - {\omega _{\text{p}}}/2} \right)/2\pi = + 8.20{\text{ MHz}}$ is an absorption dip at small ${\bar p_{RT}}$ and splits into two amplification peaks as ${\bar p_{RT}}$ is increased. In the region of ${\bar p_{RT}} > - 30{\text{ dBm}}$, three of the six transitions shift toward higher frequency, and the others toward lower frequency as ${\bar p_{RT}}$ is increased. For ${{\Delta }}/2\pi = + 0.05{\text{ MHz}}$, two transitions are observed, and the frequency of the transition observed at $\left( {{\omega _{{\text{ref}}}} - {\omega _{\text{p}}}/2} \right)/2\pi = 0.05{\text{ MHz}}$ is constant and that of the other transition becomes lower in frequency as ${\bar p_{RT}}$ is increased. For ${{\Delta }}/2\pi = - 8.10{\text{ MHz}}$, one of the three transitions is observed as an amplification peak, and the others are absorption dips. One of the absorption dips shifts toward higher frequency, and the other toward lower frequency as ${\bar p_{RT}}$ is increased.

Zoom In Zoom Out Reset image size

We numerically simulate the reflection coefficient for comparison with the experiments. The experimentally determined frequency detuning ${{\Delta }}$, Kerr nonlinearity $\chi$, and internal (external) photon loss rate ${\kappa _{\text{i}}}$ (${\kappa _{\text{e}}}$) were used for the simulations. The calculations are performed using Python and QuTiP package [30, 31]. The reflection coefficients are calculated following [28] as

$$\begin{align}{{\Gamma }} = 1 + \mathop \sum \limits_{\tilde m,\;\tilde n} {\xi _{\tilde m\tilde n}}\end{align} \tag{ 3 }$$

with

$$\begin{align}{\xi _{\tilde m\tilde n}} = \frac{{{\kappa _{\text{e}}}{{\left| {{X_{\tilde m\tilde n}}} \right|}^2}\left( {{\rho _{\tilde m\tilde m}} - {\rho _{\tilde n\tilde n}}} \right)}}{{i\left( {{\omega _{{\text{ref}}}} - {\omega _{\tilde m\tilde n}}} \right)\; + \;\frac{{{\kappa _{\text{e}}} + {\kappa _{\text{i}}}}}{2}\left( {{Y_{\tilde m\tilde m}} + {Y_{\tilde n\tilde n}}} \right)}}.\end{align} \tag{ 4 }$$

Here, ${\xi _{\tilde m\tilde n}}$ is the contribution of the transition from the state $\left| {\tilde n} \right\rangle$ to the state $\left| {\tilde m} \right\rangle$ in the reflection coefficient. The state $\left| {\tilde n} \right\rangle$ is the energy eigenstate of the KPO, and the state is labeled by the photon number of the state when the two-photon drive is absent. Here, ${\omega _{\tilde m\tilde n}}$ is the transition frequency from the state $\left| {\tilde n} \right\rangle$ to the state $\left| {\tilde m} \right\rangle$ calculated by diagonalizing the Hamiltonian equation (2). ${\rho _{\tilde n\tilde n}}$ is the $\tilde n$th diagonal element $\left\langle {\tilde n} \right|\rho \left| {\tilde n} \right\rangle$ of the steady-state density matrix obtained by solving a Gorini–Kossakowski–Sudarshan–Lindblad (GKSL) master equation with a single-photon loss. Here, ${\kappa _{{\text{total}}}}$ is used as the single-photon loss rate. ${X_{\tilde m\tilde n}} = \left\langle {\tilde m} \right|{a^\dagger }\left| {\tilde n} \right\rangle$ is the matrix element of the one-photon transition, and ${Y_{\tilde m\tilde m}} = \left\langle {\tilde m} \right|{a^\dagger }a\left| {\tilde m} \right\rangle$ is the expectation value of the photon number for the state $\left| {\tilde m} \right\rangle$. Fock states with a photon number up to 30 are considered in the calculation.

Figures 3(a)–(c) show the result of the numerical simulations, which are qualitatively in agreement with the experimental results shown in figures 2(a)–(c). For quantitative comparison in the following sections, here we convert the horizontal axes in figure 2, the power of the two-photon drive ${\bar p_{RT}}$, into the amplitude of the two-photon drive $\beta$, which can be calculated by the first-order derivative of resonance frequency with respect to the current applied to the pump line $d{\omega _{\text{r}}}/di$ (appendix A):

$$\begin{align}\beta = \frac{1}{4}\left| {\frac{{d\omega }}{{di}}} \right|\sqrt {\frac{2}{{{Z_0} \times 1000}}} {10^{\frac{{{{\bar p}_{{\text{KPO}}}}{\text{ }}}}{{20}}}}\end{align} \tag{ 5 }$$

Zoom In Zoom Out Reset image size

From the measurement results shown in figure 1(c), we obtained $d\omega /di = 2\pi \times \left( { - 8.27} \right){\text{ MHz}}/\mu {\text{A}}$ at this bias point. Another parameter, the two-photon drive power at the chip ${\bar p_{{\text{KPO}}}}$, is ${\bar p_{\text{KPO}}} = {\bar p_{{\text{RT}}}} + {R_{{\text{RT}} - {\text{KPO}}}}$, where ${R_{{\text{RT}} - {\text{KPO}}}}$ is the attenuation of the two-photon drive inside and outside the refrigerator. The attenuation ${R_{{\text{RT}} - {\text{KPO}}}}$ is estimated by comparing the measured transition frequencies with the theoretical ones, which were calculated using the Hamiltonian in the absence of a coherent drive under the rotating wave approximation. If a transition frequency is separated from others by a frequency larger than the line width, the transition frequency is extracted from the measured spectrum (figures 2(a)–(c)) by fitting it to the following equation [28]

$$\begin{align}{\Gamma _{\tilde m\tilde n}} = 1 - \frac{{{{\kappa }}_{\text{e}}^{\left( {\tilde m\tilde n} \right)}}}{{i\left( {{\omega _{{\text{ref}}}} - {\omega _{\tilde m\tilde n}}} \right) + \frac{{\kappa _{\text{e}}^{\left( {\tilde m\tilde n} \right)} + \kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}}}{2}}},\end{align} \tag{ 6 }$$

where $\kappa _{\text{e}}^{\left( {\tilde m\tilde n} \right)}$ and $\kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}$ correspond to nominal external and internal loss rates, respectively. To extract ${R_{{\text{RT}} - {\text{KPO}}}}$, we used ${R_{{\text{RT}} - {\text{KPO}}}}$ as a fitting parameter minimizing the square of the difference between the measured and calculated transition frequencies. The estimated values of ${R_{{\text{RT}} - {\text{KPO}}}}$ are −57.0 dB, −57.6 dB, and −57.6 dB with ±0.1 dB uncertainties for ${{\Delta }}/2\pi = + 8.20\;{\text{MHz}},\; + 0.05\;{\text{MHz}},\;{\text{and}}\; - 8.10\;{\text{MHz}},$ respectively. From an independent measurement using a chip with a through transmission line, we estimate that the two-photon drive is attenuated by 58 dB from the microwave source to the KPO with an accuracy of ±2 dB, which is consistent with the above result.

Figures 3(d)–(f) show the same experimental data as figure 2 plotted as a function of $\beta$ in the horizontal axis. The theoretical calculation of the transition frequencies ${\omega _{\tilde m\tilde n}}$'s is also plotted by gray dashed lines. The measured frequencies agree with the calculated transition frequencies within the error of 0.90 MHz, 0.29 MHz, 0.95 MHz for ${{\Delta }}/2\pi = + 8.20\;{\text{MHz}},\; + 0.05\;{\text{MHz}},\;{\text{and}}\; - 8.10\;{\text{MHz}},$ respectively.

The experimental results (figures 3(d)–(f)) are highly consistent with the numerical calculations (figures 3(a)–(c)). The agreement indicates that the theoretically predicted energy levels and the population distribution of the stationary state are reproduced in the experiment. From here, we interpret each transition spectrum shown in figures 3(a)–(c). Note that the two-photon drive does not change the photon number parity, whereas the one-photon drive used as the probe signal does. Therefore, transitions can occur between the states $\left| {2\tilde k} \right\rangle$ and $\left| {\widetilde {\left( {2k^{\prime} + 1} \right)}} \right\rangle$, where both $k$ and $k^{\prime}$ are non-negative integers. For the transitions $\left| {\tilde n} \right\rangle \to \left| {\tilde m} \right\rangle$ at ${\omega _{\tilde m\tilde n}}$ to be visible, both $\left| {{X_{\tilde m\tilde n}}} \right|$ and $\left( {{\rho _{\tilde m\tilde m}} - {\rho _{\tilde n\tilde n}}} \right)$ must be non-zero (equation (4)). Whether the transition at ${\omega _{\tilde m\tilde n}}$ appears as a peak or dip in the spectrum is determined by the sign of the population difference $\left( {{\rho _{\tilde m\tilde m}} - {\rho _{\tilde n\tilde n}}} \right)$. Thus, the reflection measurement provides information about the eigenenergy, population distribution, and transition matrix element.

First, we consider ${{\Delta }}/2\pi = 8.20\;{\text{MHz}}$. Figures 3(g) and (j) show the numerically calculated energy levels and their populations, respectively, for the four highest energy eigenstates of the KPO in the rotating frame as a function of $\beta$. The energy level and the state population are schematically drawn in figures 4(a) and (c). In the case of $\beta /2\pi = 0.21{\text{ MHz}} \ll \chi /2\pi$ (figure 4(a)), the states $\left| {\tilde 0} \right\rangle$ and $\left| {\tilde 2} \right\rangle$ are almost degenerate. Because $\beta$ is sufficiently smaller than ${\kappa _{{\text{total}}}}$, the population is concentrated in the state $\left| {\tilde 0} \right\rangle$. Therefore, the transition $\left| {\tilde 0} \right\rangle \to \left| {\tilde 1} \right\rangle$ ($\left| {\tilde 1} \right\rangle \to \left| {\tilde 0} \right\rangle$) should appear as an absorption dip (amplification peak). Figure 4(b) plots both experimentally obtained and numerically calculated $\Gamma$'s as a function of ${\omega _{{\text{ref}}}}$. As expected, we observe an absorption dip at ${\omega _{\tilde 1\tilde 0}}$ and a small amplification peak at ${\omega _{\tilde 0\tilde 1}}$ and their difference in magnitude comes from the difference in their transition matrix elements of the one-photon drive. Similarly, the transitions at ${\omega _{\tilde 2\tilde 1}}$ and ${\omega _{\tilde 3\tilde 0}}$ appear as absorption dips. The fact that the transition between the states $\left| {\tilde 0} \right\rangle$ and $\left| {\tilde 3} \right\rangle$ exhibit weak probe power indicates that at least one of them deviates from the Fock states. We note that the discrepancy observed for the dip at ${\omega _{\tilde 2\tilde 1}}$ is possibly due to its proximity to the peak at ${\omega _{\tilde 0\tilde 1}}$. As discussed in [28], when two transitions are close to each other, there is an interference effect, which is neglected in equations (3) and (4).

Zoom In Zoom Out Reset image size

Figures 4(c) and (d) show the schematic energy diagram and the reflection coefficient at $\beta /2\pi = 7.40\;{\text{MHz}}$, respectively. Because the population inversion (${\rho _{\tilde 0\tilde 0}} < {\rho _{\tilde 1\tilde 1}}$) is created, the transition $\left| {\tilde 0} \right\rangle \to \left| {\tilde 1} \right\rangle$ at ${\omega _{\tilde 1\tilde 0}}$ turns out to be an amplification peak, which is referred to as dressed-state amplification [32]. Because the degeneracy of the states $\left| {\tilde 0} \right\rangle$ and $\left| {\tilde 2} \right\rangle$ is lifted at $\beta /2\pi = 7.40\;{\text{MHz}}$ (figure 3(g)), the transition $\left| {\tilde 2} \right\rangle \to \left| {\tilde 1} \right\rangle$ is separated from the transition $\left| {\tilde 0} \right\rangle \to \left| {\tilde 1} \right\rangle$ and becomes visible as another amplification peak at ${\omega _{\tilde 1\tilde 2}}$. The transition $\left| {\tilde 2} \right\rangle \to \left| {\tilde 3} \right\rangle$ is separated from the transition $\left| {\tilde 0} \right\rangle \to \left| {\tilde 3} \right\rangle$ and also becomes visible as another dip at ${\omega _{\tilde 3\tilde 2}}$. Because four states $\left| {\tilde 0} \right\rangle$, $\left| {\tilde 1} \right\rangle$, $\left| {\tilde 2} \right\rangle$, and $\left| {\tilde 3} \right\rangle$ are involved, there are eight possible transitions $\left| {\tilde 0} \right\rangle \leftrightarrow \left| {\tilde 1} \right\rangle$, $\left| {\tilde 0} \right\rangle \leftrightarrow \left| {\tilde 3} \right\rangle$, $\left| {\tilde 2} \right\rangle \leftrightarrow \left| {\tilde 1} \right\rangle$, and $\left| {\tilde 2} \right\rangle \leftrightarrow \left| {\tilde 3} \right\rangle$, which appear as four absorption dips and four amplification peaks. However, the transitions $\left| {\tilde 3} \right\rangle \to \left| {\tilde 0} \right\rangle$ and $\left| {\tilde 3} \right\rangle \to \left| {\tilde 2} \right\rangle$ cannot be observed because the matrix element of the one-photon transition $\left| {{X_{\tilde 0\tilde 3}}} \right|$ and $\left| {{X_{\tilde 2\tilde 3}}} \right|$ are too small as shown in figure 3(m). This makes the number of transitions six. Signal-to-noise ratio in figure 4(d) is lower than that in figure 4(b). This is because the amplitudes of the peaks and dips in figure 4(d) are smaller than those in figure 4(b). The signal-to-noise ratio can be improved by increasing the average number of the VNA.

For $\beta /2\pi > 10{\text{ MHz}}$, some of the transitions disappear as shown in figure 3(a). The peak and dip in the transitions $\left| {\tilde 0} \right\rangle \to \left| {\tilde 1} \right\rangle$ at ${\omega _{\tilde 1\tilde 0}}$ and $\left| {\tilde 1} \right\rangle \to \left| {\tilde 0} \right\rangle$ at ${\omega _{\tilde 0\tilde 1}}$ disappear because the population difference between ${\rho _{\tilde 0\tilde 0}}$ and ${\rho _{\tilde 1\tilde 1}}$ is small (figure 3(j)). Similarly, the dip in the transition $\left| {\tilde 2} \right\rangle \to \left| {\tilde 3} \right\rangle$ at ${\omega _{\tilde 3\tilde 2}}$ disappears because the difference between ${\rho _{\tilde 2\tilde 2}}$ and ${\rho _{\tilde 3\tilde 3}}$ is small. The peak height of the transition $\left| {\tilde 2} \right\rangle \to \left| {\tilde 1} \right\rangle$ at ${\omega _{\tilde 1\tilde 2}}$ becomes small due to decreasing $\left| {{X_{12}}} \right|$ (figure 3(m)).

In the case of ${{\Delta }}/2\pi = + 0.05\;{\text{MHz}}$, the states $\left| {\tilde 0} \right\rangle$ and $\left| {\tilde 1} \right\rangle$ are almost degenerate regardless of the value of $\beta$ as shown in figure 3(h). Within the probe frequency range in this measurement (figure 3(e)), the transitions $\left| {\tilde 0} \right\rangle \leftrightarrow \left| {\tilde 1} \right\rangle$ and $\left| {\tilde 1} \right\rangle \leftrightarrow \left| {\tilde 2} \right\rangle$ can contribute to the spectrum. However, the transition $\left| {\tilde 2} \right\rangle \to \left| {\tilde 1} \right\rangle$ at ${\omega _{\tilde 1\tilde 2}}$ is not visible because the matrix element of the one-photon transition $\left| {{X_{\tilde 1\tilde 2}}} \right|$ is small (figure 3(n)). For the transitions $\left| {\tilde 0} \right\rangle \leftrightarrow \left| {\tilde 1} \right\rangle$, the frequencies ${\omega _{\tilde 1\tilde 0}}$ and ${\omega _{\tilde 0\tilde 1}}$ are the same, so the dip at ${\omega _{\tilde 1\tilde 0}}$ overlaps with the peak at ${\omega _{\tilde 0\tilde 1}}$, resulting in two distinguishable spectra. As the populations of ${\rho _{\tilde 0\tilde 0}}$ and ${\rho _{\tilde 1\tilde 1}}$ approach 0.5 (figure 3(k)), the dip depth at ${\omega _{\tilde 1\tilde 0}}$ appears to decrease (figures 3(b) and (e)). In contrast, although the matrix element of the one-photon transition $\left| {{X_{\tilde 2\tilde 1}}} \right|$ slightly decreases (figure 3(n)), the peak of $\left| {\tilde 1} \right\rangle \to \left| {\tilde 2} \right\rangle$ at ${\omega _{\tilde 2\tilde 1}}$ becomes more visible, because the population difference between ${\rho _{\tilde 1\tilde 1}}$ and ${\rho _{\tilde 2\tilde 2}}$ increases (figure 3(k)). For $\beta \gg \chi$, the steady state of the states $\left| {\tilde 0} \right\rangle$ and $\left| {\tilde 1} \right\rangle$ are the same mixed state [10, 11].

For ${{\Delta }}/2\pi = - 8.10\;{\text{MHz}}$, there are no degenerate states at $\beta = 0$ (figure 3(i)). We observe one amplification peak at ${\omega _{\tilde 0\tilde 1}}$ and two absorption dips at ${\omega _{\tilde 1\tilde 0}}$ and ${\omega _{\tilde 2\tilde 0}}$, which can be attributed to the relationship between the populations ${\rho _{\tilde 0\tilde 0}} > {\rho _{\tilde 1\tilde 1}} > {\rho _{\tilde 2\tilde 2}}$ (figure 3(l)). As the population ${\rho _{\tilde 1\tilde 1}}$ increases, the transition $\left| {\tilde 1} \right\rangle \to \left| {\tilde 2} \right\rangle$ become more visible as the dip at ${\omega _{\tilde 2\tilde 1}}$ in $\left| \Gamma \right|$ (figures 3(c) and (f)). In the same regime, although the population difference between ${\rho _{\tilde 0\tilde 0}}$ and ${\rho _{\tilde 1\tilde 1}}$ decreases, the peak at ${\omega _{\tilde 0\tilde 1}}$ become visible because $\left| {{X_{\tilde 0\tilde 1}}} \right|$ increases (figure 3(o)).

The transition frequencies described thus far have been calculated by diagonalizing the Hamiltonian equation (2), and are in good agreement with the experimental data. As shown in appendix B, however, some of the transition-frequency shifts for small $\beta$ can be calculated analytically when ${{\Delta }}\sim \chi /2$, in which the frequency of the two-photon drive is in resonance with the energy difference between the states $\left| 0 \right\rangle$ and $\left| 2 \right\rangle$. Namely, the frequency shifts of ${\omega _{\tilde 1\tilde 0}},\;{\omega _{\tilde 0\tilde 1}},\;{\omega _{\tilde 1\tilde 2}},\;{\text{and}}\;{\omega _{\tilde 2\tilde 1}}$ when ${{\Delta }} = \chi /2$ are given by $+ \left( {\frac{\chi }{2} + \frac{{3{\beta ^2}}}{\chi }} \right) - \sqrt 2 \beta ,$ $- \left( {\frac{\chi }{2} + \frac{{3{\beta ^2}}}{\chi }} \right) + \sqrt 2 \beta , + \left( {\frac{\chi }{2} + \frac{{3{\beta ^2}}}{\chi }} \right) + \sqrt 2 \beta$ and $- \left( {\frac{\chi }{2} + \frac{{3{\beta ^2}}}{\chi }} \right) - \sqrt 2 \beta$, respectively, where the term $\frac{{3{\beta ^2}}}{\chi }$ is caused by the Stark shift of the state $\left| 1 \right\rangle$ and the term $\sqrt 2 \beta$ results from the Rabi splitting originating from the states $\left| {0} \right\rangle$ and $\left| {2} \right\rangle$. Figure 5 shows the comparison of transition frequencies obtained from diagonalizing the Hamiltonian when ${{\Delta }} = 2\pi \, \times \,8.20\,{\text{MHz}}\sim \chi /2$ (gray lines) and analytical formula (dashed lines). They are in good agreement for small $\beta$. In the range of $\beta < 2\pi \; \times \;2\;{\text{MHz}}$, the difference between them is less than 0.38 MHz, which is smaller than ${\kappa _{{\text{total}}}}$. However, a discrepancy can be seen in the range of $\beta > 2\pi \; \times \;2\;{\text{MHz}}$ especially for ${\omega _{\tilde 1\tilde 0}}\;{\text{and}}\;{\omega _{\tilde 0\tilde 1}}$. This is because the contribution of the states with large photon numbers cannot be ignored when $\beta$ is large. Namely, the analytical formula is derived by considering only the Fock states with the photon number up to 2. However, as $\beta$ increases, the population of the Fock states with larger photon numbers increases and the states $\left| {\tilde 0} \right\rangle$, $\left| {\tilde 1} \right\rangle$, and $\left| {\tilde 2} \right\rangle$ can no longer be well approximated as $\left| {{ \pm _2}} \right\rangle$ and $\left| {1^{\prime}} \right\rangle$, which are analytically derived in the subspace spanned by $\left| 0 \right\rangle$, $\left| 1 \right\rangle$, and $\left| 2 \right\rangle$ (see appendix B).

Zoom In Zoom Out Reset image size

Thus far, we have investigated the $\beta$ dependence of the frequency and the amplitude of peaks and dips in $\left| \Gamma \right|$ and found good agreement between the experiment and the theory as shown in figures 3(a)–(f). However, the overall linewidth of the peaks and dips in the experiments is smaller than the theoretical as seen in the figures, where it should be noted that the color scales are all the same. To examine this difference, we investigate the nominal loss rates $\kappa _{\text{e}}^{\left( {\tilde m\tilde n} \right)}$ and $\kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}$ obtained by fitting the measured $\left| \Gamma \right|$ to ${\Gamma _{mn}}$ in equation (6).

By comparing equations (3) and (6), the nominal external loss rate $\tilde \kappa _{\text{e}}^{\left( {\tilde m\tilde n} \right)}$ and internal loss rate $\tilde \kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}$ are given by

$$\begin{align}\tilde \kappa _{\text{e}}^{\left( {\tilde m\tilde n} \right)}& = - {\kappa _{\text{e}}}{\left| {{X_{\tilde m\tilde n}}} \right|^2}\left( {{\rho _{\tilde m\tilde m}} - {\rho _{\tilde n\tilde n}}} \right),\end{align} \tag{ 7 }$$

$$\begin{align}\tilde \kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}& = \left( {{\kappa _{\text{e}}} + {\kappa _{\text{i}}}} \right)\left( {{Y_{\tilde m\tilde m}} + {Y_{\tilde n\tilde n}}} \right) + {\kappa _{\text{e}}}{\left| {{X_{\tilde m\tilde n}}} \right|^2}\left( {{\rho _{\tilde m\tilde m}} - {\rho _{\tilde n\tilde n}}} \right).\end{align} \tag{ 8 }$$

In figure 6, we plot $\kappa _{\text{e}}^{\left( {\tilde m\tilde n} \right)}$ and $\kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}$ for ${{\Delta }}/2\pi = + 8.20\;{\text{MHz}}$ as a function of $\beta$. For comparison, we also plot those predicted from the model, $\tilde \kappa _{\text{e}}^{\left( {\tilde m\tilde n} \right)}$ and $\tilde \kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}$ (equations (7) and (8)). As shown in figure 6(a), $\kappa _{\text{e}}^{\left( {\tilde m\tilde n} \right)}$ agrees well with $\tilde \kappa _{\text{e}}^{\left( {\tilde m\tilde n} \right)}$ for all $\left( {\tilde m,\;\tilde n} \right)$, except that $\kappa _{\text{e}}^{\left( {\tilde 3\tilde 0} \right)}$ differs significantly from $\tilde \kappa _{\text{e}}^{\left( {\tilde 3\tilde 0} \right)}$ in the range of $\beta /2\pi {\text{ }}$ from 0.4 MHz to 4 MHz. This deviation is due to the fact that the dips in the transitions $\left| {\tilde 0} \right\rangle \to \left| {\tilde 3} \right\rangle$ at ${\omega _{\tilde 3\tilde 0}}$ and $\left| {\tilde 2} \right\rangle \to \left| {\tilde 3} \right\rangle$ at ${\omega _{\tilde 3\tilde 2}}$ overlap with each other in this pump power range (figure 3(d)), so that there is the interference effect between the transitions as mentioned earlier in figure 4(b). Although the frequency difference between ${\omega _{\tilde 3\tilde 0}}$ and ${\omega _{\tilde 3\tilde 2}}$ is also small in the region of $\beta /2\pi < 0.4\;{\text{MHz}}$, the deviation is not seen because the population difference between ${\rho _{\tilde 2\tilde 2}}$ and ${\rho _{\tilde 0\tilde 0}}$ makes the transition $\left| {\tilde 0} \right\rangle \to \left| {\tilde 3} \right\rangle$ dominant (figure 3(j)).

Zoom In Zoom Out Reset image size

In contrast, the agreement between $\kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}$ and $\tilde \kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}$ is poor as shown in figure 6(b). Only $\kappa _{\text{i}}^{\left( {\tilde 1\tilde 0} \right)}$ and $\kappa _{\text{i}}^{\left( {\tilde 3\tilde 0} \right)}$ are consistent with the theory in the limited range of $\beta / 2\pi < 0.1\;{\text{MHz}}$. Overall, $\tilde \kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}$ is smaller than the theoretical estimation $\kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}$, except for $\kappa _{\text{i}}^{\left( {\tilde 3\tilde 0} \right)}$ in the range of $\beta /2\pi {\text{ }}$ from 0.4 MHz to 4 MHz, where we observed deviation in $\kappa _{\text{e}}^{\left( {\tilde 3\tilde 0} \right)}$ as discussed above.

The discrepancy between $\tilde \kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}$ and $\kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}$ is not well understood at this time. Because the nominal external photon loss rate, which is related to the single-photon decay to the signal line, is consistent with the theory, we believe that dephasing or losses other than the single-photon loss can cause the discrepancy between $\tilde \kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}$ and $\kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}$. We do not include the phase relaxation effect in equations (4) and (6), but we found that simply including phase relaxation in the master equation by $\mathcal{L} = {a^\dagger }a$ as the Lindblad operator does not reproduce $\tilde \kappa _{\text{i}}^{\left( {\tilde m\tilde n} \right)}$. One possibility is the 1/f spectrum of flux noise [33, 34], which is not taken into account in the simulation. The phase relaxation and noise spectra will need to be further investigated in future studies.

We performed spectroscopic measurements of a KPO under the two-photon drive induced by the magnetic flux modulation and compared the results with theoretical calculations. Since we use the continuous measurement of the reflection coefficient of the KPO, our method is simpler than previous studies [21, 22] for charge-driven SNAIL transmons with pulse operations. The transition frequencies are in good agreement with the calculations after adjusting the attenuation values for the microwave transmission line, which turned out to be consistent with the results of independent measurements. We also demonstrated that some transition frequencies can be interpreted as Rabi splitting and the Stark shift and provided a simple analytical formula for them. The magnitude of the reflection coefficients is also in good agreement with the calculation, and we showed that its behavior under the two-photon drive, such as the change from the absorption dip to amplification peak, can be explained by the difference in the population between the initial and final states of the transition. These results demonstrate that the system is well characterized by the model Hamiltonian including the two-photon drive induced by the flux modulation. We also investigated the nominal loss rates and found that the external loss is in good agreement with the theoretical calculation, while the internal loss rate shows significant deviation from the theory, which still needs clarification.

The spectroscopic results have potential applications in quantum annealing machines and gate-model quantum computers. The estimation of the magnitude of the two-photon drive at the quantum chip is critical for efficient parameter scheduling [13, 26, 27]. Additionally, the accurate estimation of KPO eigenenergies are useful for gate operations of a cat qubit [23].

We thank Y Kitagawa for his assistance in the device fabrication. We also thank T Aoki, K Koshino, and Y Kano for their fruitful discussions. A part of this work was conducted at the AIST Nano-Processing Facility supported by Nanotechnology Platform Program of the Ministry of Education, Culture, Sports, Science and Technology (MEXT), Japan. The devices were fabricated in the Superconducting Quantum Circuit Fabrication Facility (Qufab) in the National Institute of Advanced Industrial Science and Technology (AIST). This paper is based on the results obtained from a Project, JPNP16007, commissioned by the New Energy and Industrial Technology Development Organization (NEDO). This work was partially supported by JST Moonshot R&D (Grant Number JPMJMS226C).

The data cannot be made publicly available upon publication because no suitable repository exists for hosting data in this field of study. The data that support the findings of this study are available upon reasonable request from the authors.

In this appendix, we derive the formula for the amplitude of the two-photon drive applied to the KPO (equation (5)) based on its relation to the modulation amplitude of the resonance frequency.

First, we derive the relation between the amplitude of the frequency modulation $\delta \omega$ and the power of the two-photon drive. We assume the time dependence of the resonance frequency as ${\omega _{\text{r}}}\left( t \right) = {\omega _{\text{r}}} + \delta \omega \cos {\omega _p}t$, where ${\omega _{\text{r}}}$ is the static resonance frequency and ${\omega _p}$ is the frequency of the two-photon drive. Then, $\delta \omega$ is related to the derivative of the resonance frequency with respect to the bias current $\frac{{d{\omega _{\text{r}}}}}{{di}}$ as $\delta \omega = \frac{{d{\omega _{\text{r}}}{\text{ }}}}{{di}}{i_{{\text{amp}}}}$, where ${i_{{\text{amp}}}}$ is the amplitude of the ac current induced by the two-photon drive. Here, we assume that the mutual inductance between the pump line and the SQUID is frequency independent. Because the root mean square value of the current $\left| i \right|$ is $\left| i \right| = 1/\sqrt 2 {i_{{\text{amp}}}}$, the average power $\bar p$ is given by $\bar p = {Z_0}{\left| i \right|^2} = \frac{{{Z_0}i_{{\text{amp}}}^2}}{2}$, where ${Z_0}$ is a characteristic impedance of the pump line. Therefore, $\bar p$ (in dBm) can be written as $\bar p = 10{\log _{10}}\left[ {\frac{{{Z_0}i_{{\text{amp}}}^2 \times 1000}}{2}} \right]$, i.e. ${i_{{\text{amp}}}} = \sqrt {\frac{2}{{{Z_0} \times 1000}}} {10^{\frac{{\bar p{\text{ }}}}{{20}}}}$, and $\delta \omega$ is written as

$$\begin{align}\delta \omega = \left| {\frac{{d\omega }}{{di}}} \right|\sqrt {\frac{2}{{{Z_0} \times 1000}}} {\text{ }}{10^{\frac{{\bar p{\text{ }}}}{{20}}}}.\end{align} \tag{ A1 }$$

The relation between the two-photon drive $\beta$ and frequency-modulation amplitude can be determined from the Hamiltonian equation (2) in a laboratory frame,

$$\begin{align}\frac{\mathcal{H}}{\hbar } = {{{\omega }}_0}{a^\dagger }a - \frac{{{\chi }}}{{12}}{\left( {{a^\dagger } + a} \right)^4} + 2\beta \cos \left( {{\omega _{\text{p}}}t} \right){\left( {{a^\dagger } + a} \right)^2}.\end{align} \tag{ A2 }$$

By expanding and rearranging the right-hand side of this equation into normally ordered products, the coefficient of ${a^\dagger }a$ becomes ${\omega _0} - \chi + 4\beta \cos \left( {{\omega _{\text{p}}}t} \right)$. The constant terms are the resonance frequency ${\omega _{\text{r}}} = {\omega _0} - \chi$ of the KPO without the two-photon drive, and the oscillating term represents the frequency modulation caused by the two-photon drive. Therefore, the frequency-modulation amplitude $\delta \omega$ can be written as $\delta \omega = 4\beta$, which together with equation (A1) leads to equation (5).

In this appendix, we derive formulae for four of the transition frequencies as a function of the two-photon drive amplitude $\beta$, which are discussed in section 6 and shown in figure 5. We investigate the physical meaning of the frequency shift by deriving its formula analytically.

We consider the situation where ${\omega _{\text{p}}}$ matches the energy difference between $\left| 0 \right\rangle$ and $\left| 2 \right\rangle$. In this case, the two-photon drive acts as a resonant drive for $\left| 0 \right\rangle$ and $\left| 2 \right\rangle$ and as a detuned drive for the states with two different photon numbers, e.g. for $\left| 1 \right\rangle$ and $\left| 3 \right\rangle$ (figure 7(a)). It is well known that a resonant drive applied to a two-level system causes Rabi splitting and a non-resonant drive causes the Stark shift [2]; here we extend the same idea to the two-photon drive applied to the KPO.

Zoom In Zoom Out Reset image size

If the states with a photon number of less than 3 are considered, the Hamiltonian equation (2) is rewritten as

$$\begin{align} \frac{{{\mathcal{H}_{\left( {0 - 3} \right)}}}}{\hbar } = \sqrt 2 \beta \left( {\left| 0 \right\rangle \left\langle 2 \right| + \left| 2 \right\rangle \left\langle 0 \right|} \right) + \sqrt 6 \beta \left( {\left| 1 \right\rangle \left\langle 3 \right| + \left| 3 \right\rangle \left\langle 1 \right|} \right) + \frac{1}{2}\chi \left| 1 \right\rangle \left\langle 1 \right| - \frac{3}{2}\chi \left| 3 \right\rangle \left\langle 3 \right|.\end{align} \tag{ A3 }$$

We decompose ${\mathcal{H}_{{\text{RWA}}}}$ as

$$\begin{align}{\mathcal{H}_{\left( {0 - 3} \right)}} = {\mathcal{H}_{{\text{Rabi}}}} + {\mathcal{H}_{{\text{Stark}}}},\end{align} \tag{ A4 }$$

where

$$\begin{align}\frac{{{\mathcal{H}_{{\text{Rabi}}}}}}{\hbar }& = \sqrt 2 \beta \left( {\left| 0 \right\rangle \left\langle 2 \right| + \left| 2 \right\rangle \left\langle 0 \right|} \right),\end{align} \tag{ A5 }$$

$$\begin{align}\frac{{{\mathcal{H}_{{\text{Stark}}}}}}{\hbar }& = \sqrt 6 \beta \left( {\left| 1 \right\rangle \left\langle 3 \right| + \left| 3 \right\rangle \left\langle 1 \right|} \right) + \frac{1}{2}\chi \left| 1 \right\rangle \left\langle 1 \right| - \frac{3}{2}\chi \left| 3 \right\rangle \left\langle 3 \right|\end{align} \tag{ A6 }$$

${\mathcal{H}_{{\text{Rabi}}}}$ represents a resonant drive for $\left| 0 \right\rangle$ and $\left| 2 \right\rangle$ states, and ${\mathcal{H}_{{\text{Stark}}}}$ represents a detuned drive for $\left| 1 \right\rangle$ and $\left| 3 \right\rangle$ states (figure 7(a)).

First, we focus on ${\mathcal{H}_{{\text{Rabi}}}}$. This term is described as a two-level system of $\left| 0 \right\rangle$ and $\left| 2 \right\rangle$. The diagonalization of ${\mathcal{H}_{{\text{Rabi}}}}$ yields

$$\begin{align} \frac{{{\mathcal{H}_{{\text{Rabi}}}}}}{\hbar } \simeq \sqrt 2 \beta \left( {\left| {{ + _2}} \right\rangle \left\langle {{ + _2}} \right| - \left| {{ - _2}} \right\rangle \left\langle {{ - _2}} \right|} \right).\end{align} \tag{ A7 }$$

Here, we define $\left| {{ \pm _2}} \right\rangle = \frac{1}{{\sqrt 2 }}\left( {|2\rangle \pm |0\rangle } \right)$. This is analogous to the derivation of the Rabi splitting for two-level systems, which can be observed as a Mollow triplet.

Next, we focus on ${\mathcal{H}_{{\text{Stark}}}}$. This term is described as a two-level system of $\left| 1 \right\rangle$ and $\left| 3 \right\rangle$ with a detuned drive. The diagonalization of ${\mathcal{H}_{{\text{Stark}}}}$ yields

$$\begin{align}\frac{{{\mathcal{H}_{{\text{Stark}}}}}}{\hbar } = \left( { - \frac{1}{2}\chi + \chi \sqrt {1 + \frac{{6{\beta ^2}}}{{{\chi ^2}}}} } \right)\left| {1^{\prime}} \right\rangle \left\langle {1^{\prime}} \right| + \left( { - \frac{1}{2}\chi - \chi \sqrt {1 + \frac{{6{\beta ^2}}}{{{\chi ^2}}}} } \right)\left| {3^{\prime}} \right\rangle \left\langle {3^{\prime}} \right|,\end{align} \tag{ A8 }$$

where

$$\begin{align}\left| {1^{\prime}} \right\rangle = \frac{1}{{\sqrt {1 + {\delta ^2}} }}\left( {\left| 1 \right\rangle - \delta \left| 3 \right\rangle } \right),\end{align} \tag{ A9 }$$

$$\begin{align}\left| {3^{\prime}} \right\rangle = \frac{1}{{\sqrt {1 + {\delta ^2}} }}\;\left( {\delta \left| 1 \right\rangle + \left| 3 \right\rangle } \right).\end{align} \tag{ A10 }$$

Here, we define $\delta = \frac{{\chi - \sqrt {{\chi ^2} + 6{\beta ^2}} }}{{\sqrt 6 \beta }}$. If the amplitude of the two-photon drive is much smaller than the Kerr nonlinearity, i.e. ${\left( {\frac{\beta }{\chi }} \right)^2} \ll 1$, equation (A8) is written as

$$\begin{align}\frac{{{\mathcal{H}_{{\text{Stark}}}}}}{\hbar } \simeq \left( {\frac{1}{2}\chi + \frac{{3{\beta ^2}}}{\chi }} \right)\left| {1^{\prime}} \right\rangle \left\langle {1^{\prime}} \right| + \left( { - \frac{3}{2}\chi - \frac{{3{\beta ^2}}}{\chi }} \right)\left| {3^{\prime}} \right\rangle \left\langle {3^{\prime}} \right|.\end{align} \tag{ A11 }$$

Comparing equations (A6) and (A11), it can be seen that the ac Stark shift is caused by the two-photon drive. Namely, the eigenenergies of $\left| 1 \right\rangle$ and $\left| 3 \right\rangle$ are shifted $+ \frac{{3{\beta ^2}}}{\chi }$ and $- \frac{{3{\beta ^2}}}{\chi }$, respectively. From equations (A7) and (A11), ${\mathcal{H}_{{\text{RWA}}}}$ can be written as

$$\begin{align} \frac{{{\mathcal{H}_{\left( {0 - 3} \right)}}}}{\hbar } \simeq \sqrt 2 \beta \left({\left|{ + _2}\right\rangle\left\langle{+ _2}\right|-\left|{ - _2}\right\rangle \left\langle{ - _2}\right|} \right) + \left( {\frac{1}{2}\chi + \frac{{3{\beta ^2}}}{\chi}}\right)\;\left|{1^{\prime}}\right\rangle \left\langle {1^{\prime}} \right| + \left( { - \frac{3}{2}\chi- \frac{{3{\beta ^2}}}{\chi }}\right)\left|{3^{\prime}}\right\rangle \left\langle {3^{\prime}}\right|.\end{align} \tag{ A12 }$$

This rewriting from equation (A3) to equation (A12) corresponds to the rewriting of the Fock states to the dressed states shown as the orange and pink levels in figure 7(b), which are induced by the two-photon drive shown by the vertical arrows in their corresponding colors in figure 7(a).

Transition frequencies between the dressed states can be calculated from figure 7(b). Namely, frequencies for the transitions $\left| { \pm_{2} ,\;N} \right\rangle \leftrightarrow \left| {1^{\prime},\;N+1} \right\rangle$ (dark blue dashed arrows) are $\frac{{{\omega _p}}}{2} + \left( {\frac{\chi }{2} + \frac{{3{\beta ^2}}}{\chi }} \right) \pm \sqrt 2 \beta$ and those for the transitions $\left| {1^{\prime},\;N} \right\rangle \leftrightarrow \left| { \pm_2 ,\;N+1} \right\rangle$ (dark red dashed arrows) are $\frac{{{\omega _p}}}{2} - \left( {\frac{\chi }{2} + \frac{{3{\beta ^2}}}{\chi }} \right) \pm \sqrt 2 \beta$. They are plotted in figure 5 in their corresponding colors.