Optimal control of two qubits via a single cavity drive in circuit quantum electrodynamics

TL;DR

This paper demonstrates how optimal control techniques can be used to design high-fidelity, robust two-qubit gates in circuit quantum electrodynamics by driving a single cavity, with all simulations based on realistic experimental parameters.

Contribution

It introduces a method to implement high-fidelity two-qubit gates in circuit QED using a single cavity drive optimized through control techniques, addressing spurious interactions.

Findings

Achieved high-fidelity two-qubit gates with optimized pulse shapes.

Demonstrated robustness against system parameter errors.

Validated results with experimentally relevant parameters.

Abstract

Optimization of the fidelity of control operations is of critical importance in the pursuit of fault-tolerant quantum computation. We apply optimal control techniques to demonstrate that a single drive via the cavity in circuit quantum electrodynamics can implement a high-fidelity two-qubit all-microwave gate that directly entangles the qubits via the mutual qubit-cavity couplings. This is performed by driving at one of the qubits' frequencies which generates a conditional two-qubit gate, but will also generate other spurious interactions. These optimal control techniques are used to find pulse shapes that can perform this two-qubit gate with high fidelity, robust against errors in the system parameters. The simulations were all performed using experimentally relevant parameters and constraints.

Figures

Equations

$$\begin{array}{rl}H=& {\omega }_r{a}^{†}a+\sum _{i=1,2}\sum _{j_i}{\omega }_{j_i}|j_i⟩⟨j_i|+\sum _{i=1,2}{g}_i(a^{†}{c}_i+a{c}_i^{†})\\ & +[\epsilon (t)a^{†}{e}^{-i{\omega }_{d}t}+{\epsilon }^{\ast }(t)a{e}^{+i{\omega }_{d}t}],\end{array}$$

$$\begin{array}{rl}H=& {\omega }_r{a}^{†}a+\sum _{j=1,2}\frac{{\omega }_a^{(j)}}{2}{\sigma }_z^{(j)}+\sum _{j=1,2}{g}_j{\sigma }_x^{(j)}(a^{†}+a)\\ & +\epsilon (t)(a^{†}{e}^{-i{\omega }_{d}t}+a{e}^{+i{\omega }_{d}t}),\end{array}$$

$${\omega }_r\alpha +\epsilon e^{-i{\omega }_{d}t}-i\dot{\alpha }=0,$$

$$\begin{array}{rl}{H}^{\prime }=& {\omega }_r{a}^{†}a+\sum _{j=1,2}\frac{{\omega }_a^{(j)}}{2}{\sigma }_z^{(j)}+\sum _{j=1,2}{g}_j{\sigma }_x^{(j)}(a^{†}+a)\\ & +\sum _{j=1,2}{\Omega }_R^{(j)}{\sigma }_x^{(j)}\cos ({\omega }_{d}t),\end{array}$$

$$U=\mathrm{exp}\bigg{[}\frac{g_{1}}{\Delta_{1}}(a^{\dagger}\sigma_{-_{1}}-a\sigma_{+_{1}})+\frac{g_{2}}{\Delta_{2}}(a^{\dagger}\sigma_{-_{2}}-a\sigma_{+_{2}})\bigg{]}.$$

$$\begin{array}{rl}{H}^{\prime \prime }=& {\omega }_r{a}^{†}a+\sum _{j=1,2}\frac{{\omega }^{\prime }{}_a^{(j)}}{2}{\sigma }_z^{(j)}\\ & +\frac{g_1{g}_2({\Delta }_1+{\Delta }_2)}{2{\Delta }_1{\Delta }_2}({\sigma }_{+}^{(1)}{\sigma }_{-}^{(2)}+{\sigma }_{-}^{(1)}{\sigma }_{+}^{(2)})\\ & +\sum _{j=1,2}{\Omega }_R^{(j)}{\sigma }_x^{(j)}\cos ({\omega }_{d}t),\end{array}$$

$$\omega^{\prime}{}_{a}^{(j)}=\omega_{a}^{(j)}+2\frac{g_{j}^{2}}{\Delta_{j}}\bigg{(}a^{\dagger}a+\frac{1}{2}\bigg{)}.$$

$$H^{eff}=H^0+\frac{1}{2}[H^2,S^{(1)}]+H_d(t)+[H_d(t),S^{(1)}],$$

$$S^{(1)}=-\frac{J}{{\Delta }_{12}}({\sigma }_{{+}_1}{\sigma }_{{-}_2}-{\sigma }_{{-}_1}{\sigma }_{{+}_2}),$$

$$\begin{split}H^{eff}=&\omega_{r}a^{\dagger}a+\frac{\tilde{\omega}_{a}^{(1)}}{2}\sigma_{z}^{(1)}+\frac{\tilde{\omega}_{a}^{(2)}}{2}\sigma_{z}^{(2)}\\ &+\Omega_{R}^{(1)}\bigg{(}\sigma_{x}^{(1)}+\frac{J}{\Delta_{12}}\sigma_{z}^{(1)}\sigma_{x}^{(2)}\bigg{)}\mathrm{cos}(\omega_{d}t)\\ &+\Omega_{R}^{(2)}\bigg{(}\sigma_{x}^{(2)}-\frac{J}{\Delta_{12}}\sigma_{x}^{(1)}\sigma_{z}^{(2)}\bigg{)}\mathrm{cos}(\omega_{d}t),\end{split}$$

$$i\dot{U}(t)=\bigg{(}H_{0}+\sum_{j}c_{j}(t)H_{j}\bigg{)}U(t),$$

$$\mathcal{F}=\bigg{\lvert}\frac{1}{n_{s}}\mathrm{tr}[W^{\dagger}\hat{O}U(T)\hat{O}]\bigg{\rvert}^{2},$$

$$H(t_k)=H_0+\sum _j{c}_j(t_k)H_j,t_k=k\tau .$$

$$\begin{array}{rl}& \mathrm{maximize}{\mathcal{F}}_{wc}={\min }_{\delta }\mathcal{F}(\theta ,\delta )\\ & \mathrm{subject}\mathrm{to}\theta \in \Theta ,\delta \in \Delta ,\end{array}$$

$$\begin{array}{rl}& \mathrm{maximise}{\min }_i[\mathcal{F}(\theta ,{\delta }_i)+{\tilde{\theta }}^{⊺}{\nabla }_{\theta }\mathcal{F}(\theta ,{\delta }_i)]\\ & \mathrm{subject}\mathrm{to}\theta +\tilde{\theta }\in \Theta ,\tilde{\theta }\in {\tilde{\Theta }}_{trust}.\end{array}$$

$$H_0={\Delta }_r{a}^{†}a+\sum _{j=1,2}\frac{{\tilde{\Delta }}_a^{(j)}}{2}{\sigma }_{z_j},$$

$$\begin{split}H_{j}=&\frac{2g_{1}}{\omega_{r}-\omega_{d}}\bigg{(}\sigma_{x_{1}}+\frac{J}{\Delta_{12}}\sigma_{z_{1}}\sigma_{x_{2}}\bigg{)}\\ &+\frac{2g_{2}}{\omega_{r}-\omega_{d}}\bigg{(}\sigma_{x_{2}}-\frac{J}{\Delta_{12}}\sigma_{x_{1}}\sigma_{z_{2}}\bigg{)},\end{split}$$

$$U_{des}=\mathrm{exp}\bigg{(}-i\frac{\pi}{4}\sigma_{x_{1}}\sigma_{z_{2}}\bigg{)},$$

$$\tilde{U}^{m}\ket{+y}\ket{+y}=\Pi^{m}_{k=1}exp\big{(}-iH(t_{k})\tau\big{)}\ket{+y}\ket{+y},$$

$$\begin{array}{rl}& {\tilde{U}}^8|+y⟩|+y⟩\approx \alpha |{\Phi }^{+}⟩+\beta |{\Psi }^{-}⟩,\\ & F\approx |⟨{\Phi }^{-}|\left(\alpha |{\Phi }^{+}⟩+\beta |{\Psi }^{-}⟩\right){|}^2=0,\end{array}$$

$$s(t_l)=\sum _{k=0}^{N-1}{T}_{l,k}c(t_k),$$

$$\begin{split}H=&\omega_{r}a^{\dagger}a+\sum_{j=1,2}\bigg{(}\omega_{a}^{(j)}c_{j}^{\dagger}c_{j}+\frac{\delta_{j}}{2}c_{j}^{\dagger}c_{j}(c_{j}^{\dagger}c_{j}-1)\bigg{)}\\ &+\sum_{j=1,2}g_{j}(a^{\dagger}c_{j}+ac_{j}^{\dagger})\\ &+(\varepsilon(t)a^{\dagger}e^{-i\omega_{d}t}+\varepsilon^{*}(t)ae^{+i\omega_{d}t}).\end{split}$$

$$\begin{array}{rl}{H}_d=& {\epsilon }_x(t)(a^{†}{e}^{-i{\omega }_{d}t}+a{e}^{+i{\omega }_{d}t})\\ & +i{\epsilon }_y(t)(a^{†}{e}^{-i{\omega }_{d}t}-a{e}^{+i{\omega }_{d}t}).\end{array}$$

Full text

††preprint: APS/123-QED

I Introduction

Fault-tolerant architectures for quantum computation require individual gate operations to be performed with high fidelity Aharonov and Ben-Or (2008). This requires strict control of all the parameters in the system and robustness of the performed function with respect to any noise or uncertainties in the system. In superconducting circuits, experimental two-qubit gate fidelities $\mathcal{F}>0.99$ have recently been demonstrated Sheldon et al. (2016); Martinis and Geller (2014); Barends et al. (2014); Ghosh et al. (2013), at the threshold for fault tolerance with a surface code architecture Fowler et al. (2009, 2012), but well below that required for a gate-model machine Aharonov and Ben-Or (2008); Nielsen and Chuang (2010). Although some of the infidelity is due to decoherence, there are also errors due to fluctuations and inaccuracies in the driving field and measured system parameters.

The systems are also required to be scalable, so that more qubits can be added in order to realize a full quantum computer DiVincenzo (2000). A crucial challenge in scaling up is to maintain high coherence, which becomes more difficult as more controls are added. Therefore, it can be ideal to keep the circuit complexity to a minimum in order to reduce the potential avenues for decoherence. Then, the problem becomes one of how to control the system with fewer sources. Control in quantum mechanics is an open problem and equivalent to a non-linear optimization, which is difficult to parametrize and solve analytically, but is an important problem for quantum computing.

In this paper, we investigate an all-microwave gate in circuit quantum electrodynamics (cQED) using a single microwave drive applied to a single cavity containing two transmon Koch et al. (2007) qubits. The system dynamics are controlled by selecting the frequency, phase, and amplitude of the microwave control drive Rigetti and Devoret (2010); de Groot et al. (2010, 2012). Using only a single microwave drive for the control means the circuit complexity is kept to a minimum while also minimizing the external sources of noise. We use the sequential convex programming (SCP) algorithm developed in Kosut et al. (2013) to find pulses that are capable of implementing an entangling operation between two qubits in one cavity to a high degree of fidelity via the cross-resonance gate Rigetti and Devoret (2010); Chow et al. (2011). This gate was first implemented in circuit QED using microwave drives local to each qubit and utilizes the qubit-qubit coupling via the common mode of the cavity. Using a fixed system and modifying only the drive pulses, we show that it is possible to implement the desired unitary with fidelities $\mathcal{F}>0.9875$ using a two-level approximation for the qubits, and $\mathcal{F}>0.9639$ for a full multi-level transmon model. We also show that these fidelities are robust to uncertainty in the system parameters and incorporation of experimentally realistic pulse filtering, while also imposing constraints on the drive power.

Advanced pulses for control in quantum systems have been investigated previously, such as dynamical decoupling schemes which use trains of pulses to cancel out environmental noise and get rid of any dephasing on the qubit Viola and Lloyd (1998); Viola et al. (1999). DRAG (derivative removal via adiabatic gate) was designed in order to remove any leakage to the non-computational levels of the qubit, here two pulses are applied to the system with the second being the derivative of the first Motzoi et al. (2009); Chow et al. (2010). Other examples include spline-shaped pulses for the resonator-induced phase gate where the pulse shapes are used to remove unwanted effects during the gate such as photon loss and residual entanglement between the qubits and cavity after the gate implementation Cross and Gambetta (2015), and SWIPHT (sppeding up wave forms by inducing phases to harmful transitions) that purposely drives the nearest harmful transition so that it undergoes one cyclic revolution while maximizing the fidelity to the desired unitary Economou and Barnes (2015). While these methods have proven to be effective, they have not included any potential sources of error in the system which will inevitably be present in a real experiment. More recently, optimal control methods were used to find the minimum time to perform quantum operations in superconducting systems, in particular the cross-resonance gate was investigated for this purpose using the common method of dedicated qubit control drives Kirchhoff et al. (2015). Work was also performed in superconducting circuits for a single global cavity drive, searching for the best parameter regime to perform single- and two-qubit gates in the shortest possible time and maximizing the fidelity to the gates of interest for the best parameter regime Goerz et al. (2016). In these cases again potential error sources were not included in the investigation.

Quantum optimal control theory is one method for designing pulse shapes to perform required interactions Werschnik and Gross (2007); Brif et al. (2010). In particular, it has been used to develop algorithms for numerically designing pulse shapes that maximize a chosen fidelity function Khaneja et al. (2005); De Fouquieres et al. (2011); Ferrie and Moussa (2015); Machnes et al. (2015); Schirmer and De Fouquieres (2011); Kosut et al. (2013); Goerz et al. (2014); Spörl et al. (2007). It has shown success in designing pulses for single-qubit gates, two-qubit gates, and readout Ginossar et al. (2010). There have also been applications of quantum optimal control algorithms towards high fidelity gates that are robust to errors in the system Kosut et al. (2013); Dong et al. (2014); Goerz et al. (2014), but these have focused either on single-qubit gates, assuming that the control came via changing the system parameters, for example, via flux tuning, or by direct driving of the qubit of interest.

The rest of the paper is organized as follows: In Sec. II we define the Hamiltonian for two qubits in a three-dimensional (3D) cavity driven by a single microwave drive, then perform a set of transformations relevant to the system in order to derive the two-qubit interaction. We see that upon making these transformations that there will be some unwanted interaction terms that are being driven resonantly that will be harmful to the desired entangling operation. In Sec. III we describe the quantum optimal control method, and in particular we review the sequential convex programming algorithm and how it can be used to perform robust quantum optimal control. In Sec. IV we present our results and discussion, showing that using SCP we can design a pulse shape that can perform the desired entangling operation with high fidelity for both a two-level system and a multi-level system even in the presence of resonant unwanted rotations, errors in the system parameters, and pulse filtering.

II All-Microwave Gate Using a Single Cavity Drive

The cross-resonance gate is an all-microwave gate that utilizes coupling between two qubits to generate a two-qubit operation Chow et al. (2011); Rigetti and Devoret (2010). For all-microwave gates the system is set up so that the qubits are far detuned from one another, and also from the cavity (i.e., the system is in the dispersive regime). This ensures that the qubits and cavity do not interact until some external control is applied, and will extend the lifetimes of the qubits when they are not being operated on. In the case where there are microwave drives local to each qubit, a microwave resonant with the target qubit would be applied directly to the other control qubit via a direct microwave line. Due to the coupling between the two qubits via the common mode of the cavity this would activate a two-qubit operation which could be used to generate entanglement. This system ideally requires a local drive on each qubit, but it is interesting to consider a case in which only a single drive is available that globally addresses both qubits. This could be an advantage in large-scale systems, as it reduces the number of required controls.

One such simple architecture is to have a single drive coupled to the cavity and use this to perform all the control; this is relevant when the system consists of two-qubits in a 3D cavity. This has been shown to be of interest due to the long coherence times enabled by the design Paik et al. (2011). In the case of the cross-resonance gate this would mean that the microwave drive that is resonant with the target qubit will not only activate the two-qubit operation, but will also cause unwanted rotations of the target qubit.

In superconducting circuits, one of the qubits of interest, and the focus of this work, is the transmon Koch et al. (2007). This is a charge-based qubit that is designed in such a way as to be less susceptible to charge noise, thus increasing coherence times. However, the trade off is that the system is less anharmonic and is thus not a true two-level system but is in fact multi-level. The Hamiltonian for two transmons coupled to a resonator with a single-microwave control is given by:

[TABLE]

where $\omega_{r}$ is the cavity resonance frequency, $i$ denotes transmon $1$ and $2$, $\omega_{j_{i}}$ is the frequency of the $j^{\mathrm{th}}$ level of the $i^{\mathrm{th}}$ transmon $\ket{j_{i}}$, $g_{i}$ is the coupling between transmon $i$ and the cavity, $\varepsilon(t)$ is the time-dependent pulse envelope, $\omega_{d}$ is the drive frequency, $a^{(\dagger)}$ are the annihilation (creation) operators of the cavity mode photon, and $c^{(\dagger)}$ are the annihilation (creation) operators of the transmon excitations.

In order to investigate the two-qubit operation of interest, it is instructive to first look at a case in which the transmon is approximated as a two-level system. This will be less computationally expensive in numerical optimization calculations, and will give an indication as to the effectiveness of these techniques in such a system while also revealing the single- and two-qubit interactions during the operation. The Hamiltonian for two two-level systems and a microwave drive coupled to a common cavity mode is given by:

[TABLE]

where $j$ denotes qubits 1 and 2, $\omega_{a}^{(j)}$ is qubit $j$ transition frequency, $g_{j}$ is the coupling between qubit $j$ and the cavity, and $\sigma_{x/z}^{(j)}$ are the Pauli spin matrices for qubit $j$. Here we have assumed a real drive.

When the drive is far detuned from the cavity, as is the case here, the effect of the drive on the cavity is not important. In order to see the effect of the drive on the qubits, the time-dependent displacement operator $D(\alpha)=\mathrm{exp}(\alpha a^{\dagger}-\alpha^{*}a)$ is be applied Blais et al. (2007), choosing $\alpha$ to satisfy:

[TABLE]

this Hamiltonian becomes:

[TABLE]

where $\Omega_{R}^{(j)}=2g_{j}\varepsilon/(\omega_{r}-\omega_{d})$. This shows an important difference to Chow et al. (2011), where each drive term would only apply to the qubit it was localized to. In this case, the single-microwave control drives both qubits.

Since the qubits are far detuned from the cavity, with $g_{j}\ll|\Delta_{j}|=|\omega_{a}^{(j)}-\omega_{r}|$, the dispersive transformation can be applied:

[TABLE]

Expanding to second order in the small parameter $g_{j}/\Delta_{j}$ and dropping fast oscillating terms gives the effective Hamiltonian:

[TABLE]

where $\omega^{\prime}{}_{a}^{(j)}$ is the shifted qubit frequency given by:

[TABLE]

This is now the Hamiltonian for two coupled qubits with a single-microwave drive that drives both qubits.

Following Ref. Richer (2013), the final step is to perform a Schrieffer-Wolff transformation which transforms the effective Hamiltonian via:

[TABLE]

where $H^{0}$ is the unperturbed part of the Hamiltonian [given on the first line of Eq. (6)], $H^{2}$ is the small perturbation term that contains off-diagonal terms [second line of Eq. (6)], and $H_{d}(t)$ is the drive term [third line of Eq. (6)]. For this derivation, $S^{(1)}$ is given by:

[TABLE]

with $J=g_{1}g_{2}(\Delta_{1}+\Delta_{2})/2\Delta_{1}\Delta_{2}$ and $\Delta_{12}=\omega_{a_{1}}^{\prime}-\omega_{a_{2}}^{\prime}$, where $J\ll\Delta_{12}$. The final Hamiltonian is then given by:

[TABLE]

where $\tilde{\omega}_{a}^{(1)}/2=\omega^{\prime}{}_{a}^{(1)}+J^{2}/\Delta_{12}$, and $\tilde{\omega}_{a}^{(2)}=\omega^{\prime}{}_{a}^{(2)}-J^{2}/\Delta_{12}$. This Hamiltonian contains the two-qubit terms $\sigma_{z}^{(1)}\sigma_{x}^{(2)}$ and $\sigma_{x}^{(1)}\sigma_{z}^{(2)}$. To activate one of these terms a drive at the correct frequency must be applied. In order to perform the $\sigma_{x}^{(1)}\sigma_{z}^{(2)}$ operation, for example, a microwave drive at $\omega_{d}=\tilde{\omega}_{a}^{(1)}$ must be applied. This $\sigma_{x}^{(1)}\sigma_{z}^{(2)}$ is two single-qubit rotations away from being a controlled-not (CNOT) gate.

According to Eq. (10), choosing $\omega_{d}=\tilde{\omega}_{a}^{(1)}$ will not only perform the $\sigma_{x}^{(1)}\sigma_{z}^{(2)}$ operation, this will also drive a single-qubit rotation on qubit 1 while also driving some off-resonant rotations of qubit 2 and the $\sigma_{z}^{(1)}\sigma_{x}^{(2)}$ term. However, performing the relevant single- or two-qubit operation is not as simple as just choosing the correct frequency and driving, the pulse shape also plays an important part. For example, choosing to drive at $\omega_{d}=\tilde{\omega}_{a}^{(1)}$ using the pulse shape shown in Fig. 1 would achieve a fidelity of $\mathcal{F}=0.1461$ for a $\pi/2$ rotation on the $\sigma_{x_{1}}\sigma_{z_{2}}$ operation, and $\mathcal{F}=0.3230$ for a $\pi/2$ rotation on the $\sigma_{x_{1}}$ operation.

In very few cases, solutions for the problem of finding pulse shapes to perform required operations can be found analytically. In most cases, however, such as here, this is not possible, as the solutions must be found numerically. There may also exist many solutions to the problem; this makes finding a solution non-trivial.

III Quantum Optimal Control

The aim of optimal control theory is to maximize or minimize a function subject to certain constraints and bounds. For quantum information processing the function we wish to maximize is the fidelity of the system unitary evolution after some time $T$, $U(T)$, with respect to a desired target unitary operator, $W\in\mathbb{C}^{n_{s}\times n_{s}}$ where $n_{s}$ is the size of the system. In quantum systems the evolution of the system can be described by:

[TABLE]

with $H_{0}$ the drift (i.e. not controlled) Hamiltonian, $j$ the index of the control, $c_{j}(t)$ the control function at time $t\in[0,T]$, and $H_{j}$ the control Hamiltonian that is under the influence of control $j$. The fidelity function used here is:

[TABLE]

where $\hat{O}$ is the projector onto the subspace of interest. With this definition we have $\mathcal{F}\in[0,1]$, when $\mathcal{F}=1$ there is no measurable distinction between $W$ and $U(T)$, and no leakage out of the relevant subspace.

In general, maximizing the fidelity function is performed numerically and can be very computationally expensive, particularly when using smooth analytic control functions, due to the large size of the control space (i.e. the dimension of the control space). For the purposes of optimization it can be more efficient to approximate the controls as a series of piecewise constant amplitudes over $N$ uniform time intervals $\tau=T/N$, giving the Hamiltonian for time $t_{k}=k\tau$ in Eq. ((11)) as:

[TABLE]

In this case, and with $H_{j}$ time-independent, $U(T)$ is now given by $U(T)=\prod_{k=1}^{N}U_{k}$ with U_{k}=\mathrm{exp}\big{(}-iH(t_{k})\tau\big{)}.

The most successful methods for finding optimal controls have been local optimizers based on gradient ascent Khaneja et al. (2005); Brif et al. (2010); Werschnik and Gross (2007). For these algorithms an initial control pulse is chosen, the fidelity and the gradient are calculated and the optimizer is run which uses the gradient at each iteration of the optimization to determine in which direction in the landscape to travel. These algorithms have been shown to converge on a solution in fewer iterations than other algorithms and are computationally less expensive, however they do not guarantee a global maximum as this depends on where you start in the control landscape.

When applying optimal control to quantum information problems the pulses that are developed must be robust to any errors in the system, such as errors in the measured system parameters. In order to perform numerical optimal control simulations the system parameters must be entered into the simulator, but this requires the parameters to be known with absolute precision which will not be the case in experiment. The estimated parameters will come with some error range, which must be included in the numerical calculation. As discussed in Kosut et al. (2013), a method for doing this is to sample points $\delta$ from the parameter range $\Delta$, calculate the fidelity for each of these points and then to maximize the worst-case fidelity, $\mathcal{F}_{wc}$:

[TABLE]

where $\theta$ are the controls and $\Theta$ is the set of allowed controls subject to any constraints imposed. Maximizing the minimum fidelity of the range, rather than the average fidelity of the range, places a more strict requirement on the optimizer and ensures that the control found by the optimizer does not include points with very low fidelity that are balanced by high fidelity terms.

Sequential convex programming

The method proposed for solving the optimization problem stated above is sequential convex programming Schittkowski and Zillober (1995); Kosut et al. (2013). This is a gradient based local optimizer and is therefore efficient at locating local optimal solutions, provided the initial guess is a good one. The algorithm is initialized with a control $\theta\in\Theta$, which is typically a convex set, or well approximated by one. Points $\delta_{i}$ are then sampled from the error range $\Delta$, and a convex trust region $\tilde{\Theta}_{trust}$ is initialized. The trust region is chosen so that the linearized fidelity $\mathcal{F}(\theta,\delta_{i})+\tilde{\theta}^{\intercal}\nabla_{\theta}\mathcal{F}(\theta,\delta_{i})$, where $\tilde{\theta}\in\tilde{\Theta}_{trust}$, used in the optimization step retains sufficient accuracy.

With the initial points set, the fidelities $\mathcal{F}(\theta,\delta_{i})$ and gradients $\nabla_{\theta}\mathcal{F}(\theta,\delta_{i})$ are calculated for each points $\delta_{i}$ selected from the error range. The linearised fidelity is then used to solve for the increment $\tilde{\theta}$ in the convex optimisation:

[TABLE]

If $\mathrm{min}_{i}\mathcal{F}(\theta+\tilde{\theta},\delta_{i})>\mathrm{min}_{i}\mathcal{F}(\theta,\delta_{i})$ then replace $\theta$ by $\theta+\tilde{\theta}$, increase the trust-region $\tilde{\Theta}_{trust}$ and repeat the process of calculating fidelities, gradients and solving for the increment $\tilde{\theta}$. If, however, the $\mathrm{min}_{i}\mathcal{F}(\theta+\tilde{\theta},\delta_{i})<\mathrm{min}_{i}\mathcal{F}(\theta,\delta_{i})$ then decrease the trust region $\tilde{\Theta}_{trust}$ and repeat the optimisation step with the same $\theta$. This process is repeated until some stopping criteria are satisfied, such as the number of iterations reaching the maximum number imposed or the gradient is below some threshold such that it is flat and has thus found a local maximum.

IV Cross Resonance Using SCP

To perform quantum optimal control the effective Hamiltonian in eqn. (10) must be cast in the form of eqn. (13). This investigation focuses on having all the system parameters fixed and the control variable as the pulse shape $\varepsilon(t)$. To begin the investigation, we first look at the Hamiltonian in the two-level approximation. Moving to a frame rotating at the drive frequency, the drift Hamiltonian is given by:

[TABLE]

where $\Delta_{r}=\omega_{r}-\omega_{a}$ and $\tilde{\Delta}_{a_{j}}=\tilde{\omega}_{a}^{(j)}-\omega_{d}$. The control Hamiltonian is given by:

[TABLE]

and thus $c(t)=\varepsilon(t)$. In these simulations the system parameters are set as $\omega_{r}/2\pi=6.44$ GHz, $\omega_{a}^{(1)}/2\pi=4.50$ GHz, $\omega_{a}^{(2)}/2\pi=4.85$ GHz, $g_{1}/2\pi=g_{2}/2\pi=133$ MHz; these are parameters that have been used in a previous cQED experiment Leek et al. (2009). The desired unitary is:

[TABLE]

which, when applied to the state $\ket{+y}\ket{+y}$ [where $\ket{+y}=(\ket{0}+i\ket{1})/\sqrt{2}$], gives the maximally entangled Bell state $\ket{\Phi^{-}}=(\ket{00}-\ket{11})/\sqrt{2}$. In order to perform this, we choose $\omega_{d}=\tilde{\omega}_{a_{1}}$.

To ensure realistic pulses are produced, constraints on the maximum and minimum amplitudes have been imposed. By doing this the control space becomes fragmented, as there will be areas that are off limits to the optimizer. This means that the control space may contain many local maxima that the optimizer may become stuck in, as the optimization problem is very sensitive to initial conditions.

In Chow et al. (2011), it was stated the pulse used to perform the cross-resonance gate was a slow Gaussian turn on with a flat top. Therefore, as an initial pulse guess for the SCP algorithm we have chosen to use a pulse that has the form of a Gaussian turn on, a flat top, and then a Gaussian turn off, as shown in Fig. 1.

From these initial conditions, with $\mathcal{F}=0.1461$, a pulse sequence can be generated that is able to perform the desired unitary with a fidelity of $\mathcal{F}=0.9945$. Figure 2 shows the changes that have been made from the original pulse sequence, shaded area, to the optimal pulse shape, outlined area. Due to a good choice of initial guess the optimiser hugs this shape and rapidly finds a solution with high fidelity, even with the constraint on maximum amplitude and few pulse pixels to optimize over. Figure 3 shows, with respect to the state:

[TABLE]

how the entanglement (given by $2|ad-bc|$ for an arbitrary two-qubit state $\ket{\psi}=a\ket{00}+b\ket{01}+c\ket{10}+d\ket{11}$) and the fidelity (given by $\mathcal{F}=|\bra{\Phi_{-}}\tilde{U}^{m}\ket{+y}\ket{+y}|^{2}$) change with each successive $c(t_{k})$. We see that although the fidelity fluctuates, the entanglement monotonically increases with each pulse. For example, at $k=8$, the entanglement is higher than 0.6 but the fidelity is nearly 0 because the outcome state is given by

[TABLE]

where $\ket{\Phi^{+}}=(\ket{00}+\ket{11})/\sqrt{2}$, $\ket{\Psi^{-}}=(\ket{01}-\ket{10})/\sqrt{2}$, and $\alpha,\beta\in\mathbb{C}$. This suggests that the pulse shape continually performs the desired two-qubit operations, but has optimised to produce single-qubit rotations at each step such that at the end the qubits will be rotated into the correct basis and that the function performed at the end is effectively just the desired two-qubit operation.

For the purposes of our optimization, the assumption of piecewise constant controls has been made. This is to make the optimization less computationally expensive as there will be fewer parameters to optimize; the control space dimension will be given by the number of piecewise constant parts in the pulse. It has been suggested that piecewise constant controls are not as powerful as analytic controls, and that analytic controls should always be used where possible Machnes et al. (2015). However, the assumption was made that the piecewise constant controls were used as an approximation of a continuous function and thus would need many time slices to accurately portray such a function. This then made the optimization very computationally expensive. Here instead we choose our control amplitudes, not as an approximation of a continuous function, but as a stand alone pulse shape. State-of-the-art arbitrary waveform generators (AWGs) can generate pulses that change amplitude at least on the order of $1$ ns. We found that the optimal choice of pulse with fewer piecewise constant parts was a pulse over $200$ ns that has 16 individual constant pulses that can be optimized, so that each constant pulse is over an interval of $12.5$ ns. This is well within the capability of many AWGs, and allows us to carry out a resource-efficient optimization. Use of pulses with shorter time steps, to approximate a continuous function, does not necessarily improve the fidelity. In Fig. 4, a pulse sequence generated is shown where each pulse length is 2 ns, getting closer to the limit of the AWGs, which implements the desired unitary with $\mathcal{F}=0.9940$. This optimisation performed the calculations around 12 times more slowly than with fewer piecewise constant parts (3 hours for the larger control space compared with 15 minutes for the smaller one).

Optimal pulses with filtering effect

In general the pulse shown in figure 2 will not be the pulse that reaches the cavity and qubits. The microwave pulse will be (mostly low pass) filtered by control hardware; we can model this by discretizing the control further and using a linear transfer function to approximate the filtering effect Motzoi et al. (2011). In this case, each new piecewise constant pixels $s(t_{l})$ are given by:

[TABLE]

where $T_{l,k}$ is the transfer function as given in Motzoi et al. (2011) for a Gaussian filter.

This is applied to the pulse in Fig. 2, the fidelity calculated drops to $\mathcal{F}=0.8303$, but upon re-optimizing with this pulse as the initial guess a fidelity of $\mathcal{F}=0.9947$ is achieved. Figure 5 shows the optimal pulse before re-optimizing (colored area) and the new pulse generated after the filter effect is taken into account (outlined area). Thus, pulses that are robust to any filtering effect may be produced by the algorithm, which is ideal for experimental implementation.

Robust pulses for errors in system parameters

As mentioned in Sec. III, the values of the system parameters are going to have some error associated with them, as any measured system parameter is given with some error range. A sample from the error ranges must be taken and the SCP algorithm can be performed while optimizing for the minimum fidelity in the range. This is similar to the approach taken in Dong et al. (2014), except that the more stringent condition of maximizing the minimum fidelity is made, as in Kosut et al. (2013), rather than optimizing for the average fidelity. This method is first performed for an error of $\pm 1\%$ in the transition frequency of qubit 2, $\omega_{a}^{(2)}/2\pi=4.85\pm 0.05$ GHz. 11 points are sampled from the error range, which has been shown to be adequate to cover the range Dong et al. (2014). Figure 6 shows the range of fidelities for each point sampled from the error range of $\omega_{a}^{(2)}$. This shows that even in the presence of error in one parameter the SCP algorithm has been able to find a solution that can produce the desired unitary with $\mathcal{F}>0.986$ for the whole range. When filtering is included a solution is found which gives the range of fidelities in figure 7.

Since there are two qubits, potential errors in both of the qubits parameters, $\omega_{a}^{(1)}$ and $\omega_{a}^{(2)}$, must be accounted for. This proves to be a greater challenge since the drive is on resonance with the dressed qubit-1 transition frequency and thus will cause the drive to become slightly off resonant if there are errors. Reducing the error range for $\omega_{a}^{(1)}$ to $\pm~{}0.1\%$ achieves a solution with the range of fidelities shown in Fig. 8, where all points have fidelity $\mathcal{F}>0.9875$. Using the robustness method can be effective in designing a pulse to be robust to errors in some range, but that it becomes increasingly difficult as more error parameters are introduced. This can perhaps be solved by choosing different start points for the optimiser, but thus far the flat-top-like Gaussian has proven to be the optimal start point for the Hamiltonian and desired unitary of interest. Future work will look at a scheme to over come this.

Multi-level Transmon Model

We now proceed to look at the case in which a full model of the transmons is incorporated. This is more computationally expensive to calculate, but is important for capturing leakage out of the two-level logical basis. In the deep transmon limit, with $E_{J}/E_{C}=100$, the anharmonicities for the transmons are $\delta_{1}/2\pi=-160$ MHz for transmon 1, and $\delta_{2}/2\pi=-170$ MHz for transmon 2. In this limit, the transmons can be approximated as Duffing oscillators, the Hamiltonian for two Duffing oscillators coupled to a common cavity mode with a single drive is given by:

[TABLE]

The first two lines of this equation can be diagonalised to find the Hamiltonian for the cavity + two transmons in the dressed basis, where the frequencies of each component will now include dependencies on all the other parts. For the purposes of optimal control this will now form the drift Hamiltonian. The drive term can then be transformed to form the new operators in the dressed basis, and can be used as the new drive term $H_{c}$ for the simulations.

For the multi-level simulations the drive term has been cast into two parts to include complex control as it has been shown that using both quadratures can be useful in suppressing leakage Motzoi et al. (2009). The drive term in this case becomes:

[TABLE]

In the new transmon limit it is that the control must be discretized more in order to reach a good fidelity. In this case each piecewise control amplitude is now $2$ ns long, well within the capabilities of current AWGs. In the case where there are no errors in the system there are many points in time that perform well with fidelities $\mathcal{F}>0.9999$, due to this extra discretization of the control. Figure 9 shows the fidelity of the optimised pulses, where each initial pulse has taken the form of a flat-top Gaussian with dt=2 ns for each piecewise amplitude, against time. The fidelities converge to $\mathcal{F}>0.999$ for all times $T>100\mathrm{ns}$.

We now investigate the robustness of one of these pulses to parameters variation. In figure 10 we plot the pulse fidelity $\mathcal{F}$ against $\omega_{a}^{(2)}$. While the pulse performs well for the specific chosen parameters, it can clearly be seen that the fidelity is highly dependent on the frequency $\omega_{a}^{(2)}$, hence not robust to variation in this parameter.

As in the previous section, the robust methods are used to find a pulse robust to errors in $\omega_{a}^{(2)}$ first. In this case, the algorithm finds a solution that gives a results of $F\approx 0.9937$ for all values in the range $\omega_{a}^{(2)}/2\pi\in[4.80,4.90]$ GHz, shown in Fig. 11, for a time of $199$ ns. Performing again the robustness method with variations in the two-qubit parameters, a solution is found that is able to perform the desired entangling gate with a fidelity of $\mathcal{F}_{i}\approx 0.9639$ for all parameters in the range $\omega_{a}^{(2)}/2\pi\in[4.80,4.90]$ GHz and $\omega_{a}^{(1)}/2\pi\in[4.495,4.505]$ GHz, for a time of $199$ ns.

In the multilevel case, the same feature appears when attempting to a find a solution that is robust to fluctuations as in the two qubit parameters, i.e., that for just errors in $\omega_{a}^{(2)}$ the algorithm is able to find a solution with $\mathcal{F}>0.99$ for an error of $\pm 1\%$, but that if we wish to include $\omega_{a}^{(1)}$ it is more difficult to account for this. Nonetheless, with an error of $\pm 0.1\%$ in the parameter $\omega_{a}^{(1)}$ and $\pm 1\%$ in $\omega_{a}^{(2)}$ a pulse is found that achieves $F>0.96$.

One of the causes of the discrepancy between the fidelities of the two-level case and the multi-level case is the anharmonicity of the transmons we have simulated. Currently, we are operating deep in the transmon regime with $E_{j}/E_{c}=100$ and so one of the limiting factors is down to leakage out of the computational subspace as fluctuations in the qubits $\omega_{01}$ transition bring them even closer to the $\omega_{12}$ transition.

For this paper, we have chosen not to include errors in the coupling strengths between the transmons and the cavity $g_{j}$. Preliminary tests with errors in coupling strengths have shown minimal effect on the output given without including these errors up to $10\%$, while this range is feasible in experiments this would merely add to the parameter range selection in the simulations. Errors specifically in the cavity frequency have also not been included since the drive and qubits resonances are far off resonant from the cavity.

This paper focuses on achieving high-fidelity controls on time scales much shorter than decoherence times and therefore we do not include decoherence in our simulations. Obviously, it is important to correct control errors in this regime if they represent the largest source of infidelity. For example, a state-of-the-art circuit with two transmons with $50$-$\mu s$ coherence times acted on with an entangling unitary operation that requires $200$ ns will see a probability of corruption of the operation due to decoherence estimated at $0.4\%$. Since typical operations errors due unoptimized controls will be larger than this we can focus on optimizing without including decoherence. Usually this is the only relevant regime in which we will gain by optimizing. However, since some forms of decoherence can be tackled actively with dynamical decoupling schemes Viola and Lloyd (1998); Viola et al. (1999), it would be interesting to consider optimizing for Hamiltonian control errors and external decoherence together in the future.

V Conclusion

In conclusion, we have shown that robust quantum control can produce pulse shapes that achieve a desired unitary with high fidelity for a realistic quantum system. In particular, we have shown that in a system where a single-microwave drive coupled to a cavity containing two transmon qubits is chosen on resonance with a qubit, modifying the shape of the driving microwave pulse can produce a desired two-qubit interaction while mitigating the unwanted rotation of the qubit that is also on resonance with the drive. This can be done even in the presence of filtering on the control and errors in the system parameters with a modest amount of resources, and can be achieved even when realistic constraints are placed on the pulses.

We have seen that including constraints on the pulse opens up more areas of local maxima in the control space, and in this case we found that there were many. These “traps” may be what is limiting the range of robustness in the two-qubit frequencies. Future work will look at how to make the error range for $\omega_{a}^{(1)}$ larger, potentially by combining these methods with a non-local optimizer in order to circumvent local traps.

We have shown that it is still possible to achieve high fidelity control with reduced circuit complexity, by increasing the complexity of the control. This shows a trade-off between the circuit complexity and pulse complexity, and that as quantum computers grow we are likely to require more complex pulse shapes if we wish to keep the circuit complexity down. Future work will study the limitations of current algorithms to solving these problems.

Acknowledgements

E.G. acknowledges financial support from the EPSRC grants (Grants No. EP/L026082/1 and No. EP/L02263X/1). This work was partially supported by the KIST Institutional Program (Project No. 2E26680-16-P025). The data underlying this work are available without restriction. Details of the data and how to request access are available from the University of Surrey publications repository doi:10.15126/surreydata.00813890