CNOT gates

For the CNOT gates, we choose the operation time tf = 500ns, which is comparable to that of the fastest C-phase gate in the experiment [72]. And we choose the same number of control parameters kmaxfor both control pulses ΩX(t)and ΩY(t)in Eq. (3.27) and in Eq.

(3.28), and vary kmax= 7to kmax= 13for the J1 optimization. After the J1 optimization (100 random initial guesses), we show the optimized J1 values versus the optimization iterations (NM algorithm) in Figure 3.12. One can see that if kmaxis too small (kmax≤ 9), not all initial guesses can achieve the control parameter sets with the lowest J1 values after the optimization; if kmax is too large (kmax = 13), more optimization iterations are needed to achieve the control parameter sets with the lowest J1 values for some samples.

Therefore, we choose kmax= 11for the CNOT gate optimization. However, the lowest J1

values after the optimization are around 4 × 10⁻⁶ and these values can not be improved to arbitrarily small (to the machine limit) even if we increase kmax or use smaller time-step

for simulation. To exclude the approximation issue, we also do the J1 optimization by the ideal Hamiltonian in Eq. (3.1) without two approximations SWA and RWA, and observe the same result as the J1 optimization with SWA and RWA. Therefore, we think that the controllability with only ΩX(t)and ΩY(t)controls is not enough to fully control the system.

Even so, J1 is still over three orders of magnitude smaller than the threshold of surface codes 10⁻². Next, we see the distribution of the maximum pulse strength Ω^Max_X and Ω^Max_Y (the maximum values of |ΩX(t)| and |ΩY(t)| within the gate operation time tf) for the ensemble of J1optimized control parameter sets in Figure 3.13, and we observe that only a small portion of the ensemble with Ω^Max_X < 1mTand Ω^Max_Y < 1mT. Therefore, we increase the ensemble size of the J1 optimization to add more samples, satisfying Ω^Max_X < 1mTand Ω^Max_Y < 1mT, after the optimization. Before implementing the second step of the two-step optimization, we show the corresponding hJ2,U −i in Eq. (3.36) versus the corresponding Ω^Max_X and Ω^Max_Y for the J1 optimized ensemble (1100 samples) in Figure. 3.14. One can see that most lower values of hJ2,U −iappear in the region with Ω^Max_X > 1mTand Ω^Max_Y > 1mT, and it implies that stronger pulse strength has benet for suppressing the electrical noise.

However, for the maximum pulse strength constraint in the realistic system, we need to choose the initial guesses for the second step optimization from the J1 optimized control parameter sets in the region with Ω^Max_X < 1mTand Ω^Max_Y < 1mT(49 samples).

After the second step optimization, we rst lter out the optimized control parameter sets with Ω^Max_X ≥ 1mTor Ω^Max_Y ≥ 1mT, and then nd an optimal control parameter set in the remaining sets with J1 = 4.53 × 10⁻⁶and hJ2,U −i = 2.78 × 10⁻⁵ and with Ω^Max_X = 0.89mT and Ω^Max_Y = 0.98mT. To see the robustness against the uncertainty αt0 for this optimal control parameter set, we simulate the ensemble average indelity hIi with the realistic Hamiltonian in Eq. (3.22) (standard deviation of the electrical noise σU − = 3GHz; the cuto frequency for ltering eects, ω0/2π = 425.4MHz) without SWA and RWA as shown in blue circle-line in Figure 3.15. At the point αt0 = 0, J2,t0 = 0, our predicted ensemble average indelity should be hIi ∼= J1+ hJ2,U −i = 4.53 × 10⁻⁶+ 2.78 × 10⁻⁵ ∼= 3 × 10⁻⁵, but it contradicts with hIi ∼= 10⁻⁴ simulated by the realistic Hamiltonian in Figure 3.15.

This is because the ltering eects degrade J1 from 4.53×10⁻⁶to 10⁻⁴ as shown in Figure 3.11. Therefore, we use the optimal control parameter set after two-step optimization as the initial guess for the ne-tuning optimization. The hIi versus αt0 after the ne-tuning optimization is shown in the red diamond-line in Figure 3.15, and at the point αt0 = 0, hIi

10³ 10⁴ 10⁵ 10^-6

10^-5 10^-4 10^-3 10^-2 10^-1

Figure 3.12.: The optimized J1 values versus optimization iterations after the J1 optimiza-tion for the CNOT gate.

0 0.5 1 1.5 2 2.5 3

0 0.5 1 1.5 2 2.5 3

Figure 3.13.: The distribution of the maximum pulse strength Ω^Max_X and Ω^Max_Y after the J1

optimization for the CNOT gate.

Figure 3.14.: The hJ2,U −i topography of the CNOT gate with σU −= 3GHz versus max-imum pulse strength Ω^Max_X and Ω^Max_Y for the J1 optimized ensemble (1100 samples).

recovers to our original estimation ∼ 3 × 10⁻⁵. Besides, one can see that the robustness curves (hIi versus αt0) after the two-step optimization (blue circle-line) and after the

ne-tuning optimization (red diamond-line) coincide for larger αt0. This is because, for large αt0, hJ2,U −i and other higher-order noise contributions are larger than J1 and thus dominate in hIi, and hJ2,U −i is not sensitive to the ltering eects as shown in Figure 3.11 and other higher-order noise contributions could also be insensitive to the ltering eects, which results in the overlap of the two curves for large αt0. After the ne-tuning optimization, we obtain a new optimal control pulse, which can recover the J1degradation from the ltering eects, and the pulse shift from the original optimal control pulse after the two-step optimization is ∼ 10⁻³mT as shown in Figure 3.16.

Next, we compare the performance of our optimal CNOT gates with the maximum pulse strength Ω^Max_X and Ω^Max_Y smaller than 1mT (Ω<1mT) and smaller than 1.5mT (Ω<1.5mT), and the C-phase gate (simulation for the ideal C-phase realized in the experiment [72]) in Figure 3.17. To see the ability to suppress the static and stochastic electrical noise βU − with standard deviation σU − = 3GHz, let us take αt0 = 0. At αt0 = 0, the ensemble average indelity hIi of the optimal CNOT gate of Ω<1mT (red diamond-line)

-0.15 -0.1 -0.05 0 0.05 0.1 0.15 10^-5

10^-4 10^-3 10^-2 10^-1

Figure 3.15.: hIi of the CNOT gate versus αt0 after the two-step optimization and after the ne-tuning optimization.

0 100 200 300 400 500

-3 -2 -1 0 1 2 3 10^-3

Figure 3.16.: Pulse shift after the ne-tuning optimization for the CNOT gate.

-0.15 -0.1 -0.05 0 0.05 0.1 0.15 10^-5

10^-4 10^-3 10^-2 10^-1 10⁰

Threshold of surface codes

Figure 3.17.: Robust performance against uncertainty αt0 in t0 for the optimal CNOT gates of Ω<1mT(red diamond-line) and Ω<1.5mT(yellow square-line), and the C-phase gate (blue circle-line).

is improved near two orders of magnitude compared with that of the C-phase gate (blue circle-line). If the maximum pulse strength Ω^Max_X and Ω^Max_Y is relaxed to be smaller than 1.5mT (Ω<1.5mT), hIi of the optimal CNOT gate (yellow square-line) is improved more than two orders of magnitude. For the robust performance against the uncertainty αt0 in t0, the C-phase gate (blue circle-line) can be robust only to αt0/t0. 1% for the threshold of the surface codes (hIi . 10⁻²). Our optimal CNOT gate of Ω<1mT (red diamond-line) can be robust to αt0/t₀ ∼ 10%, and that of Ω<1.5mT (yellow square-line) further robust to α_t0/t₀∼ 15%. The corresponding optimal pulses after the ne-tuning optimization for the CNOT gate of Ω<1mT and the CNOT gate of Ω<1.5mT are shown in Figure 3.18.

3.4.2. Single-qubit gates

In this subsection, we demonstrate the performance of two single-qubit gates I2 ⊗ X₁ (Identity gate for dot2 qubit and X-gate for dot1 qubit) and H2⊗ I₁ (Hadamard gate for dot2 qubit and Identity gate for dot1 qubit). We suitably choose the number of control parameters kmax = 8 for both ΩX(t) and ΩY(t) of these two gates and choose the gate operation time tf = 200ns and tf = 250ns for I2⊗ X₁ and H2⊗ I₁, respectively. The

0 100 200 300 400 500 -1

-0.5 0 0.5 1

0 100 200 300 400 500

-1.5 -1 -0.5 0 0.5 1 1.5

Figure 3.18.: The optimal control pulses of the CNOT gate with Ω<1mT and Ω<1.5mT.

hJ_{2,U −}i topography of I2 ⊗ X₁ gate after the J1 optimization is shown in Figure 3.19, and we compare it with that of the CNOT gate in Figure 3.14. One can see that the height dierence in the hJ2,U −i topography of the CNOT gate is only around one order of magnitude, while around two orders of magnitude for I2⊗X₁gate. Furthermore, the lowest hJ_{2,U −}iarea (deep blue area) for I2⊗ X₁ gate is closer to the area with Ω^Max_X < 1mTand Ω^Max_Y < 1mT than that for the CNOT gate. That is, the hJ2,U −i topography of I2⊗ X₁ gate around the area Ω^Max_X < 1mTand Ω^Max_Y < 1mTis more steep than that of the CNOT gate. Therefore, for the second step optimization of I2⊗ X₁ gate, all initial guesses in the area Ω^Max_X < 1mTand Ω^Max_Y < 1mTow into the area ΩX,max > 1mTor ΩY,max > 1mT more easily than that of the CNOT gate does.

So, we should add an extra cost function, the uence (a measure of the eld energy) [43],

F ≡ ˆ _tf

0

|Ω_X(t)|²dt + ˆ _tf

0

|Ω_Y(t)|²dt (3.37)

in the second step optimization and in the ne-tuning optimization to modify the cost function topography in the control parameter space, and the total cost function becomes

J1+ hJ2,U −i + ξ · F , (3.38)

Figure 3.19.: hJ2,U −i topography of I2 ⊗ X₁ gate with σU − = 3GHz versus maximum pulse strength Ω^Max_X and Ω^Max_Y for the J1 optimized ensemble (2389 samples).

where the constant parameter ξ determines the contribution ratio of the uence F in the total cost function. If ξ is too small, F doesn't work and the steep hJ2,U −itopography still exists, while if ξ is too large, F dominates all the topography and thus hJ2,U −i can not be suppressed. And we nd the optimal ξ = 10⁻⁶. The performance of the optimal single-qubit I2⊗ X₁ gate and H2⊗ I₁gate is shown in Figure 3.20, and the corresponding optimal pulses after the ne-tuning optimization are shown in Figure 3.21. Both optimal gates can suppress the static and stochastic electrical noise with σU − = 3GHz to hIi ∼= 10⁻⁵ (at α_t0 = 0), and can be robust to αt0/t₀ more than 15% for the threshold of surface codes, hIi . 10⁻².

-0.15 -0.1 -0.05 0 0.05 0.1 0.15 10^-5

10^-4 10^-3 10^-2 10^-1

Threshold of surface codes

Figure 3.20.: Robust performance against uncertainty αt0 in t0 for the optimal single-qubit I2⊗ X₁ gate and H2⊗ I₁ gate.

0 50 100 150 200

-1 0 1

0 50 100 150 200 250

-1 0 1

Figure 3.21.: The optimal control pulses of the single-qubit I2⊗ X₁ gate and H2⊗ I₁ gate.

4. Conclusion

Our two-step optimization method can provide robust control pulses of high-delity quan-tum gates for stochastic time-varying noise and systematic error. Besides, our method is quite general, and can be applied to dierent system models, noise models, and noise CFs (PSDs). We apply our robust control method to the realistic system, quantum-dot electron spin qubits in isotopically puried silicon. We use the realistic system parameters from the experiment, characterize the noise model and noise strength from the experiment, and also consider experiment constraints such as the power limitation of the on-chip ESR line and the nite bandwidth of waveform generators, and nally demonstrate the high-delity and robust single-qubit gates and CNOT gates for this realistic system by our robust control method. Therefore, our method will make essential steps toward constructing high-delity and robust quantum gates for FTQC in realistic quantum computing systems.

When our optimal pulses are applied to the qubits in the laboratory, the gate delity could degrade from our prediction for some unknown factors without taking into account in our simulation model. Therefore, some closed-loop optimization methods [85, 86] imple-mented in the laboratory are suggested to calibrate our optimal pulses for recovering the

delity degradation from these unknown factors. The cost function for the closed-loop op-timization in the laboratory is just the ensemble average indelity hIi which is obtained via many repetitions of indelity measurement in the experiment. If the uctuations of these unknown factors are small, then hIi ∼= J1+ hJ₂i. If these unknown factors doesn't alter our original hJ2i topography simulation and change only J1, the function of the closed-loop optimization in the laboratory corresponds to that of our ne-tuning optimization. If these unknown factors destroy our original hJ2i topography simulation largely, our opti-mal control parameter set is no longer on the at bottom of the realistic hJ2itopography.

Thus the function of the closed-loop optimization in the laboratory is equivalent to that of our second step of the two-step optimization. To estimate how much time is required

for implementing the closed-loop optimization in the laboratory, we rst assume the cor-responding ne-tuning optimization via NM algorithm is implemented in the laboratory, and, in the case of our optimal CNOT gate (22 control parameters) in Sec. 3.4, around 1000optimization iterations are required for NM algorithm, where the average number of cost function calls per each iteration is ∼ 3 [87]. So we need to input the cost function to NM algorithm ∼ 3000 times, and the cost function in the closed-loop optimization is just the ensemble average indelity hIi. Each hIi, obtained by randomized benchmarking [88, 89], can be performed in 2 seconds in the laboratory of the superconducting qubit system [86], and thus we can perform the corresponding ne-tuning optimization of our optimal CNOT gate in the laboratory in ∼ 1.7 hours by ∼ 3000 hIi measurements. How-ever, the optimization iterations of the second-step optimization is over 500 times as that of the ne-tuning optimization. Therefore, around ∼ 35 days are required for the closed-loop optimization in the laboratory when the unknown factors change the original hJ2i topography simulation largely, and in this case we think the practical way is to character-ize these unknown factors in experiment, and then input the detailed information of these unknown factors to simulate the precise hJ2i topography for our two-step optimization via classical computers. Another improvement way is to use more ecient optimal algorithms other than NM algorithm to reduce the total number of measurements for the closed-loop optimization in the laboratory [90].

To conclude, the optimal control theory enables us to construct robust and high-delity gate pulses against noise and uncertainty in qubit systems. The optimized pulses af-ter closed-loop calibration and optimization in the experiments can implement desired quantum gates with target performance. Several advanced experiments using the optimal control pulses have been demonstrated [86, 91, 92, 93, 94]. Thus the optimal control the-ory is practical and applicable experimentally and can provide an essential input into the realization of large-scale FTQC.

.

A. : Derivation of Eqs. (2.12)-(2.14)

We present the derivation of Eqs. (2.12)-(2.14) and discuss the role of the extra term  in Eq. (2.12). Substituting the total system propagator in the Dyson expansion U(tf) = U_I(t_f) · (I + Ψ₁+ Ψ₂+ · · · )into the indelity denition I of Eq. (2.2), we obtain

I = J₁

−2

4ⁿRe{Tr[U_T^†U_I(t_f)]^?· Tr[U_T^†U_I(t_f) · (Ψ1+ Ψ2+ · · · )]}

−1 4ⁿ

Tr[U_T^†U_I(t_f) · (Ψ1+ Ψ2+ · · · )]

2

. (A.1)

The rst term J1on the right hand side of Eq. (A.1) is the gate indelity for the ideal system dened in Eq. (2.13). Then we dene the error shift matrix U of the ideal propagator U_I(t_f) at time tf from the target gate UT up to a global phase φ as

UI(tf) = e^iφUT(I + U). (A.2)

Note that when the gate indelity J1 for the ideal system is made small, the matrix elements of U also become small. Substituting the expression of UI(t_f)of Eq. (A.2) back to Eq. (A.1), we obtain

I = J₁+ {− 1

2ⁿ⁻¹Re[Tr(Ψ₁)]}

+J2+ (U, Ψj) + O( ˜H^m_N, m ≥ 3), (A.3)

where J2 is dened in Eq. (2.14),

(U_, Ψ_j) =

− 1

2ⁿ⁻¹Re{Tr[U_(Ψ₁+ Ψ₂+ · · · )]}

− 2

4ⁿRe{Tr[U]^?· Tr[Ψ₁+ Ψ2+ · · · ]}

− 2

4ⁿRe{Tr[U(Ψ1+ Ψ2+ · · · )] · Tr [Ψ1+ Ψ2+ · · · ]^?}

− 2

4ⁿRe{Tr[U_]^?· Tr[U_(Ψ₁+ Ψ₂+ · · · )]}

− 1

4ⁿ|Tr[U_(Ψ₁+ Ψ₂+ · · · )]|², (A.4) and O( ˜H^m_N, m ≥ 3)denotes other higher-order terms without containing U. The rst-order noise term, −Re[Tr(Ψ1)]/2ⁿ⁻¹, in Eq. (A.3) actually vanishes for Tr(Ψ1) = −i´_tf

0 Tr[HN(t⁰)]dt⁰ is a purely imaginary number, where the noise Hamiltonian HN(t⁰)is Hermitian [with βj(t) being real] and thus Tr[HN(t⁰)]is a real number. This result of no rst-order noise contri-bution in I is similar to that in Ref. [42]. This is also the reason why there is no rst-order noise contribution in ensemble average hIi of Eq. (2.15). Equations (2.12)-(2.14) can then be easily obtained from Eq. (A.3) with the identication of  = (U, Ψ_j).

We discuss below the property and the role of  = (U, Ψ_j)in Eq. (2.12) or in Eq. (A.3).

The extra contribution  = (U, Ψj) to the gate indelity, with the detailed form shown in Eq. (A.4), is related to the error shift matrix U and all Dyson expansion terms Ψj. As noted earlier, if J1 is small, then the matrix elements of U are also small. Moreover, if the noise strength is not too strong such that |Ψj+1|  |Ψ_j|, then the extra contribution

 = (U_, Ψ_j)is also small. Therefore when running optimization for a low noise strength, for which the higher-order term O( ˜H^m_N, m ≥ 3)becomes negligible (see Appendix B), the extra contribution  can be omitted as J1 is minimized to a small number. Consequently, one can focus on the optimization of only J1+ hJ₂i.

The advantage of introducing J1and  in our method is to enable more degrees of freedom in control parameters for optimization. There are actually no J1and  contributions in the gate indelity expression of the robust control method of SUPCODE [30, 31] and the lter-transfer-function method [49, 50]. In these methods, J1, or, equivalently, the error shift matrix U, is set exactly to 0 by imposing some constraints on the control parameters. In contrast, our method can tolerate some error of U and thus have more degrees of freedom

in control parameters as long as J1 and the extra contribution hi in gate indelity hIi are made just smaller than hJ2i. This advantage of having more degrees of freedom for optimization plays an important role in nding better control pulses as the number of qubits, the number of controls, and the number of noise sources increase.

B. : Estimation of higher-order contributions

Here we estimate the contributions of higher-order terms O( ˜H^m_N, m ≥ 3)and discuss when they can be neglected. We express the higher-order terms as O( ˜H^m_N, m ≥ 3) = P

p≥3Jp, where Jp denotes the p-th order noise term of the gate indelity. Detailed forms of the

rst two lowest-order terms in O( ˜H^m_N, m ≥ 3) are is the q-th order Dyson expansion term. To make an estimation of the magnitude of Ψq, we take the Z-noise model for the single-qubit gate operations in Sec. 2.3.1.1 as an example.

Substituting the noise Hamiltonian ˜H_N(t) = βZ(t)RZ(t)with RZ(t) = U_I^†(t)[ω0Z/2]UI(t) in the interaction picture into Ψq, we obtain

Ψ_q = (−i)^q The time integral contribution {´_tf

0 ω₀dt₁´_t1

0 ω₀dt₂· · ·´_tq−1

0 ω₀dt_q}can be estimated to be

about ∼ (ω0tf)^q/q!. By combining the above estimations, the magnitude of |Ψq,jk| is of the order of ∼ (ω0tfσZZ)^q/q!. Then substituting the estimated value of |Ψq,jk| into J₂ in Eq. (2.14), J3 in Eq. (B.1), and J4 in Eq. (B.2), we obtain the magnitude ra-tio J3/J₂ ∼ (ω₀t_fσ_ZZ)/3 and J4/J₂ ∼ (ω₀t_fσ_ZZ)²/12. The single-qubit gate operation time in Sec. 2.3.1.1 is ω0t_f = 20. If we choose the noise uctuation σZZ = 10⁻³, then the ratio J3/J₂ ∼ (6 × 10⁻³) and J4/J₂ ∼ (3 × 10⁻⁵), and thus the higher-order terms O( ˜H^m_N, m ≥ 3)can be safely neglected. If, however, σZZ ∼ 10⁻¹, then ω0t_fσ_ZZ ∼ 2. In this case, J3/J2 ∼ 2/3 and J4/J2 ∼ 1/3, so the higher-order terms O( ˜H^m_N, m ≥ 3) can not be neglected. Comparing our estimation with the results of the full-Hamiltonian simulation, one nds that the ensemble average of the gate indelity hIi of the IDG strategy scales as the second power of σZZ (because hJ2i dominates) for small σZZ until σZZ ∼ 10⁻¹ in Figure 2.3(a) for low-frequency noise γZZ = 10⁻⁷ω₀and in Figure 2.4(a) for high-frequency noise γZZ = 10⁻¹ω₀. This is consistent with our estimation. In other words, if σZZ is con-siderably smaller than 10⁻¹, O( ˜H^m_N, m ≥ 3) can be ignored. Therefore, even in the case where the full-Hamiltonian simulation is not available, we can use this estimation method to determine the criterion for neglecting the higher-order terms O( ˜H_N^m, m ≥ 3).

Bibliography

[1] A. G. Fowler, M. Mariantoni, J. M. Martinis, and A. N. Cleland, Surface codes:

Towards practical large-scale quantum computation, Phys. Rev. A, vol. 86, p. 032324, Sep 2012.

[2] K. Khodjasteh and L. Viola, Dynamically error-corrected gates for universal quantum computation, Phys. Rev. Lett., vol. 102, p. 080501, Feb 2009.

[3] K. Khodjasteh and L. Viola, Dynamical quantum error correction of unitary opera-tions with bounded controls, Phys. Rev. A, vol. 80, p. 032314, Sep 2009.

[4] K. Khodjasteh, D. A. Lidar, and L. Viola, Arbitrarily accurate dynamical control in open quantum systems, Phys. Rev. Lett., vol. 104, p. 090501, Mar 2010.

[5] K. Khodjasteh, H. Bluhm, and L. Viola, Automated synthesis of dynamically cor-rected quantum gates, Phys. Rev. A, vol. 86, p. 042329, Oct 2012.

[6] J. R. West, D. A. Lidar, B. H. Fong, and M. F. Gyure, High delity quantum gates via dynamical decoupling, Phys. Rev. Lett., vol. 105, p. 230503, Dec 2010.

[7] A. M. Souza, G. A. Álvarez, and D. Suter, Experimental protection of quantum gates against decoherence and control errors, Phys. Rev. A, vol. 86, p. 050301, Nov 2012.

[8] F. F. Fanchini, J. E. M. Hornos, and R. d. J. Napolitano, Continuously decoupling a hadamard quantum gate from independent classes of errors, Phys. Rev. A, vol. 76, p. 032319, Sep 2007.

[9] F. F. Fanchini, R. d. J. Napolitano, B. Çakmak, and A. O. Caldeira, Protecting the

√swap operation from general and residual errors by continuous dynamical decou-pling, Phys. Rev. A, vol. 91, p. 042325, Apr 2015.

[10] A. Z. Chaudhry and J. Gong, Decoherence control: Universal protection of two-qubit states and two-qubit gates using continuous driving elds, Phys. Rev. A, vol. 85, p. 012315, Jan 2012.

[11] X. Xu, Z. Wang, C. Duan, P. Huang, P. Wang, Y. Wang, N. Xu, X. Kong, F. Shi, X. Rong, and J. Du, Coherence-protected quantum gate by continuous dynamical decoupling in diamond, Phys. Rev. Lett., vol. 109, p. 070502, Aug 2012.

[12] J. Clausen, G. Bensky, and G. Kurizki, Bath-optimized minimal-energy protection of quantum operations from decoherence, Phys. Rev. Lett., vol. 104, p. 040401, Jan 2010.

[13] J. Clausen, G. Bensky, and G. Kurizki, Task-optimized control of open quantum systems, Phys. Rev. A, vol. 85, p. 052105, May 2012.

[14] M. Wenin and W. Pötz, Minimization of environment-induced decoherence in quan-tum subsystems and application to solid-state-based quanquan-tum gates, Phys. Rev. B, vol. 78, p. 165118, Oct 2008.

[15] P. Rebentrost, I. Serban, T. Schulte-Herbrüggen, and F. K. Wilhelm, Optimal

[15] P. Rebentrost, I. Serban, T. Schulte-Herbrüggen, and F. K. Wilhelm, Optimal