On The Power Of Coherently Controlled Quantum Adiabatic Evolutions

A major challenge facing adiabatic quantum computing is that algorithm design and error correction can be difficult for adiabatic quantum computing. Recent work has considered addressing his challenge by using coherently controlled adiabatic evolutions in the place of classically controlled evolution. An important question remains: what is the relative power of controlled adiabatic evolution to traditional adiabatic evolutions? We address this by showing that coherent control and measurement provides a way to average different adiabatic evolutions in ways that cause their diabatic errors to cancel, allowing for adiabatic evolutions to combine the best characteristics of existing adiabatic optimizations strategies that are mutually exclusive in conventional adiabatic QIP. This result shows that coherent control and measurement can provide advantages for adiabatic state preparation. We also provide upper bounds on the complexity of simulating such evolutions on a circuit based quantum computer and provide sufficiency conditions for the equivalence of controlled adiabatic evolutions to adiabatic quantum computing.

The quantum adiabatic theorem is an essential tool for quantum information processing and quantum control [1][2][3][4][5][6]. It states that the evolution generated by a slowly varying Hamiltonian (relative to the minimum eigenvalue gap) maps eigenstates of the initial Hamiltonian to eigenstates of the final Hamiltonian [7]. This process provides a simple and error robust method for state preparation that is used extensively in quantum simulation, adiabatic quantum computing as well as pulse design. A drawback of adiabatic evolution is that it is often much slower than competing state preparation methods. Finding ways of reducing "diabatic" errors (which result from using a finite evolution time) is vitally important for practical applications of adiabatic state preparation.
Two major strategies have been proposed to minimize the error in adiabatic evolutions: local adiabatic evolution and boundary cancellation methods. Local adiabatic evolution [8,9] (LAE) minimizes the time to reach the adiabatic regime by choosing the evolution speed such that the adiabatic condition is satisfied at each instant throughout the evolution. In a typical scenario of LAE, the rate at which the Hamiltonian changes is fast in the beginning and the end of evolution, when the distance between the ground state and the first excited state is large, and small in the middle around the minimal gap. This approach optimizes the scaling of the evolution time with the size of the system and works best to reduce diabatic errors in the short time or "Landau-Zener" regime (so called because the Landau-Zener formula provides a better approximation to the resultant state than adiabatic perturbation theory does).
Boundary cancellation methods minimize the error in the adiabatic approximation once the adiabatic condition is met [10,11]. These methods polynomially improve the error scaling, relative to LAE, by setting the first n − 1 derivatives of the Hamiltonian to zero at the boundaries. This strategy tends to lead to taking the Hamiltonian to be slowly varying near the beginning and end of the evolution, which typically is where the eigenvalue gap is largest. Since the Hamiltonian will often vary slowly when the gap is large, it forces the evolution to speed up around the minimal gap, which retards the convergence to the adiabatic regime (the regime where adiabatic perturbation theory applies).
These two approaches are typically at odds: LAE says that you should move quickly when the gap is large to minimize the error, which is often at the beginning and end of the adiabatic passage [8,9]; whereas boundary cancellation methods show that it is often best to move very slowly at the beginning and end of the evolution [10,12]. The question is: can these two objectives be simultaneously satisfied and if so how ? We consider a model of adiabatic quantum computation that can achieve both goals. Our hybrid model for adiabatic computation uses a small control register that the user has universal control over, and a larger adiabatic system that is coherently controlled by the smaller register. These generalizations grant us two freedoms: (a) the adiabatic subsystem can evolve under a superposition of different adiabatic evolutions and (b) measurement can be used on the control qubits without exciting the system out of the groundstate. These freedoms allow us to escape the constraints of unitarity and implement a wider class of operations including linear combinations of unitaries [13], which we use to increase the resilience of the evolution to diabatic errors. This model also subsumes those of [14,15].
Unlike the previous methods, we do not search for a single optimal adiabatic evolution. Instead we take two (or more) evolutions that generate errors that are oriented in opposite directions as in Figure 2 and then use the non-deterministic circuit in Figure 1 to suppress these errors by performing an appropriate weighted average of the evolutions. We then show that a linear combination of adiabatic evolutions can asymptotically decrease the error in the adiabatic approximation. The resultant averaged adiabatic evolution can have the benefits of both LAE and boundary cancellation methods: the convergence to the adiabatic regime is comparable to LAE, while the error scaling in the adiabatic regime is comparable to that of boundary cancellation methods.
In the following section we review the adiabatic theorem. We then provide the gadget that we use to cancel the leading order diabatic errors in Section II. We illustrate the utility of this method in Section III where we apply the gadget to approximately cancel the dominant diabatic transition. We provide methods in Section V that simultaneously suppress every transition, assuming that the adiabatic paths obey a particular symmetry condition. Finally, we discuss how our techniques can combine the best features of local-adiabatic evolution and boundary cancellation methods in Section IV and then discuss the implementation of our model of coherently controlled adiabatic evolution using a quantum computer in Section VII.

I. REVIEW OF THE ADIABATIC THEOREM
It is not possible to provide a closed form solution to the Schrödinger equation for the case of time-dependent Hamiltonians in general. It is customary in such cases to express the time evolution operator, which is the formal solution to as A wide array of approximation methods exist for U (t) including the Magnus expansion [16], Dyson series [17], Floquet theory [18], the Landau-Zener formula [19] and the adiabatic approximation. The adiabatic approximation is widely used to approximate quantum dynamics in cases where rate of change of the time-dependent Hamiltonian is slow relative to an appropriate power of the minimum eigenvalue gap. In essence, the approximation states that if you prepare a system in an eigenstate of the Hamiltonian and evolve sufficiently slowly then the quantum system will remain in the eigenstate throughout the evolution. This lack of excitation throughout the process makes it analogous to reversible adiabatic processes in statistical mechanics. This analogy is not exact since the change in Von-Neumann entropy is also zero for any unitary process and so the "adiabatic" moniker persists for largely historical reasons.
Since the adiabatic approximation requires slow evolution, it is useful to consider how the approximation error scales as the speed of the transition from the initial to the final Hamiltonian decreases. This makes it natural to parameterize time via the variable s where and T is the total time for the adiabatic passage. While an adiabatic evolution occurs on t ∈ [0, T ], s ∈ [0, 1] regardless of the actual duration of the evolution. This means that if the Hamiltonian is re-parameterized as H(s), then we can increase T to make the evolution slower without fundamentally changing the form of the evolution. We need to introduce some further notation before we can discuss the adiabatic approximation in greater detail. Firstly we define |n(s) to be the instantaneous eigenvectors of the time-dependent Hamiltonian: and we make no assumptions about the ordering of E n (i.e. we do not assume that E 0 ≤ E 1 ). Also for notational simplicity, we define |g(s) := |E 0 (s) . We refer to this state as |g(s) because it will represent the ground state in many practical examples of adiabatic QIP. The eigenvalue gaps will also be key to our analysis and so we use the following notation for them: The adiabatic approximation is often expressed in many different ways. The simplest of these states that In general the adiabatic approximation holds if for integers a, b, c and d that depend on the properties of the Hamiltonian [20][21][22][23]. A common misconception is that the adiabatic approximation holds if although this criteria is appropriate for sufficiently slow evolutions under smoothly varying Hamiltonians [20,21].
We refer to such results as zeroth order adiabatic theorems, because they provide an estimate of the error that is correct to zeroth order in T −1 , meaning that they simply tell you that the error is zero if the adiabatic process is infinitely slow. In order to show that we can combine different adiabatic evolutions together to cancel the error, we need to have a sharper adiabatic condition that approximates the error to at least O(1/T ). It is necessary for us to use a first order adiabatic approximation, which provides us with the error in the adiabatic approximation correct to O(1/T 2 ) [20]: This result can easily be found by using path integral methods using techniques also discussed in [24][25][26] and upper bounds on the magnitude of the sum of all O(1/T 2 ) terms are given in [20]. Equation (9) tells us something surprising: the leading order contribution to the error in the adiabatic approximation does not depend on the minimum gap but rather the eigenvalue gap at the beginning and the end of the evolution, which motivates takingḢ(s) = 0 or equivalently ṅ(s)|g(s) = 0 on the boundary as per boundary cancellation methods [10]. The apparent contradiction posed by (9) is easily resolved. Adiabatic conditions like (7) give criteria for the convergence of the adiabatic perturbation series of U (T, 0) in powers of 1/T and equations such as (9) give a truncated expression for the power series. This means that after a critical evolution speed, the error in the adiabatic approximation no longer depends on the minimum gap; whereas, the error depends crucially on the minimum gap before this point. We refer to the regime where the minimum gap dictates the error as the Landau-Zener regime and the regime where it does not as the adiabatic regime.
Similarly, the first order adiabatic theorem relies on several conditions outlined in [20]. First, the Hamiltonian must be twice differentiable and 3-times piecewise differentiable with all such derivatives upper bounded by a constant. Second, the system must be already in the adiabatic regime (i.e. the Θ(1/T ) contribution to the error is much greater than the sum of all O(1/T 2 ) contributions). Third, we require that the norm of the Hamiltonian be upper-bounded by a constant for all times during the evolution. These criteria guarantee the validity of (9).
A common way to reduce errors in both the Landau-Zener regime, as well as the adiabatic regime, is to change the path used in the adiabatic evolution. The most frequently used adiabatic path, known as linear interpolation, is where H 0 is the initial Hamiltonian and H 1 is the final Hamiltonian. There are of course many ways that one could imagine transitioning from the initial Hamiltonian to the final Hamiltonian. Each of these ways represents a particular "adiabatic path" and (10) is known as the linear adiabatic path. More generally we could consider a path of the form where f (0) = g(0) = 0 and f (1) = g(1) = 1. Such paths can be extremely important for adiabatic quantum computing because they allow the evolution to slow down through, or even avoid, parts of the evolution that contribute substantially to the error; however, here we assume the simple case of g(s) = f (s). We do not require that the range of f is [0, 1] here. In fact, some of the adiabatic paths that we consider will attain negative values and values greater than 1.
Other examples of non-linear paths include local-adiabatic evolution, which seeks to minimize the error in the Landau-Zener regime by choosing the evolution speed to be smallest near the minimum gap. Boundary cancellation methods on the other hand choose paths that minimize the error in the adiabatic regime by choosing the evolution speed to be zero at s = 0 and s = 1. These two strategies are seemingly orthogonal. At present there is no known method that combines the best features of local adiabatic paths and the paths yielded by boundary cancellation methods. Our work provides a way to achieve this, thereby illustrating that controlled adiabatic evolution affords greater power than conventional adiabatic evolution.

II. CONTROLLED ADIABATIC EVOLUTION USING A SMALL NUMBER OF ANCILLAS
The central idea behind our approach is to use a gadget that was recently proposed in [27] to non-deterministically implement the weighted average of two or more adiabatic evolutions. This idea of using controlled adiabatic evolutions and measurement has been recently explored by Itay Hen [14] and is also used in holonomic quantum computing [15]; however, these results do not consider using coherent control and measurement to suppress diabatic errors. The gadget that we use for this averaging process is given in Figure 1. The circuit in Figure 1 probabilistically implements linear combinations of unitary operations, as seen through the following argument We then see that if the ancilla register is measured to be 0 then the circuit performs a weighted combination of U A and U B on the state |ψ otherwise the circuit implements the difference between the two operators.
The generalization to cases where multiple U A and U B are used is trivial, it simply involves increasing the number of qubits used to control the overall rotation [27]. Such circuits, or variants thereof, are also used in [28,29]. For the case of adiabatic evolution, we know that to zeroth order where T = max{T A , T B }. This means that, to leading order, both U A and U B generate the same evolution up to a global phase and hence we expect the success probability to be high if the phases picked up by |g under both evolutions are comparable. Rather than choosing different paths that apply the same phase to |g(s) , we counter-rotate the evolution of each eigenstate by including an additional phase to each unitary. This affords us much greater freedom to choose adiabatic paths for U B and U A . In particular, we choose these phases such that We see from the choices of phases in (15) that (13) gives the failure probability of the linear combination O(1/T ) in the limit of large T . This means that the failure probability will typically be extremely small for adiabatic processes and even if a failure is observed, the gadget in Figure 1 informs the user that a failure has occurred and the state preparation process can be repeated until success is obtained. We see from numerical experiments that the failure probability of these circuits has a near-negligible impact on the cost of coherently controlled adiabatic state preparation in the adiabatic regime. Generalization of these ideas to cases where more than two unitary evolutions are averaged is straight forward and is discussed in detail in [13]. We present the two-unitary case explicitly here since the majority of our results focus on averaging two different adiabatic evolutions.

III. A GENERAL METHOD FOR CANCELING A SINGLE TRANSITION
Our first approach is a generalization of the strategy employed by Wiebe and Babcock in [30], which suppresses the dominant transition in the adiabatic passage for adiabatic paths satisfying by choosing the evolution time T appropriately. Our strategy is to suppress a single transition, not by choosing a single time and requiring a symmetry condition as per [30], but by interfering the adiabatic evolution with a dual evolution as suggested in Figure 2. This allows such errors to be suppressed for any evolution time and any primary path. We also provide a method for suppressing the two most significant diabatic transitions in Appendix A. We wish to choose, for fixed H A , an adiabatic path that quadratically suppresses the transition |g(0) → |e(1) where |e(s) is any given instantaneous eigenstate of H(s). From (9) and (15), we see that if we combine U A (T A , 0) with U B (T B , 0) and achieve a successful measurement outcome then we obtain a result proportional to So to leading order, the linear combination will give the correct result. Then using, (9) it is clear that This transition can therefore be canceled, to O(1/T 2 ), by choosing θ, T B and f B such that the weighted average of the diabatic transitions to the state |e is zero: where ė(s)|g(s) = e(s)|Ḣ(s)|g(s) γe,g(s) . Thus it is reasonable to expect that this condition can be met by choosing θ and f B properly. The remaining question is: how can this be done in practice? We provide two strategies for finding f B for any fixed f A such that these errors cancel to leading order.

A. Partially Anti-Symmetric Combination
Our first method chooses the paths f A and f B to satisfy an anti-symmetric condition on the derivatives at the beginning and end of the evolution. This approach is most useful in cases where it is desirable for f B to be as similar to f A as possible. When optimizing these paths it is important to note that although f A and f B are arbitrary interpolation functions that describe the adiabatic paths, they are constrained to obey Furthermore, let us choose f B such that its derivatives are symmetric with f A at s = 0 and anti-symmetric at s = 1 Then using (22), (19) simplifies to Equation (23) can be satisfied for any f A and T A by setting This solution reduces to that of [30] in the limit as T A → 0; however, a non-trivial secondary path will always be needed if the symmetry condition demanded by [30] is not held. An important consequence of taking the derivatives to be negative at s = 1 is that there exists s such that f (x) > 1 for all x ∈ (s , 1). This is a consequence of the fact that H(s) is twice differentiable and hence f B is continuous; from which the result directly follows from the mean value theorem. Thus f b does not monotonically approach 1 as s → 1, but rather it overshoots the value and then reverses direction to end the evolution at s = 1. Such reversals of direction are analogous to the backwards time steps used in Trotter-Suzuki methods and, although non-traditional, are not necessarily problematic for adiabatic evolution.
We see from this discussion that controlled adiabatic paths can be used to suppress diabatic errors in ways that are impossible using traditional adiabatic optimization strategies. In particular, for any optimization strategy, such as local adiabatic evolution, we can always find a second path to add to the primary path to suppress a chosen transition to one order higher. These ideas can also be generalized to suppress more than one transition; however, finding a closed form solution is difficult in such cases. We discuss generalizing this method to simultaneously suppress two diabatic transitions in Appendix A. A drawback of this approach is that it cannot be used for arbitrary small times because (24) forces a difference at least π 1 0 γg,e(ξ) between evolution times. We address this issue below by providing a method that does not require a shift in time, but requires a more substantial deformation to the primary adiabatic path.

B. Completely Anti-Symmetric Combination
An alternative approach is to set the derivatives for the second path to be completely antisymmetriċ Pluging (26) into (19) we obtain .
The error is suppressed when (25) holds and In other words, the gap integrals for both paths must be equivalent modulo 2π. This removes the difficulty with offsetting one of times. However, in this case, the path f B both begins and ends the evolution by moving backwards. Alternatively, we can modify one boundary from each path. This backwards motion at s = 0 means that there exists δ > 0 such that the range of f B (s) is within [−δ, 1 + δ]. Again, this use of backward evolution is atypical of conventional approaches to adiabatic evolution where the additional evolution time/speed required by backwards evolution would tend to be detrimental. In contrast, such backwards evolutions can lead to substantial reductions in the cost for coherently controlled adiabatic evolution.

C. Interpolation
There are many ways that these requirements can be satisfied by a dual path to f A . The way that we satisfy these requirements is by adding a smooth polynomial continuation of f A about s = 1 that allows the derivative to loop around and attain the opposite value. This interpolation must have piecewise continuous third derivatives in order to guarantee that the O(1/T 2 ) terms will remain sub-dominant in the limit of large T . This naturally leads to a quartic interpolation that takes the following form for a partially anti-symmetric combination where ∆ is a free parameter that controls when f B switches from the original adiabatic path f A to the polynomial interpolation. The parameters are then set by requiring We could also have equally well choose the backwards evolution to start at s = 0 rather than s = 1. Although seemingly arbitrary, this choice can have a substantial impact on the error depending on whether the gap is larger at s = 0 and s = 1. We also make use of this fact later in Section V where we exploit this fact to suppress every transition simultaneously for Hamiltonians that satisfy a certain symmetry property. The case of fully anti-symmetric boundaries is similar except now two polynomial interpolations are needed: where ∆ is a free parameter that controls how rapidly f B switches from the original adiabatic path, f A , to the polynomial interpolation. The parameters are then set by requiring In particular, the coefficients in (31) can then be found by substituting (31) into (32). It then follows that for any fixed path f A , we can choose f B such that the dominant transition is suppressed to O(1/T 2 ). This opens the possibility that our error suppression methods may allow adiabatic state preparation to be performed using less evolution time (or equivalently, fewer gates on a quantum computer) than existing methods. However, local adiabatic evolution is known to be optimal for performing adiabatic Grover's search [8,9], so we cannot expect that the algorithm will outperform all existing adiabatic algorithms in every time regime. We will see below, that although our method does not outperform local adiabatic evolution for short times it can come very close to matching its performance while giving substantially reduced error for slow evolutions.

IV. COMPARISON TO LOCAL ADIABATIC EVOLUTION
We focus our numerical results on the case of adiabatic Grover's search. The Hamiltonian for adiabatic Grover's search is and sufficiently slow evolution of this Hamiltonian causes the initial eigenstate |+ ⊗n to transition to the marked state |0 as per Grover's search. Local adiabatic evolution is known to be optimal for adiabatic Grover's search [8,9] meaning that the quadratic speedup over classical algorithms is attained for the adiabatic path. The path for local adiabatic evolution, for cases where the search space is N -dimensional, is Our goal is to compare the cost of performing this adiabatic quantum algorithm using local adiabatic evolution to the cost incurred by using our methods. We choose f A to be the path given by local adiabatic evolution whereas f B is taken to be a continuation of the local adiabatic path that satisfies (22) or (26) in all of the following numerical examples.
There are several ways that the cost of an adiabatic algorithm can be measured. The most straight forward method is to compare the time required for the evolution. Although this cost metric is appropriate in cases where the norm of the Hamiltonian is fixed, it is not appropriate for comparing different adiabatic evolutions because the energy required for both paths may differ substantially. Since there is a duality between energy and time in quantum mechanics, a fast evolution that requires a lot of energy may be dynamically equivalent to a slow evolution that requires little energy. Thus we need to consider not just the time but also the energy. For this reason, we use the following cost metric for the case where we combine j evolutions (where j ≥ 1): where P (0) is the success probability of the gadget which is given by (13). Here we implicitly assume that the cost of the rotations and the control logic is negligible and that each of the evolutions can be implemented in parallel. These assumptions may not hold in general, but they are appropriate for quantum computer simulations of such adiabatic evolutions because the query complexity of such evolutions depends on the maximum evolution time chosen rather than the total evolution time. This point will be made clear in Section VII. We see in Figure 3 that including the second adiabatic path with partially anti-symmetric boundary conditions (as per Section III A) to LAE yields comparable performance to LAE for short evolutions and also provides the improved scaling of boundary cancellation methods in the adiabatic regime. In particular, the second path follows the interpolation strategy of (29): it follows LAE (i.e. (34)) until s = 0.8 and then smoothly transitions to a fourth-order polynomial. Unlike the method of [30], this approach does not only give superior error scaling over a discrete set of points; although, this technique does enforce a minimum evolution time as discussed in Section III A. An important drawback of this method is that there is a manifest lack of symmetry in the derivatives in the adiabatic regime for this method. This means that the adiabatic interference effects that appear in the LAE and boundary cancellation paths will not appear here [30]. Note that if the Search Hamiltonian did not have symmetric derivatives or spectrum, the adiabatic interference effects would not appear and so they are an artifact of having a highly structured test Hamiltonian. Figure 4 tells a similar story. In that case we use fully anti-symmetric boundary conditions and add a second path that interpolates between LAE and polynomial evolution as per (31) with ∆ = 0.2. This also corresponds to evolution under LAE for 80% of the time. Unlike the case in Figure 3 adiabatic interference patterns are again visible in the adiabatic regime because the two polynomials used to create the fully anti-symmetric boundary conditions between the two paths at s = 0 and s = 1 ensure that the derivatives are the same at the boundary; thereby allowing such interference effects to emerge again which causes the error to be substantially reduced on a discrete subset of points as per [30]. As a consequence we can clearly see that our method substantially outperforms both methods for a range of evolutions with cost ranging from [50, 100] due in part to the presence of adiabatic interference effects that are absent from the boundary cancellation method.
In both of the cases considered, our methods are less effective at suppressing errors in the adiabatic regime than boundary cancellation methods. This is because the O(1/T 2 ) terms in the error in the adiabatic approximation also depend onḢ(s). Such terms are zero for boundary cancellation methods and so we generically expect from the triangle inequality that boundary cancellation will lead to less error in this regime. An important point to note is that although these test cases do not outperform LAE for fast evolutions or boundary cancellation methods for slow evolutions, they can outperform both methods for evolutions that operate at an intermediate speed. This implies that these methods are not just a compromise between the two approaches: they also provide superior scaling in a region that is badly addressed by existing adiabatic optimization methods.

V. SUPPRESSING EVERY TRANSITION FOR SYMMETRIC H
We now consider suppressing errors for Hamiltonians where H(0) and H(1) have the same spectra. Although restrictive, this condition is satisfied in many natural problems [3,8,30]. This symmetry is very useful because it guarantees that two adiabatic interpolations exist between H 0 and H 1 such that the amplitudes for every state orthogonal to |g(1) that arise due to f A are equal and opposite to those that arise under f B . This means that the linear combination will simultaneously suppress diabatic leakage into every state. In contrast, the methods discussed in Section III A and Section III B only guarantee this for a single (but arbitrarily chosen) transition. First let us assume that the following conditions are met for f A (s) and f B (s) for all states |n = |g . We will see that these conditions can always be met if the spectrum of H(s) is symmetric about s = 1/2. Such conditions do not naturally arise for all adiabatic passages but there are many examples where such Hamiltonians are natural. A natural example is the Search Hamiltonian; however, such an application is trivial because the quantum dynamics occurs within a two-dimensional subspace. Other examples occur in adiabatic gates [14,31,32] and holonomic quantum computing [15,33].
It is then clear that if we take θ = π/4 and T A = T B then every transition will be suppressed from O(1/T ) to O(1/T 2 ) under these assumptions. The question remaining is: when can we make these conditions hold? A natural case that covers a wide range of adiabatic protocols is the case where the eigenvalue gap is symmetric. That is to say that γ g,n (s) = γ g,n (1 − s) for all s and n = g. It is difficult to find a second adiabatic path that satisfies the conditions in Section III C, for the choice f A = f because the anti-symmetry required by (37), (37) necessitates the use of adiabatic paths similar to those in Section III B. Such paths will typically violate (36) because including the reversal near s = 0 and s = 1 will change the gap integral.
A better approach is to modify both paths. It is easy to see by substitution that if we let f B be given by (29) and (30) and then choose f A to be the time reversed version of this path (i.e. f A (s) = 1 − f B (1 − s)). It is then easy to see by substitution that under there assumptions The assumption that γ g,n (s) = γ g,n (1 − s) then directly implies that γ g,n (f A (s)) = γ g,n (f B (1 − s)), which gives us the desired result of This fact becomes immediately obvious in light of Figure 5, where the spectrum for a Search Hamiltonian with f A and f B chosen to be time reverses of each other is given. After substituting (41) into (39) and using the assumption that γ g,n (1) = γ g,n (0) we see that Thus for any adiabatic path parameterized by f (s) and any evolution time T , we can always choose two paths f A and f B that both incur diabatic errors that cancel to O(1/T 2 ). This fact is demonstrated numerically for a Hamiltonian that satisfies these requirements in Figure 6. In contrast, without these assumptions, only the dominant transition is suppressed to O(1/T 2 ) which implies that the diabatic errors scale as O(1/T ) for sufficiently long evolutions. This result is much stronger than that of [30], which only leads to suppression of all diabatic errors if the gap integrals for each transition are rational multiples of each other and, even then, will only work at specially chosen values of T . This precludes the technique's use for almost all Hamiltonians. In contrast, coherent control of the adiabatic path allows all of the transitions to be suppressed for a wide class of adiabatic protocols and this result holds for any T . This clearly shows demonstrates that coherently controlled adiabatic evolution allows us to circumvent the limitations of existing adiabatic optimization schemes. We will also see this below, where we show that these methods can be used in concert with boundary cancellation methods.

VI. INCORPORATING BOUNDARY CANCELLATION
The phase cancellation method allows one to reduce the error in adiabatic evolution to O 1 T m+1 by setting first m derivatives of the Hamiltonian zero at the boundaries as shown in [10,12]. Our approach can incorporate these results in a natural way. When we set first m derivatives to zero on the boundary, the remaining error is [30] ( (43) This expression is analogous to (9), as can be seen by substituting into (9). Equation (43) implicitly assumes that ṅ(s)|n(s) = 0, which can always be assumed to be true because the phases of |n(s) are arbitrary. The previous methods can then be used after making these substitutions and using a higher-order polynomial to perform the interpolation. For example, we can generalize the method of Section V with the following modification to the conditions required for both f A and f B : and picking T A = T B with θ = π/4. This enables exponentially accurate adiabatic approximations if m is chosen as a function of T . Alternatively, cancellation of the leading order transition to O(1/T m+1 ) can be obtained by using the exact same ideas within the methods of Section III A and Section III B after using the substitutions in (44) and conditions similar to (45) for the polynomial interpolations. It should also be noted that in systems where the adiabatically transported subspace is a one-or two-dimensional subspace then such approaches also can be used to make the overall error scaling O(1/T m+1 ). This raises the possibility that controlled adiabatic evolution can be combined with boundary cancellation methods to substantially reduce the cost of performing a high-accuracy adiabatic state preparation.

VII. COSTING CONTROLLED ADIABATIC EVOLUTIONS
There are two types of costly resources in coherently controlled adiabatic evolutions. The first is the cost of evolving the adiabatic register, denoted |ψ in Figure 1. The second is the cost of performing the required rotations on the control register. In this section we will provide a complete cost analysis of this model under the assumption that it is being simulated using a circuit based quantum computer that is further equipped with oracles that compute the necessary properties of the Hamiltonian. We will then conclude that coherently controlled adiabatic using sparse, rowcomputable Hamiltonians evolution is polynomially equivalent to the circuit model. Other appropriate cost models, such as bounding the energy and time required to implement the controlled Hamiltonians using k-local Hamiltonians will not be discussed here. Since we show that coherently controlled adiabatic evolution is polynomially equivalent to the circuit model, it will immediately follow that it is also polynomially equivalent to adiabatic quantum computation using local Hamiltonians.
The first important result that we need to show this is an upper bound on the number of oracle queries needed to simulate a time-dependent Hamiltonian within fixed error on a quantum computer. We will use this result to upper bound the query complexity of performing the controlled adiabatic evolutions. In order to understand the theorem, we will define a smoothness classification for Hamiltonians: , for all t ∈ I and p ∈ {0, 1, · · · , P }.
Now with Definition 1 we can state the following corollary, which gives the query complexity of the simulaiton algorithm. The number of one-and two-qubit needed for the simulation is at most proportional to the number of queries made.
Corollary 1 (Cor. 6 of [34]). Let {H µ : R → C 2 n ×2 n ; µ = 1, . . . , M } be a set of time-dependent Hermitian operators that is Λ-2k-smooth on I = (t 0 , t 0 + ∆t) \ {t 1 , . . . , t L }, where t 0 < t 1 < · · · < t L < t 0 + ∆t, with the additional conditions where N T is the number of bits used to represent t and n H is the number of qubits used to encode the matrix elements of H. Then the query complexity for simulating evolution generated by H(t) = µ T † µ H µ (t)T µ , for fixed set of unitary basis changing operators T µ , within an error of using time-ordered Trotter-Suzuki formulas with error O(∆t 2k+1 ) is where C is the number of oracle calls needed to simulate a one-sparse Hamiltonian, and the number of basis change operations is at most N queries /(3Cd 2 ).
Note that we need to use a result for simulating piecewise smooth Hamiltonians because the interpolations used in our method will typically cause H(f j (s)) to be non-analytic at either one or two points. The result of Corollary 1 is then useful because it provides the cost of performing such a simulation despite such complications. Note that in the cases we consider L = 1 for partially anti-symmetric boundary conditions and L = 2 for completely antisymmetric boundary conditions. Also, for simplicity we cite a method that does not use adaptively chosen timesteps. Such adaptive methods are given in [34] and lead to similar scaling where Λ is replaced by the time average of the instantaneous values of Λ.
There are two types of oracles that are required by this corollary. Firstly an oracle is required that outputs the location of the j th non-zero matrix element in a given row, where where L j (i) gives the j th non-zero element in row i. We cost a single query to O 1 as n queries to a single qubit oracle because each query to this oracle yields n bits, and it is more realistic to cost the algorithm by the number of qubits output if the dimension of the Hilbert space is large. The corollary also requires an oracle for matrix elements of H(s) O 2 : |i |k |s |0 → |i |k |s |H(s) i,k .
We also use (for convenience) a new oracle, O f , whose role is to prepare a quantum state encoding the particular value of f A (s) or f B (s) that is needed in a given timestep. In general, if we wish to find the value of f p (s) we use the oracle in the following way: This oracle is crucial to our approach because it allows us to remove the multiple controls used in Figure 1. For example, Our cost analysis of the controlled adiabatic evolution follows by converting the controlled evolution in Figure 1 into the evolution of a single larger Hamiltonian. This larger Hamiltonian can then be simulated by conventional means (such as a Trotter-Suzuki decomposition as per Corollary 1).
Theorem 2. Assume that we wish to simulate a coherently controlled adiabatic evolution that uses the controlled evolutions {T e −i 1 0 H(f1(s))T1 , . . . , T e −i 1 0 H(fp(s))Tp } such that H(s) is a Hamiltonian satisfying H (s) ≤ Γ and each H(f j ) is Λ-2k-smooth and all remaining assumptions of Corollary 1 are held for ∆t = max j T j . The query complexity of performing the simulation within error at most using k th -order time-ordered Trotter-Suzuki formulas and oracles that yield one bit per query obeys where z n is the number of times that n → 2 log 2 n must be iterated before achieving a value that is less than or equal to 6 and N f (s) ≤ log 2 (Γ) .
Proof. To see this, note that the controlled unitary evolutions in Figure 1 produce a time-evolution operator of the following block-diagonal form: where U 1 , . . . , U p are the p controlled adiabatic evolutions. By expanding out the unitaries as time-ordered operator exponentials, we see that the ideal time evolution operator is of the form Consider the Hamiltonian H = j |j j| ⊗ H j . It is easy to see using Taylor expansion and the fact that each of the terms in H commute that Here represents the direct sum operation. Thus we have that Now let H(s) = j |j j| ⊗ H j (s). It then follows from the definition of the ordered-operator exponential and the block-diagonal structure of (55) that It then follows that the controlled evolutions in (53) can be expressed as a simulation of a single time-dependent Hamiltonian by taking H j (s) → H(f j (s))T j /T in (56). For simplicity, let us take T = max j T j .
Next we need to find the properties of the dilated Hamiltonian H(s) that describes the controlled evolutions in the controlled adiabatic evolution. Firstly, assuming each H j is the sum of M Hamiltonians that can be efficiently transformed into d-sparse matrices, it follows that H can be expressed as a similar sum. Similarly, since H(s) = j |j j| ⊗ H(f j (s))T j /T , it follows from the fact that H(s) has a direct sum structure and the assumption that each H(f j (s)) is Λ-2k-smooth that for any non-negative integer q ≤ 2k Hence, for T = max j T j , H(s) is at most Λ-2k-smooth. We then have from (46) that the cost of simulating the effective Hamiltonian H(s) using k th -order Trotter-Suzuki formulas is at most The remaining issue is the calculation of C. In order to compute C we need to first show that we can simulate a query to the Hamiltonian oracles for H using those for H. We specifically require two oracles: one that computes the locations of the i th (potentially) non-zero matrix element in any row of H and another that evaluates that matrix element at a fixed value of s.
The oracle for finding the column index for a specified element in row x of H can be constructed as follows. The oracle O 1 has the property that where y(q) is the column index of the q th element in row x. Then for any j we can construct the corresponding oracle by exploiting the block diagonal structure of H via O 1 (q) |j |x |0 = |j |x |y(q) + 2 n (j − 1) = |j |x |j |y(q) .
The oracle O 1 (q) can therefore be enacted using one query to O 1 and a polynomial size arithmetic circuit. The second oracle O 2 (q) gives for a specific value of s that is specified, the value of H(s). Specifically, after taking into account the block diagonal structure of H, we need the oracle to be of the form This oracle can be implemented using one query to O f and one query to O 2 .
In [34], it is assumed that the time is provided to the oracles via classical control. Here, we assume that the time is provided via a quantum register so we must add the cost of preparing this register to the cost, C, of simulating a one-sparse matrix. Lemma 9 of [34] gives us that, the query complexity (costed at 1/per bit of output yielded by O 1 and O 2 ) is For each one-sparse Hamiltonian that appears in the Trotter-Suzuki decomposition, the time register must be initialized once [34]. This causes an additional source of error and if we are to fit it within our error budget, we must reduce the error in other parts of the simulation algorithm. There are three sources of error in the simulation algorithm: Trotter-Suzuki error, error due to finite n H and error due to finite n T (we have neglected errors in synthesizing single qubit operations etc). Each of these three sources of error is chosen to be at most /3 in [34]. Therefore, if we reduce the error tolerance in the Trotter-Suzuki approximation to /6 and allow an error tolerance of /6 for approximating f j (s) then the overall error will remain at most . Thus the overall complexity becomes The error in e −iH(s)T is at most ∆H(s) T [35], where ∆H(s) is the error in implementing the Hamiltonian. By Taylor's theorem this is at most Γ max j |∆f j (s)| max j T j , where ∆f j (s) is the error incurred by approximating f j (s) to a finite number of digits. Let us define the number of digits used to express f j as n fj . Then the error in f j is ∆f j ≤ 2 −n f j . Hence it suffices to choose Thus since we have to both compute the value of f j (s) to n fj bits of precision using queries to O 2 and then uncompute it, equations (64) and (62) give as claimed.
We therefore see from Theorem 2 that this model of adiabatic computation can be efficiently simulated using the posited oracles under reasonable smoothness assumptions. This naturally leads to the following corollary: O f can be efficiently computed for each j by assumption and hence for any fixed j and s the state |f j (s) can be prepared efficiently. Furthermore, because p ∈ O(poly(n)), it follows that the state p j=1 a j |j |f j (s) can be efficiently prepared. Thus O f can be efficiently simulated in the circuit model as well. This implies that quantum computers can efficiently simulate this class of coherently controlled adiabatic evolutions.
Local Hamiltonians are a subset of d-sparse Hamiltonians. Therefore the class of adiabatic evolutions considered includes a set of Hamiltonians that generate a family of evolutions that are polynomially equivalent to the circuit model [36]. Thus if we ignore the control register then the controlled adiabatic evolution can be reduced to a universal adiabatic quantum computer. Thus our model of computation is polynomially equivalent to both the circuit model and adiabatic quantum computation.
We now see that controlled adiabatic quantum computation using piecewise smooth, sparse, bounded, rowcomputable Hamiltonians is not an exponentially more powerful model of computation than traditional adiabatic computation. Apart from showing that this is a reasonable model of quantum computation, it also shows that the maximum evolution time used is a reasonable metric for the cost of the evolution (once made dimensionless by multiplying by a characteristic energy of the system). For most of the adiabatic paths considered, the contribution of the derivatives of the Hamiltonian to Λ is negligible and thus in practice it suffices to ignore their contributions. Also, because this algorithm scales near-linearly with the evolution time, this analysis clearly shows that our model can only potentially provide sub-polynomial speedups over circuit based quantum computation for fixed d and M .

VIII. CONCLUSION
The big question that our work addresses is: does coherent control over an adiabatic state preparation protocol yield any concrete benefit? We find that it does. Specifically, it allows us to combine the best features of local adiabatic evolution and boundary cancellation methods, which are two optimization strategies that are traditionally at odds with each other. These combined strategies provide better error scaling than any known adiabatic optimization technique for evolution times that are close to the transition between the Landau-Zener regime and the adiabatic regime. We finally provide an explicit quantum simulation algorithm for simulating our protocol, and in turn traditional adiabatic algorithms, that explicitly gives the cost of simulating the controlled adiabatic evolution using a quantum computer and find that this cost scales near-linearly in the evolution time.
These results only begin to scratch the surface of what is possible within this paradigm. Our approach explicitly uses linear combinations of unitary operations that are nearly unitary. This raises an interesting question of whether truly non-unitary processes will be of use in optimizing adiabatic passage. Progress towards this goal has already been reported in [37]. Also, techniques similar to ours may be of value in phase randomization protocols similar to [38]. Ultimately, these ideas may even lead to more natural ways of performing error correction or suppression in a coherently controlled adiabatic quantum computer. These are just a few examples of the many avenues of research that are opened by this work.
Hence, by choosing f A (s) and f B (s) as in Section III A and using (25), we get rid of terms containing derivatives at the end for both levels. This approach trivially generalizes to higher-dimensional systems.
Now we need to fix T A and T B in order to remove the end derivatives parts by requiring that evolutions gain opposite phases. We can rewrite already simplified (A1), (A2) as 1 0 γ A g,e1 (ξ)dξ 1 0 γ B g,e1 (ξ)dξ 1 0 γ A g,e2 (ξ)dξ 1 0 γ B g,e2 (ξ)dξ T A T B = (2n + 1)π (2m + 1)π . (A7) This system of equation has a solution, unless the determinant of the matrix equals zero. Note, that with this approach we get only a discrete set of times T A and T B for which the error vanishes in contrast to many of our prior methods. Error suppression can also be achieved for arbitrary time if we use more than 2 evolutions. A 2-level inspired solution uses 4 unitaries U A , U B , U C and U D where U A and U C are given by Hamiltonian H(f A (s)) and U B and U D by H(f B (s)). We pick the functions f A , f B based on Section III A. The goal is then to find times T A -T D and weights a-d such that following equations are satisfied First, the normalization condtion a + b + c + d = 1 (A10) must hold. Second, we choose b and T B such that they cancel the error on first level from evolution by U A . This is exactly the same problem we solved for a 2 level system, hence the proper b and T B are T B = 1 0 γ A g,e1 (ξ) dξT A + (2n + 1) π 1 0 γ B g,e1 (ξ) dξ The same can be done for U C and U D : This suppresses the first transition out of |g (typically the transition to the first excited state) to O(1/T 2 ). In addition, the errrors from the derivatives at the end on the second excited state cancel as well.
Therefore we can set the value of T C and ratio of a and c and we are still free to choose arbitrary T A . (ξ)dξT C +(2m+1)π) × cos 1 2 1 0 γ A g,e2 (ξ)dξT C + 1 0 γ B g,e2 (ξ)dξ 1 0 γ B g,e1 (ξ)dξ 1 0 γ A g,e1 (ξ) dξT C + (2m + 1) π . (A16) We can ensure that both sides of the equation pick the same phases (up to 2kπ) by setting T C and then we only need a value of k for which the cosines would have the same sign. That is possible, unless This procedure can also be used to find paths that cancel multiple transitions in higher-dimensional systems but a closed form may not necessarily exist, unlike the present case.