Demonstration of a non-Abelian geometric controlled-Not gate in a superconducting circuit

Holonomies, arising from non-Abelian geometric transformations of quantum states in Hilbert space, offer a promising way for quantum computation. These holonomies are not commutable and thus can be used for the realization of a universal set of quantum logic gates, where the global geometric feature may result in some noise-resilient advantages. Here we report the first on-chip realization of a non-Abelian geometric controlled-Not gate in a superconducting circuit, which is a building block for constructing a holonomic quantum computer. The conditional dynamics is achieved in an all-to-all connected architecture involving multiple frequency-tunable superconducting qubits controllably coupled to a resonator; a holonomic gate between any two qubits can be implemented by tuning their frequencies on resonance with the resonator and applying a two-tone drive to one of them. This gate represents an important step towards the all-geometric realization of scalable quantum computation on a superconducting platform.


INTRODUCTION
When a nondegenerate quantum system makes a cyclic evolution in the Hilbert space, it will pick up a phase, which, in general, is contributed by both the dynamical and geometric effects. The dynamical part is the time integral of the energy, while the geometric one depends upon the area enclosed by the loop that the quantum state traverses in the Hilbert space. This effect, discovered by Berry in cyclic and adiabatic evolutions [1], has been generalized to nonadiabatic [2] and nondegenerate [3] cases. If a system has degenerate energy levels, the cyclic evolution of the corresponding degenerate subspaces will produce a matrixvalued quantum state transformation that is path-dependent and referred to as non-Abelian geometric phase or holonomy [3]. The Berry phase and holonony depend upon the global geometry of the associated loops and have intrinsic resistance to certain kinds of small errors, suggesting quantum gates based on geometric operations have practical advantages as compared to dynamical gates [4][5][6][7]. In particular, it was shown that all of the elementary one-and two-qubit gates needed for accomplishing any quantum computation task could be achieved with Berry phase and holonomic transformations, offering a possibility for implementations of geometric quantum computation [8,9].
The conditional Berry phase was first observed in nuclear magnetic resonance systems [10]. However, the relatively long operation time associated with an adiabatic evolution represents an unfavorable condition for the implementation of geometric quantum computation with such controlled phase gates. As such, geometric effects without the adiabatic restriction are highly desirable for the implementation of quantum logic gates that are robust against noises [11][12][13][14][15][16]. So far, nonadiabatic geometric controlled-phase gates have been realized in ion traps [17][18][19][20] and superconducting circuits [21][22][23]. On the other hand, Sjöqvist et al. have proposed an approach for realizing a universal set of elementary gates based on nonadiabatic holonomies [24], whose robustness against noises has been analyzed [25,26]. Following this approach, a universal gate set involving two noncommutable single-qubit gates and a two-qubit controlled-Not (CNOT) gate have been experimentally realized with nuclear magnetic resonance [27] and solid-state spins [28,29]. Several groups have demonstrated holonomic single-qubit gates in superconducting circuits [30][31][32][33], which represent a promising platform for quantum computation [34]. Recently, Egger et al.
reported a holonomic operation for producing entangled states in a superconducting circuit [35]. However, a non-Abelian geometric entangling gate necessary for constructing a universal holonomic gate set has not been implemented in such scalable systems. More recently, Han et al. reported a universal set of time-optimal geometric gates with superconducting qubits [36], where single-qubit gates were realized using non-Abelian geometric phase, but the two-qubit gate was based on Abelian geometric phase.
In this paper, we propose and experimentally demonstrate a scheme for realizing non-adiabatic, non-Abelian geometric CNOT gate for two qubits, one acting as the control qubit and the other as the target qubit. The two qubits are strongly coupled to a resonator, so that the energy levels of the target qubit depend on the state of the control qubit. This conditional energy-level shift enables the target qubit to be resonantly driven by classical fields, conditional on the state of the control qubit. With suitable setting of the parameters, these classical fields can drive the degenerate subspace spanned by the two basis states of the target qubit to undergo a conditional cyclic evolution, realizing a CNOT gate between these two qubits. We realize this holonomic gate in a superconducting multi-qubit processor, where any two qubits can be selectively coupled to a common resonator but effectively decoupled from other qubits through frequency tuning. This flexibility enables direct implementation of holonomic gates between any pair of qubits on the chip, without the restriction of nearest-neighbor couplings. The measured process fidelity for the CNOT gate is above 0.9. With further improvements in the device design and fabrication, as confirmed by our numerical simulations, the gate fidelity can be significantly increased. Our scheme is applicable to other spin-boson systems, such as cavity QED and ion traps [37].

THEORETICAL MODEL
The system under consideration is composed of two qutrits coupled to a resonator. Each qutrit has three basis states, as shown in Fig. 1a, with |g and | f serving as two logic states of a qubit, and |e , lying between |g and | f , used as an auxiliary state for realizing the controlled logic operation. For simplicity, we will refer to the qutrits as qubits. As will be shown, the control qubit (Q 1 ) remains in its computational space, while the target qubit (Q 2 ) has a probability of being populated in the auxiliary level |e during the gate operation. The transition |g ←→ |e of each qubit resonantly interacts with the resonator, while | f state is effectively decoupled from the resonator. In the interaction picture, the Hamiltonian describing the qubitresonator interaction is given by where a and a † are the photonic annihilation and creation operators for the resonator, λ j is the coupling strength between the jth qubit and the resonator with angular frequency ω r . We here have set the energy of the ground state |g for each qubit to be 0. To realize the CNOT gate, the transition |g 2 ←→ |e 2 of Q 2 is driven by a classical field with angular frequency (ω r − λ 2 ), and |e 2 ←→ | f 2 is driven by a classical field with angular frequency ω f ,2 − ω r + λ 2 , wherehω f ,2 is the energy of Q 2 's state | f 2 (Fig. 1a). The interaction between the second qubit Revision on 2021.04.24 (ω r − λ 2 ) and ω f ,2 − ω r + λ 2 , respectively. Here ω r is the angular frequency of the resonator that is strongly coupled to |g 2 ←→ |e 2 with the coupling strength λ 2 , andhω f ,2 is the energy spacing between | f 2 and |g 2 . (b) Dressed states and energy levels with Q 1 initially in | f 1 . When being initially in | f 1 , Q 1 is effectively decoupled from the resonator due to the large detuning. The strong coupling between Q 2 and the resonator results in dressed states ψ ± n , whose energy levels are nonlinearly dependent on the coupling strength. The two driving fields are on resonance with the transitions |g 2 0 ←→ ψ − 1 and ψ − 1 ←→ | f 2 0 , respectively, but highly detuned from other transitions. (c) Dressed states and energy levels with Q 1 initially in |g 1 . If initially in |g 1 , Q 1 , together with Q 2 , is strongly coupled to the resonator, resulting in three dressed states Φ ± 1 and Φ 0 1 in the single-excitation subspace. The driving fields are highly detuned from transitions of |g 1 g 2 0 and |g 1 f 2 0 to these dressed states.  Fig. 2. Pulse sequence. Before the gate operation, both qubits are initialized to their ground state at the corresponding idle frequencies, where single-qubit rotations are performed to prepare them in a product state. Then a Z pulse is applied to Q 1 , tuning |g 1 ←→ |e 1 close to the resonator's frequency; Q 2 is subjected to a Z pulse, which brings |g 2 ←→ |e 2 to the resonator's frequency, and a driving pulse involving two frequency components respectively on resonance with the transitions |g 2 0 ←→ ψ − 1 and ψ − 1 ←→ | f 2 0 . After the CNOT gate, realized with these pulses, both qubits are tuned back to their idle frequencies for quantum state tomography. and the driving fields is described by where Ω ge and Ω e f denote the Rabi frequencies of the two fields driving |g 2 ←→ |e 2 and |e 2 ←→ | f 2 , respectively. We here have assumed that the phases of the fields driving the transitions |g ←→ |e and |e ←→ | f are 0 and π, respectively.
Each matrix element is characterized by two colorbars, one for the real part and the other for the imaginary part. The black wire frames denote the matrix elements of the ideal output states.
The strong couplings between the qubits and the resonator produce dressed states, whose energy levels depend on the total excitation number as well as on the number of qubits being initially populated in |g . When the control qubit is in the state | f 1 , it does not interact with the resonator, and the coupling between the target qubit and the resonator is described by the Jaynes-Cummings model, whose eigenstates are given by Here the second symbol in each ket denotes the photon number in the resonator. The eigenenergies of the dressed states ψ ± n areh nω r ± √ nλ 2 . We here consider the case that the resonator is initially in the vacuum state |0 . Consequently, the classical fields resonantly couple the states |g 2 0 and | f 2 0 to the single-excitation dressed state ψ − 1 , respectively, as sketched in Fig. 1b. We suppose that Ω ge and Ω e f are much smaller than λ 2 , so that the classical fields cannot drive the transitions from ψ − 1 to ψ ± 2 due to the large detunings. However, these offresonant couplings shift the energy levels of ψ − 1 by −2hδ 1 , with δ 1 = 2Ω 2 ge /λ 2 (see Supplemental Material). Furthermore, offresonant coupling to |h 1 |g 2 0 and |e 1 ψ ± 2 shifts the energy level of | f 1 ψ ± 1 by an amount of −hδ 2 , where δ 2 = 9λ 2 1 /4α 1 (see Supplemental Material), and |h 1 is the fourth level of Q 1 and α 1 is its anharmonicity (α j = 2ω e,j − ω f ,j , j = 1, 2). To compensate for these shifts, the angular frequency of the field driving |g 2 ←→ |e 2 should be set to ω d, With this setting and performing the transformation exp (iH int t/h), the system dynamics associated with Q 1 's state | f 1 can be described by the effective Hamiltonian where tan φ 2 = Ω e f /Ω ge .
When Q 1 is initially in the state |g 1 , it is also strongly coupled to the resonator, and there are three dressed states in the single-excitation subspace: where tan θ = λ 2 /λ 1 . The corresponding eigenenergies arē hω e andh ω e ± λ 2 1 + λ 2 2 , as shown in Fig. 1c. When and Ω e f , the qubits cannot make any transition between each of these singleexcitation dressed states and the state |g 1 g 2 0 or |g 1 f 2 0 as each of these transitions is highly detuned from the driving fields. As a consequence, Q 2 is not affected by the driving fields when Q 1 is initially in the state |g 1 . Therefore, the system dynamics is described by the effective Hamiltonian of Eq. (5). The evolution of the initial basis states |c 1 d 2 0 (c, d = g, f ) are given by When Ω e f /Ω ge remains unchanged during the interaction, the evolution satisfies the parallel-transport condition ψ cd (t)| H eff ψ c d (t) = 0, and hence is purely geometric. If the Rabi frequencies of the driving fields and the interaction time are appropriately chosen so that T 0 Ωdt = π, the degenerate qubit subspace undergoes a cyclic evolution. Consequently, the qubits return to the computational space {|g 1 g 2 , |g 1 f 2 , | f 1 g 2 , | f 1 f 2 } with the resonator left in the vacuum state |0 after the time T. With this setting, the evolution operator of the qubits in the computational basis is which is a non-Abelian holonomy. For φ = π/2, i.e., Ω ge = Ω e f , this corresponds to a CNOT gate, which flips the state of the target qubit conditional on the control qubit being in the state | f 1 . 4. Measured process matrix for the realized CNOT gate. The process matrix is measured by preparing a set of 36 distinct input product states in the computational basis {|g 1 g 2 , |g 1 f 2 , | f 1 g 2 , | f 1 f 2 } and reconstructing the density matrices for these states and for the output states produced by the CNOT gate. The |e -state populations of Q 1 and Q 2 , averaged over the 36 output states, are 2.2% and 2.8%, respectively.

EXPERIMENTAL IMPLEMENTATION
The experiment is performed in a superconducting circuit involving five frequency-tunable qubits, labeled from Q 1 to Q 5 , coupled to a resonator with a fixed frequency ω r /2π =5.584 GHz [21,38,39]. In our experiment, Q 1 and Q 2 , whose anharmonicities are 2π × 242 MHz and 2π × 249 MHz, are used as the control and target qubits, respectively. The on-resonance coupling strengths of the g-e transitions of Q 1 and Q 2 to the resonator are respectively λ 1 =2π × 20.8 MHz and λ 2 =2π × 19.9 MHz. The energy relaxation time T 1 and pure Gaussian dephasing time T * 2 for the basis state | f of Q 1 (Q 2 ) are 13.0 (10.7) µs and 2.1 (1.5) µs, while those for the intermediate state |e are 23.9 (15.9) µs and 2.7 (2.1) µs, respectively. The other qubits are on far off-resonance with the resonator so that their interactions with the resonator are effectively switched off throughout the gate operation. We note that during the gate operation, the two qubits have a probability of being populated in | f 1 e 2 , which is significantly coupled to |e 1 f 2 via virtual photon exchange as the two qubits almost have the same anharmonicity α 2π × 240 MHz. To suppress this coupling, Q 1 should be detuned from Q 2 by an amount much larger than λ 1 λ 2 /α. This detuning slightly changes the energy level configuration of the dressed states associated with Q 1 's initial state |g 1 , but does not affect the gate dynamics.
As shown in Fig. 2, the experiment starts with the initialization of Q 1 and Q 2 to the ground state |g at their idle frequencies 5.47 GHz and 5.43 GHz, respectively, which is followed by the preparation of each qubit in one of the six Except |g , each of the other single-qubit states is produced by a g-e π/2-or π-pulse followed by a e-f π-pulse. After these effective single-qubit rotations, these two qubits are prepared in a product state. We then apply square Z pulses to both qubits, tuning their |g ←→ |e transition frequencies to 5.58 GHz and 5.584 GHz and thus switching on their interactions with the resonator. Accompanying these Z pulses, a driving pulse composed of two components with frequencies of 5.565 GHz and 5.369 GHz is applied to Q 2 , resonantly connecting the computational states |g 2 0 and | f 2 0 to the dressed state ψ − 1 . The Rabi frequencies of these driving fields are Ω ge = Ω e f = 2π × 2.2 MHz. Since the resonator is initially in the vacuum state, the system dynamics is governed by the effective Hamiltonian (5) and the time evolution given by Eq. (10). After a duration of τ = 205 ns, the CNOT gate is realized.
One of the most important features of the CNOT gate is that it can convert a two-qubit product state into an entangled state. In particular, when the control qubit is initially in the superposition state (|g 1 + | f 1 ) / √ 2 and the target state in |g 2 , they will evolve to the maximally entangled state (|g 1 g 2 + | f 1 f 2 ) / √ 2 after this gate. We measure this output state by quantum state tomography. This is realized by subsequently biasing each of the two qubits back to its idle frequency right after the gate operation, applying an e-f π-pulse to each qubit, and measuring its state along one of the three orthogonal (X, Y, and Z) axes of the corresponding Bloch sphere with respect to the basis {|g , |e }. The Z measurement is directly realized by state readout, while the X (Y) measurement realized by the combination of a g-e π/2-pulse and state readout. The reconstructed output twoqubit density matrix is displayed in Fig. 3, which has a fidelity of 0.935 ± 0.016 to the ideal maximally entangled state, and a concurrence of 0.888 ± 0.029. To fully characterize the performance of the implemented CNOT gate, we prepare a full set of 36 distinct two-qubit input states before the two-qubit gates, and measure these states and the corresponding output states. With these measured results, the process matrix for the gate operation is reconstructed. The measured process matrix, χ meas , is presented in Fig. 4. The gate fidelity, defined as F = tr (χ id χ meas ), is 0.905 ± 0.008, where χ id is the ideal process matrix. The measured fidelity is in well agreement with the numerical simulation based on the Lindblad master equation, which yields a fidelity of 0.908. One of the error sources is the transitions from |g 1 g 2 0 and |g 1 f 2 0 to Φ 0 1 and Φ ± 1 and the transition from ψ − 1 to ψ − 2 induced by the drive, which cause quantum information leakage to the noncomputational space. Such a leakage error can be mitigated through the improvement of the qubit's nonlinearity or by balancing the drive amplitude and the gate operation time provided the qubits' coherence is bettered, which allows the gate fidelity to be increased by about 6.5% (see Supplemental Material). On the other hand, the qubits' energy relaxation and their dephasings contribute about 1.8% and 1.6% of the error, respectively. Our further numerical simulations show that the CNOT gate with a fidelity above 99% can be obtained with sufficiently large qubit's nonlinearity α j and qubit-resonator coupling strength λ j . For instance, with the parameters λ j /2π = 110 MHz, α j /2π = −3.69 GHz [40,41], Ω ge /2π = Ω e f /2π = 5.9 MHz, T 1 = 60 µs, and T * 2,j = 86 µs, we find a CNOT gate with the operation time about 87 ns and the fidelity of 0.991, which is at the surface code threshold for fault tolerance [42][43][44]. We note this gate is robust against the frequency fluctuations of the driving fields. Suppose that the angular frequencies of these drives deviate from the desired values by an amount of δω = 2π × 100 kHz. The infidelity incurred by this deviation is about [π(δω) 2 /(8Ω 2 ge/e f )] 2 0.1%.

CONCLUSION
In conclusion, we have proposed and demonstrated a scheme for implementing a non-Abelian geometric gate between two superconducting qubits, whose ground and second excited states act as the computational basis states. The conditional dynamics is realized by resonantly driving the transitions between the basis states of the target qubit to the single-excitation dressed states formed by this qubit and the resonator. This entangling gate, together with the previously demonstrated non-Abelian geometric single-qubit gates [30][31][32][33], constitutes a universal set of holonomic gates for realizing quantum computation with superconducting qubits. The method can be directly applied to other systems composed of qubits coupled to a bosonic mode, including cavity QED and ion traps.

DEVICE PARAMETERS
The device used in our experiment is similar to the ones reported in Refs. [1][2][3], it possesses five superconducting Xmon qutrits, whose basis states are denoted as |g , |e and | f , as shown in Fig. 1 of the main text. These qutrits are capacitively connected to a bus resonator, whose rare frequency is measured as 5.584 GHz when all the qutrits are staying in the ground state |g at their respective idle frequencies ω j /2π. Each qutrit's frequency is flexibly adjusted and thus can be controlled to couple to the bus resonator in a dispersive or resonant manner, that induces virtual-photon-mediated [2,4] or real-photon-mediated qubit-qubit coupling [5]. The coupling strength λ j between each qutrit and the resonator through the |g ←→ |e transition is measured through qutrit-resonator vacuum Rabi swap, while keeping the qutrit's |e ←→ | f transition decoupled as the qutrit's anharmonicity α j is much larger than λ j . The device is kept inside a dilution refrigerator with a base temperature below 20 mK. We pick up two qutrits, which are labelled as Q 1 and Q 2 , for our implementation. The related parameters including the qutrit states' coherence times and readout fidelities are characterized and listed in Table S1. As the computation information is encoded in |g and | f (|e as the auxiliary state), we thus refer to the qutrits as qubits. (k = e, f ) are respectively the energy relaxation time, the Ramsey dephasing time and the spin-echo dephasing time of Q j 's state |k measured at the idle point. In addition, α j is the qubit's anharmonicity and λ j is the coupling strength between Q j and the bus resonator. The probability of detecting Q j in state |k when it is prepared in state |k is F |k j . The I-Q data to differentiate these basis states are plotted in Fig. S1.
Note that the two fields also cannot drive the transition from |ψ + 1 ←→ |ψ ± 2 , for which the energy gaps areh[ω r − (1 ∓ √ 2)λ 2 ], largely detuned from the two fields by and respectively. Though the different off-resonant couplings lead to respective energy shifts to |ψ + 1 , such energy shifts are symmetric and thus neutralize to keep |ψ + 1 almost constant. Note also that off-resonant couplings from | f 1 |ψ + 1 to |h 1 |g 2 0 and |e 1 |ψ ± 2 through the resonator photon lead to energy shifts, which areh  S4. The experimental sequence, which consists of three steps: initial state preparation, gate operation and quantum state tomography. The initial state is created by applying microwave pulses with a Gaussian envelop at the idle points. Then, in the second step, qubit frequencies are tuned by rectangular waves to be near (Q 1 ) or on resonance (Q 2 ) with the resonator. Q 2 is subjected to a two-tone microwave pulse with a flattop envelop during the interaction time which lasts about 205 ns. In the third step, tomographic operations are executed before the two-qubit joint readout.  Figure S4 shows the experimental sequence, which is divide into three steps. Firstly, two successive microwave pulses are imposed on each qubit at their idle points to prepare the initial state. The first pulse with the frequency of ω j /2π realizes the |g ↔ |e rotation while the second pulse with the frequency of (ω j − α j )/2π is a flip operation between |e and | f state, known as a e-f π rotation. After the initial state preparation, rectangular pulses are applied to open the qubitresonator interaction for a time of about 205 ns. The control qubit Q 1 is biased to an optimized point close to ω r , while the target qubit Q 2 stays on resonance with the resonator when a two-tone microwave pulse with the angular frequencies ω r − λ 2 − δ 1 − δ 2 and ω f ,2 − ω r + λ 2 + δ 1 + δ 2 are applied on Q 2 . Finally, the qubits are brought back to their idle points for quantum state tomography. To extract the density matrix, we use three tomographic operations {I, X/2, Y/2}

EXPERIMENTAL SEQUENCE OF THE HOLONOMIC GATE
which are executed at the {|g , |e } space after an e-f π rotation for each qubit, as can be seen in the third step of the sequence. For each tomographic operation, we perform the two-qubit joint readout by applying a two-tone measurement pulse to the transmission line, yielding the probabilities of the two-qubit basic states {|g , |e , | f } ⊗2 . As only probabilities of |g and | f state are used for post analysis, the extracted density matrices have a trace value of smaller than 1, which indicates a leakage to the |e state. In Fig. S5 we plot the measured leakage probability in |e state for each qubit.

QUANTUM PROCESS TOMOGRAPHY
Quantum process tomography is executed by performing the state tomography for totally 36 input and corresponding output states after the gate sequence are applied. The input states are two-qubit product states which are produced by applying rotation pulses to each qubit. The mean fidelity characterized by quantum state tomography for all input states and output states are about 0.983 ± 0.003 and 0.915 ± 0.008, respectively, as shown in Fig. S6. The χ-matrix can be extracted from these input and output states by utilizing the least square optimization method with the Hermitian and positive semidefinite constraints [9]. Note that we did not apply the constraint of unit trace for both the calculation of density matrix and χ-matrix considering the leakage to non-computational states.

GATE ERROR ANALYSIS
We have performed numerical simulation to quantify the errors of our gate. The infidelity of our gate mainly comes from the imperfect decoupling between the microwave drive and the qubit and also the decoupling between the qubit and the resonator. For example, when the control qubit Q 1 is in |g state, the detuning between the microwave drive and the dressed state energy level is not large enough to decouple them, which induces a small transition from |g 1 , g 2 , 0 and |g 1 , f 2 , 0 to |Φ 0 1 and |Φ ± 1 , leading to a leakage error. Lowering the drive amplitude can effectively reduce this leakage error, but will extend the evolution time and as a result increase the decoherence error. In addition, the nonlinearities α j need to be larger to better decouple the control qubit Q 1 from the resonator and the target qubit when Q 1 is prepared in | f state during the gate operation.
Considering all these factors, in Table S2, we numerically calculate the χ-fidelities for different nonlinearities and driving amplitudes. When decoherence is neglected, increasing the nonlinearity or decreasing the driving amplitude can both improve the fidelity, as shown in the last column of Table S2. The CNOT gate with χ-fidelity larger than 0.99 can be realized by use of qubits with The gate fidelities for different nonlinearities and driving amplitudes are obtained by optimizing the evolution time and qubit frequencies. The cases with decoherence adopt the T 1 values listed in Table S1 and the pure dephasing times of about 40 µs. The dephasing time used here is estimated from the exponential fit of the Ramsey measurement data before 200 ns. The first column presents numerical data considering parameters of our experimental device, which shows a good agreement with the experimental results. The limitation of nonlinearity, restriction of driving amplitude and decoherence contribute gate errors of about 4.2%, 2.3% and 3.4% respectively.
good coherence, provided the nonlinearity reaches 2π × 1.0 GHz and the driving amplitude reduces to 2π × 1.5 MHz. However, small driving amplitude requires long evolution time, which leads to more decoherence error. For a gate time of about 200 ns, the decoherence contributes about 3.4% of the total error, as shown in the table. To achieve a short evolution time, both the nonlinearity and coupling strength need to be enlarged. Our further numerical simulations show that, given parameters (α j /2π = −3.69 GHz, λ j /2π = 110 MHz, Ω ge,e f /2π = 5.9 MHz, T 1,j =60 µs, T * 2,j =86 µs) which are accessible in recent superconducting qubits [10][11][12], the CNOT gate with an optimized operation time of 87 ns yields an χ-fidelity of 0.991, indicating the potential of our scheme in high-fidelity quantum operations. Fig. S7. Numerical process fidelities of the CNOT gate by varying the anharmonicity α j and coupling strength λ j . Here we have set α ≡ −α j , λ ≡ λ j (j = 1, 2). The simulation results are obtained by optimizing the evolution time and qubit frequencies without considering decoherence. The red dot shows the position of our current device.

ADDITIONAL NUMERICAL SIMULATION
We have further performed numerical simulation by sweeping both the anharmonicity α j and qubit-resonator coupling strength λ j . For simplicity, the drive amplitude is fixed to be Ω ge,e f /2π =2.0 MHz, similar to that used in our experiment. Fig. S7 plots the numerical process fidelities χ of the CNOT gate in parameter space, where the yellow region meets the threshold for surface code. 6