Access the full text.
Sign up today, get DeepDyve free for 14 days.
Loading next page...
/lp/association-for-computing-machinery/quantum-linear-system-solver-based-on-time-optimal-adiabatic-quantum-SrKDnzw86h
References
References for this paper are not available at this time. We will be adding them shortly, thank you for your patience.
- Publisher
- Association for Computing Machinery
- Copyright
- Copyright © 2022 Copyright held by the owner/author(s).
- ISSN
- 2643-6809
- eISSN
- 2643-6817
- DOI
- 10.1145/3498331
- Publisher site
- See Article on Publisher Site
Abstract
Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Computing and Quantum Approximate Optimization Algorithm DONG AN, Department of Mathematics, University of California, Berkeley LIN LIN, Department of Mathematics and Challenge Institute of Quantum Computation, University of California, Berkeley, and Computational Research Division, Lawrence Berkeley National Laboratory We demonstrate that with an optimally tuned scheduling function, adiabatic quantum computing (AQC) can readily solve a quantum linear system problem (QLSP) with O(κ poly(log(κ/ϵ))) runtime, where κ is the condition number, and ϵ is the target accuracy. This is near optimal with respect to both κ and ϵ, and is achieved without relying on complicated amplitude amplification procedures that are difficult to implement. Our method is applicable to general non-Hermitian matrices, and the cost as well as the number of qubits can be reduced when restricted to Hermitian matrices, and further to Hermitian positive definite matrices. The success of the time-optimal AQC implies that the quantum approximate optimization algorithm (QAOA) with an optimal control protocol can also achieve the same complexity in terms of the runtime. Numerical results indicate that QAOA can yield the lowest runtime compared to the time-optimal AQC, vanilla AQC, and the recently proposed randomization method. CCS Concepts: • Theory of computation→ Quantum computation theory;• Mathematics of comput- ing→ Numerical analysis; Additional Key Words and Phrases: Quantum linear system problem, adiabatic quantum computing, quantum approximate optimization algorithm ACM Reference format: Dong An and Lin Lin. 2022. Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Com- puting and Quantum Approximate Optimization Algorithm. ACM Trans. Quantum Comput. 3, 2, Article 5 (February 2022), 28 pages. https://doi.org/10.1145/3498331 1 INTRODUCTION Linear system solvers are used ubiquitously in scientific computing. Quantum algorithms for solv- ing large systems of linear equations, also called the quantum linear system problem (QLSP), This work was partially supported by the Department of Energy under Grant No. DE-SC0017867, the Quantum Algorithm Teams Program under Grant No. DE-AC02-05CH11231 (L.L.), by a Google Quantum Research Award, and by the NSF Quantum Leap Challenge Institute (QLCI) program through grant number OMA-2016245 (D. A. and L. L.). Authors’ addresses: D. An, Department of Mathematics, University of California, Berkeley, California 94720; email: dong_an@berkeley.edu; L. Lin, Department of Mathematics and Challenge Institute of Quantum Computation, University of California, Berkeley, California 94720 and Computational Research Division, Lawrence Berkeley National Laboratory, Berkeley, California 94720; email: linlin@math.berkeley.edu. Author current affliation: Dong An’s, Joint Center for Quantum Information and Computer Science, University of Maryland. This work is licensed under a Creative Commons Attribution International 4.0 License. © 2022 Copyright held by the owner/author(s). 2643-6817/2022/02-ART5 $15.00 https://doi.org/10.1145/3498331 ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. 5:2 D. An and L. Lin have received much attention recently [8, 10–12, 16, 17, 30, 33, 34]. The goal of QLSP is to efficiently −1 −1 N×N N compute |x = A |b/A |b on a quantum computer, where A ∈ C ,and |b∈ C is a normalized vector (for simplicity we assume N = 2 ,and A = 1). The ground-breaking Harrow, Hassidim, and Lloyd (HHL) algorithm obtains |x with cost O(poly(n)κ /ϵ),where −1 κ = AA is the condition number of A,and ϵ is the target accuracy. On the other hand, the best classical iterative algorithm is achieved by the conjugate gradient method, where the cost is at least O(N κ log(1/ϵ)), with the additional assumptions that A should be Hermitian positive defi- nite and a matrix-vector product can be done with O(N ) cost [29]. The complexity of direct meth- ods based on the Gaussian elimination procedure removes the dependence on κ, but the depen- dence on N is typically super-linear even for sparse matrices [21]. Therefore, the HHL algorithm can potentially be exponentially faster than classical algorithms with respect to N . The undesirable dependence with respect to ϵ is due to the usage of the quantum phase estimation (QPE)algo- rithm. Recent progresses based on linear combination of unitaries (LCU)[12]and quantum signal processing (QSP)[16, 22] have further improved the scaling toO(κ poly(log(κ/ϵ))) under different query models, without using QPE. However, the O(κ ) scaling can be rather intrinsic to the methods, at least before complex techniques such as variable time amplitude amplification (VTAA) algorithm [2] are applied. The VTAA algorithm is a generalization of the conventional amplitude amplification algorithm, and allows to quadratically amplify the success probability of quantum algorithms in which differ- ent branches stop at different time. In [ 2], VTAA was first used to successfully improve the com- plexity of HHL algorithm to O(κ/ϵ ).In[12], the authors further combine VTAA algorithm and a low-precision phase estimate to improve the complexity of LCU to O(κ poly(log(κ/ϵ))),which is near-optimal with respect to both κ and ϵ. It is worth noting that the VTAA algorithm is a com- plicated procedure and can be difficult to implement. Thus, it remains of great interest to obtain alternative algorithms to solve QLSP with near-optimal complexity scaling without resorting to VTAA. Some of the alternative routes for solving QLSP are provided by the adiabatic quantum computing (AQC)[1, 19] and a closely related method called the randomization method (RM)[6, 30]. The key idea of both AQC and RM is to solve QLSP as an eigenvalue problem with re- spect to a transformed matrix. Assume that a Hamiltonian simulation can be efficiently performed on a quantum computer, it is shown that the runtime of RM scales as O(κ log(κ)/ϵ) [30], which achieves near-optimal complexity with respect to κ without using VTAA algorithm as a subrou- tine. The key idea of the RM is to approximately follow the adiabatic path based on the quantum Zeno effect (QZE) using a Monte Carlo method. Although RM is inspired by AQC, the runtime complexity of the (vanilla) AQC is at least O(κ /ϵ) [1, 7, 30]. Therefore, the RM is found to be at least quadratically faster than AQC with respect to κ. In this article, we find that with a simple modification of the scheduling function to traverse the adiabatic path, the gap between AQC and RM can be fully closed, along with the following two aspects. (1) We propose a family of rescheduled AQC algorithms called AQC(p). Assuming κ (or its upper bound) is known, we demonstrate that for any matrix A (possibly non-Hermitian or dense), when 1 < p < 2, the runtime complexity of AQC(p) can be only O(κ/ϵ). Thus, AQC(p) removes a logarithmic factor with respect to κ compared to RM. (2) We propose another resched- uled algorithm called AQC(exp), of which the runtime is O(κ poly(log(κ/ϵ))). The main benefit of AQC(exp) is the improved dependence with respect to the accuracy ϵ, and this is the near-optimal complexity (up to logarithmic factors) with respect to both κ and ϵ. The scheduling function of AQC(exp) is also universal because we do not even need the knowledge of an upper bound of κ.Ex- isting works along this line [15, 25] only suggest that runtime complexity is O(κ poly(log(κ/ϵ))), which improves the dependence with respect to ϵ at the expense of a much weaker dependence ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Computing 5:3 on κ. Our main technical contribution is to again improve the dependence on κ. Since the cost of any generic QLSP solver can not be less than O(κ) [17], our result achieves the near-optimal complexity up to logarithmic factors. We remark that in the AQC based algorithm, only the total runtime T depends on κ. Beyond the runtime complexity, we also discuss gate-efficient approaches to implement our AQC(p) and AQC(exp) methods. In particular, assume that we are given access to the same query models as those in [12]: the sparse input model of a d-sparse matrix A and the prepare oracle of the state |b. We demonstrate that, when the adiabatic dynamics is simulated using the truncated Dyson series method [23], the query complexity of the AQC(p) method scales O(dκ/ϵ log(dκ/ϵ)), and that of the AQC(exp) method scales O(dκ poly log(dκ/ϵ)). Both algorithms scale almost lin- early in terms of κ, and the AQC(exp) method can achieve near-optimal scaling in both κ and ϵ. Furthermore, the asymptotic scaling of the AQC(exp) method is the same as that of LCU with VTAA method [12, Theorem 5]. However, the AQC(exp) method avoids the usage of complex VTAA rou- tine, which significantly simplifies its practical implementation. The quantum approximate optimization algorithm (QAOA)[14], as a quantum variational algorithm, has received much attention recently thanks to the feasibility of being implemented on near-term quantum devices. Due to the natural connection between AQC and QAOA, our re- sult immediately suggests that the time-complexity for solving QLSP with QAOA is also at most O(κ poly(log(κ/ϵ))), which is also confirmed by numerical results. We also remark that both QAOA and AQC schemes prepare an approximate solution to the QLSP in a pure state, while RM prepares a mixed state. Moreover, all methods above can be efficiently implemented on gate-based comput- ers and are much simpler than those using the VTAA algorithm as a subroutine. 2 QUANTUM LINEAR SYSTEM PROBLEM AND VANILLA AQC N×N Assume A ∈ C is an invertible matrix with condition number κ and A = 1. Let |b∈ C be a normalized vector. Given a target error ϵ, the goal of QLSP is to prepare a normalized state |x ,whichis an ϵ-approximation of the normalized solution of the linear system |x = −1 −1 A |b/A |b , in the sense that|x x |−|xx| ≤ ϵ. 2 a a 2 For simplicity, we first assume A is Hermitian and positive definite and will discuss the general- ization to non-Hermitian case later. The first step to design an AQC-based algorithm for solving QLSP is to transform the QLSP to an equivalent eigenvalue problem. Here, we follow the procedure introduced in [30]. Let Q = I −|bb|. We introduce 0 Q H = σ ⊗ Q = , 0 x b Q 0 then H is a Hermitian matrix and the null space of H is Null(H ) = span{|b,|b}.Here, |b = 0 0 0 |0,b := (b, 0) ,|b = |1,b := (0,b) . The dimension of H is 2N and one ancilla qubit is needed to enlarge the matrix block. We also define 0 AQ H = σ ⊗ (AQ ) + σ ⊗ (Q A) = . 1 + b − b Q A 0 Here, σ = (σ ± iσ ).Notethatif |x satisfies A|x∝ |b,wehave Q A|x = Q |b = 0. Then ± x y b b Null(H ) = span{|x ,|b} with |x = |0, x. Since Q is a projection operator, the gap between 0 1 b and the rest of the eigenvalues of H is 1. The gap between 0 and the rest of the eigenvalues of H 0 1 is bounded from below by 1/κ (see Appendix A). QLSP can be solved if we can prepare the zero-energy state |x of H , which can be achieved by the AQC approach. Let H (f (s)) = (1− f (s))H + f (s)H , 0 ≤ s ≤ 1. The function f :[0, 1] → [0, 1] 0 1 ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. 5:4 D. An and L. Lin is called a scheduling function, and is a strictly increasing mapping with f (0) = 0, f (1) = 1. The simplest choice is f (s) = s, which gives the “vanilla AQC”. We sometimes omit the s-dependence as H (f ) to emphasize the dependence on f . Note that for any s,|b is always in Null(H (f (s))),and there exists a state |x (s) = |0, x (s), such that Null(H (f (s))) = {|x (s),|b}. In particular, |x (0) = |b and |x (1) = |x ; and therefore, |x (s) is the desired adiabatic path. Let P (s) be the projection ¯ ¯ to the subspace Null(H (f (s))), which is a rank-2 projection operator P (s) = |x (s)x (s)|+|bb|. Furthermore, the eigenvalue 0 is separated from the rest of the eigenvalues of H (f (s)) by a gap Δ(f (s)) ≥ Δ (f (s)) := 1− f (s) + f (s)/κ. (1) We refer to the Appendix A for the derivation. Consider the adiabatic evolution i∂ ψ (s) = H (f (s)) ψ (s) , |ψ (0) = |b, (2) s T T T where 0 ≤ s ≤ 1, and the parameter T is called the runtime of AQC. The quantum adiabatic theorem [19, Theorem 3] states that for any 0 ≤ s ≤ 1, |1−ψ (s)|P (s)|ψ (s)| ≤ η (s), (3) T 0 T where (1) 2 (1) (1) s (2) ⎧ ⎫ ⎪ H (s ) ⎪ H (0) H (s) 1 H (s ) 2 2 2 η(s) = C + + + ds . (4) 2 2 2 3 ⎪ ⎪ T Δ (0) T Δ (f (s)) T Δ (f (s )) Δ (f (s )) ⎩ ⎭ (k ) The derivatives of H are taken with respect to s, i.e., H (s) := H (f (s)), k = 1, 2. Throughout ds the article, we shall use a generic symbol C to denote constants independent of s, Δ,T . Intuitively, the quantum adiabatic theorem in Equation (3) says that, if the initial state is an eigen- state corresponding to the eigenvalue 0, then for large enough T the state |ψ (s) will almost stay in the eigenspace of H (s) corresponding to the eigenvalue 0, where there is a double degeneracy and only one of the eigenstate |x (s) is on the desired adiabatic path. However, such degeneracy will not break the effectiveness of AQC for the following reason. Note that b|ψ (0) = 0, and ¯ ¯ H (f (s))|b = 0 for all 0 ≤ s ≤ 1, so the Schrödinger dynamics Equation (2)implies b|ψ (s) = 0, which prevents any transition of |ψ (s) to |b. Therefore, the adiabatic path will stay along |x (s). Using b|ψ (s) = 0, we have P (s) |ψ (s) = |x (s)x (s)|ψ (s). Therefore, the estimate in Equa- T 0 T T tion (3)becomes 2 2 1−|ψ (s)|x (s)| ≤ η (s). This also implies that (see appendix) |ψ (s)ψ (s)|−|x (s)x (s)| ≤ η(s). T T 2 Therefore η(1) can be an upper bound of the distance of the density matrix. (1) (2) If we simply assume H ,H are bounded by constants, and use the worst case bound 2 2 −1 that Δ ≥ κ , we arrive at the conclusion that in order to have η(1) ≤ ϵ, the runtime of vanilla AQC is T κ /ϵ. 3AQC(P)METHOD Our goal is to reduce the runtime by choosing a proper scheduling function. The key observation is that the accuracy of AQC depends not only on the gap Δ(f (s)) but also on the derivatives of H (f (s)), as revealed in the estimate in Equation (4). Therefore, it is possible to improve the accuracy ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Computing 5:5 if a proper time schedule allows the Hamiltonian H (f (s)) to slow down when the gap is close to 0. We consider the following schedule [1, 19]: f (s) = c Δ (f (s)), f (0) = 0, p > 0. (5) 1 −p Here, Δ is defined in Equation ( 1), and c = Δ (u)du is a normalization constant chosen so ∗ p that f (1) = 1. When 1 < p ≤ 2, Equation (5) can be explicitly solved as 1−p p−1 f (s) = 1− 1+ s(κ − 1) . (6) κ − 1 −1 Note that as s → 1, Δ (f (s)) → κ , and therefore the dynamics of f (s) slows down as f → 1and −1 the gap decreases towards κ . We refer to the adiabatic dynamics Equation (2) with the schedule in Equation (5) as the AQC(p) scheme. Our main result is given in Theorem 1 (See Appendix D for the proof). N×N Theorem 1. Let A ∈ C be a Hermitian positive definite matrix with condition number κ.For any choice of 1 < p < 2, the error of the AQC(p) scheme satisfies |ψ (1)ψ (1)|−|x x | ≤ Cκ/T. (7) T T 2 Therefore, in order to prepare an ϵ−approximation of the solution of QLSP, it suffices to choose the run- time T = O(κ/ϵ). Furthermore, when p = 1, 2, the bound for the runtime becomes T = O(κ log(κ)/ϵ). The runtime complexity of the AQC(p) method with respect to κ is only O(κ).Comparedto Reference [30], AQC(p) further removes the log(κ) dependence when 1 < p < 2, and hence reaches the optimal complexity with respect to κ. Interestingly, though not explicitly mentioned in [30], the success of RM for solving QLSP relies on a proper choice of the scheduling function, which approximately corresponds to AQC(p=1). It is this scheduling function, rather than the QZE or its Monte Carlo approximation per se that achieves the desired O(κ logκ) scaling with respect to κ. Furthermore, the scheduling function in RM is similar to the choice of the schedule in the AQC(p=1) scheme. The speedup of AQC(p) versus the vanilla AQC is closely related to the quadratic speedup of the optimal time complexity of AQC for Grover’s search [1, 19, 27, 28], in which the optimal time scheduling reduces the runtime fromT ∼O(N ) (i.e., no speedup compared to classical algorithms) toT ∼O( N ) (i.e., Grover speedup). In fact, the choice of the scheduling function in Reference [28] corresponds to AQC(p=2) and that in Reference [19] corresponds to AQC(1<p<2). 4AQC(EXP)METHOD Although AQC(p) achieves the optimal runtime complexity with respect to κ, the dependence on −1 ϵ is still O(ϵ ), which limits the method from achieving high accuracy. It turns out that when T is sufficiently large, the dependence on ϵ could be improved to O(poly log(1/ϵ)), by choosing an alternative scheduling function. The basic observation is as follows. In the AQC(p) method, the adiabatic error bound we consider, i.e., Equation (4), is the so-called instantaneous adiabatic error bound, which holds true for all s ∈ [0, 1]. However, when using AQC for solving QLSP, it suffices only to focus on the error bound at the final time s = 1. It turns out that this allows us to obtain a tighter error bound. In fact, such an error bound can be exponentially small with respect to the runtime [1, 15, 25, 32]. Roughly speaking, with an additional assumption for the Hamiltonian H (f (s)) that the derivatives of any order vanish at s = 0, 1, the adiabatic error can be bounded by c exp(−c T ) for some positive 1 2 constants c ,c , α. Furthermore, it is proved in [15] that if the target eigenvalue is simple, then 1 2 −1 3 −1 c = O(Δ ) and c = O(Δ ).Notethat Δ ≥ κ for QLSP, thus, according to this bound, to 1 2 ∗ ∗ ∗ obtain an ϵ-approximation, it suffices to choose T = O(κ poly(log(κ/ϵ))).Thisisanexponential speedup with respect to ϵ, but the dependence on the condition number becomes cubic again. ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. 5:6 D. An and L. Lin However, it is possible to reduce the runtime if the change of the Hamiltonian is slow when the gap is small, as we have already seen in the AQC(p) method. For QLSP, the gap monotonically decreases, and the smallest gap occurs uniquely at the final time, where the Hamiltonian H (s) can be set to vary slowly by requiring its derivatives to vanish at the boundary. We consider the following schedule −1 f (s) = c exp − ds , (8) s (1− s ) where c = exp (−1/(s (1− s ))) ds is a normalization constant such that f (1) = 1. This sched- (k ) (k ) ule can assure that H (0) = H (1) = 0 for all k ≥ 1. We refer to the adiabatic dynamics Equation (2) with the schedule in Equation (8) as the AQC(exp) scheme. Our main result is given in Theorem 2 (see Appendix E for the proof). N×N Theorem 2. Let A ∈ C be a Hermitian positive definite matrix with condition number κ.Then for large enough T > 0, the final time error |ψ (1)ψ (1)|−|x x | of the AQC(exp) scheme is T T 2 bounded by κ log κ C log(κ) exp −C . (9) Therefore, for any κ > e, 0 < ϵ < 1, in order to prepare an ϵ−approximation of the solution of QLSP, log κ 2 4 it suffices to choose the runtime T = O(κ log (κ) log ( )). Compared with RM and AQC(p), although the log(κ) dependence reoccurs, AQC(exp) achieves an exponential speedup over RM and AQC(p) with respect to ϵ (and hence giving its name), and thus is more suitable for preparing the solution of QLSP with high fidelity. Furthermore, the time scheduling of AQC(exp) is universal and AQC(exp) does not require knowledge on the bound of κ. We remark that the performance of the AQC(exp) method is sensitive to the perturbations in the scheduling function, which can affect the final error in the AQC(exp) method. This is similar to the finite precision effect in the adiabatic Grover search reported in [ 18]. Therefore, the sched- uling function should be computed to sufficient accuracy on classical computers using numerical quadrature, and implemented accurately as a control protocol on quantum computers. 5 GATE-BASED IMPLEMENTATION OF AQC We briefly discuss how to implement AQC(p) and AQC(exp) on a gate-based quantum computer. Since |ψ (s) = T exp(−iT H (f (s ))ds ) |ψ (0),where T is the time-ordering operator, it is T T sufficient to implement an efficient time-dependent Hamiltonian simulation of H (f (s)). One straightforward approach is to use the Trotter splitting method. The lowest order approxi- mation takes the form T exp −iT H (f (s )) ds ≈ exp (−iThH (f (s ))) m=1 (10) ≈ exp (−iTh(1− f (s ))H ) exp (−iThf (s )H ), m 0 m 1 m=1 where h = s/M, s = mh.Itisprovedin[31] that the error of such an approximation is O(poly (log(N ))T /M), which indicates that to achieve an ϵ-approximation, it suffices to choose M = 2 −iτH −iτH 0 1 O(poly(log(N ))T /ϵ). On a quantum computer, the operations e and e require a time- independent Hamiltonian simulation process, which can be implemented via techniques such as LCU and QSP [4, 22]. For a d-sparse matrix A, according to [5], the query complexity is f = O(д) if f = O(дpoly log(д)) for a single step. Here, the O means that we neglect the log log ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Computing 5:7 Table 1. Computational Costs of AQC(p) and AQC(exp) via a Time-dependent Hamiltonian Simulation using the Truncated Dyson Expansion [23] AQC(p) AQC(exp) Queries O(dκ/ϵ log(dκ/ϵ)) O(dκ poly(log(dκ/ϵ))) Qubits O(n + log(dκ/ϵ)) O(n + log(dκ/ϵ)) Primitive gates O(ndκ/ϵ poly(log(dκ/ϵ))) O(ndκ poly(log(dκ/ϵ))) factors. Note that the total sum of the simulation time of single steps is exactly T regardless of the choice of M, and the total query complexity is O(dT log(dT/ϵ)). Using Theorems 1 and 2, the query complexity of AQC(p) and AQC(exp) is O(dκ/ϵ log(dκ/ϵ)) and O(dκ poly(log(dκ/ϵ))), respectively. Nevertheless, M scales asO(T ) with respect to the runtimeT , which implies that the number of time slices should be at least O(κ ). Therefore the gate complexity scales superlinearly with respect to κ. The scaling of the Trotter expansion can be improved using high order Trotter– Suzuki formula as well as the recently developed commutator-based error analysis [13], but we will not pursue this direction here. There is an efficient way to directly perform the time evolution of H (f (s)) without using the splitting strategy, following the algorithm proposed by Low and Wiebe in [23], where the time- dependent Hamiltonian simulation is performed based on a truncated Dyson expansion. A detailed discussion on how to implement this algorithm in a gate-efficient way is presented in [ 20,Appendix C], and here we summarize the basic idea as follows. Assume that we are given two input query models: P that gives the locations and values of the nonzero entries of the matrix A,and P that A B produces the quantum state |b. Here, the input query models are the same as those in [12]. Then, one can construct the block-encoding representations of the matrix A [16] and the matrix Q with O(1) additional primitive gates. Next, the block-encodings of A and Q can be applied to build the block-encodings of H and H , and then the HAM-T model, which is a block-encoding of the select 0 1 oracle of the time-dependent Hamiltonian H (s) evaluated at different time steps and serves as the input model in the truncated Dyson series method [23]. Finally, after the construction of HAM- T, the adiabatic dynamics can be simulated following the procedure for solving time-dependent Schrödinger equations discussed in [23]. The costs of AQC(p) and AQC(exp) are summarized in Table 1, where for both AQC(p) and AQC(exp), almost linear dependence with respect to κ is achieved. The almost linear dependence 1−δ on κ cannot be expected to be improved to O(κ ) with any δ > 0[17]. Thus, both AQC(p) and AQC(exp) are almost optimal with respect to κ, and AQC(exp) further achieves an exponential speedup with respect to ϵ. 6 QAOA FOR SOLVING QLSP The QAOA scheme [14] considers the following parameterized wavefunction −iγ H −iβ H −iγ H −iβ H P 1 P 0 1 1 1 0 |ψ := e e ··· e e |ψ . (11) θ i Here, θ denotes the set of 2P adjustable real parameters {β ,γ } ,and |ψ is an initial wavefunc- i i i i=1 tion. The goal of QAOA is to choose |ψ and to tune θ,sothat |ψ approximates a target state. In i θ the context of QLSP, we may choose |ψ = |b, and each step of the QAOA ansatz in Equation (11) can be efficiently implemented using the quantum singular value transformation [ 16]. More specif- ically, as discussed in Section 5 and in [20], the block-encodings of H and H can be efficiently 0 1 constructed via the input models for the matrix A and the vector |b. Then, the quantum singular −iβH −iγH 0 1 value transformation can be directly applied to simulate e and e . According to [16, Corol- lary 62], the cost of each single simulation scales linearly in time and logarithmically in precision, ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. 5:8 D. An and L. Lin and hence the total complexity of implementing a QAOA ansatz scales linearly in total runtime of QAOA, defined to be T := (|β | + |γ |), and logarithmically in precision. Notice that with a i i i=1 sufficiently large P, the optimal Trotter splitting method becomes a special form of Equation (11). Hence, Theorem 2 implies that with an optimal choice of{β ,γ } , the QAOA runtimeT is at most i i i=1 O(κ poly(log(κ/ϵ))). We remark that the validity of such an upper bound requires a sufficiently large P and an optimal choice of θ. On the other hand, our numerical results suggest that the same scaling can be achieved with a much smaller P. For a given P,the optimal θ maximizes the fidelity as max F := |ψ |x | . θ θ However, the maximization of the fidelity requires the knowledge of the exact solution |x , which is not practical. We may instead solve the following minimization problem minψ |H |ψ . (12) θ θ Since the null space of H is of dimension 2, the unconstrained minimizer |ψ seems possi- 1 θ ble to only have a small overlap with |x . However, this is not a problem due to the choice of the initial state |ψ = |b. Notice that by the variational principle the minimizer |ψ maxi- i θ −iβH −iγH 0 ¯ 1 ¯ ¯ mizes ψ |P (1)|ψ . Using the fact that e |b = e |b = |b for any β,γ,weobtain θ 0 θ ¯ ¯ ¯ b|ψ = b|b = 0, which means the QAOA ansatz prevents the transition to |b, similar to AQC. Then,ψ |P (1)|ψ = ψ |x x |ψ = F , so the minimizer of Equation (12) indeed maximizes the θ 0 θ θ θ θ fidelity. For every choice of θ, we evaluate the expectation valueψ |H |ψ . Then, the next θ is adjusted θ θ on a classical computer towards minimizing the objective function. The process is repeated till con- vergence. Efficient classical algorithms for the optimization of parameters in QAOA are currently an active topic of research, including methods using gradient optimization [24, 36], Pontryagin’s maximum principle [3, 35], reinforcement learning [9, 26], to name a few. Algorithm 1 describes the procedure using QAOA to solve QLSP. ALGORITHM 1: QAOA for solving QLSP (0) 2P 1: Initial parameters θ = {β ,γ } . i i i=1 2: for k = 0, 1,... do (k ) 3: Perform Hamiltonian simulation to obtain ψ . (k ) (k ) (k ) 2 4: Measure O (θ ) = ψ |H |ψ . θ 1 θ (k ) 2 5: If O (θ ) < ϵ/κ ,exitthe loop. (k+1) 6: Choose θ using a classical optimization method. 7: end for Compared to AQC(p) and AQC(exp), the QAOA method has the following two potential ad- vantages. The first advantage is that QAOA provides the possibility of going beyond AQC-based algorithms. Notice that the Trotter splitting method is a special form of the QAOA ansatz in Equa- tion (11). If the angles{β ,γ } have been properly optimized (which is a very strong assumption i i i=1 and will be further discussed later), the total QAOA runtime T will be by definition comparable to or even shorter than the runtime of AQC with the best scheduling function (after discretization). Second, one way of implementing AQC(p) and AQC(exp) using an operator splitting method re- quires the time interval to be explicitly split into a large number of intervals, while numerical results indicate that the number of intervals P in QAOA can be much smaller. This could reduce ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Computing 5:9 the depth of the quantum circuit. Compared to AQC, QAOA has the additional advantage that it only consists of 2P time-independent Hamiltonian simulation problem, once θ is known. Despite the potential advantages, several severe caveats of using QAOA for QLSP arise when we consider beyond the time complexity. The first is that classical optimization of the angles {β ,γ } i i i=1 can be difficult. Commonly used classical optimization algorithms, such as the gradient descent method, are likely to be stuck at local optimizers and thus result in sub-optimal performance. The cost for the classical optimization is also hard to known apriori. The optimization may require many iterations, which can diminish the gain of the runtime reduction. The second is related to (k ) the accurate computation of the objective function O (θ ). Note that the minimal spectrum gap of −1 2 H isO(κ ). In order to obtain an ϵ-approximation, the precision of measuring O (θ ) = ψ |H |ψ 1 θ θ 2 2 4 4 should be at leastO(ϵ /κ ). HenceO(κ /ϵ ) repeated measurements can be needed to achieve the desired accuracy. 7 GENERALIZATION TO NON-HERMITIAN MATRICES Now we discuss the case when A is not Hermitian positive definite. First, we still assume that A is Hermitian (but not necessarily positive definite). In this case, we adopt the family of Hamiltonians introduced in [30], which overcomes the difficulty brought by the indefiniteness of A at the expense of enlarging the Hilbert space to dimension 4N (so two ancilla qubits are needed to enlarge the matrix block). Here, we define H = σ ⊗ (σ ⊗ I )Q + σ ⊗ Q (σ ⊗ I ) , 0 + z N − z N +,b +,b where Q = I −|+,b+,b|,and |± = (|0± |1). The null space of H is Null(H ) = +,b 2N 0 0 span{|0,−,b,|1,+,b}.Wealsodefine H = σ ⊗ (σ ⊗ A)Q + σ ⊗ Q (σ ⊗ A) . 1 + x +,b − +,b x Note that Null(H ) = span{|0,+, x,|1,+,b}. Therefore, the solution of the QLSP can be obtained if we can prepare the zero-energy state |0,+, x of H . The family of Hamiltonians for AQC(p) is still given by H (f (s)) = (1− f (s))H + f (s)H , 0 ≤ 0 1 s ≤ 1. Similar to the case of Hermitian positive definite matrices, there is a double degeneracy of the eigenvalue 0, and we aim at preparing one of the eigenstate via time-optimal adiabatic evolution. More precisely, for any s, |1,+,b is always in Null(H (f (s))), and there exists a state |x (s) with |x (0) = |0,−,b,|x (1) = |0,+, x, such that Null(H (f (s))) = {|x (s),|1,+,b}. Such degeneracy will not influence the adiabatic computation starting with |0,−,b for the same reason we discussed for Hermitian positive definite case (also discussed in [ 30]), and the error of AQC(p) is still bounded by η(s) given in Equation (4). Furthermore, the eigenvalue 0 is separated from the rest of the eigenvalues of H (f (s)) by a gap 2 2 2 2 Δ(f (s)) ≥ (1− f (s)) + (f (s)/κ) [30]. For technical simplicity, note that (1− f ) + (f/κ) ≥ (1− f + f/κ)/ 2 for all 0 ≤ f ≤ 1, we define the lower bound of the gap to be Δ (f ) = (1− f + f/κ)/ 2, which is exactly proportional to that for the Hermitian positive definite case. Therefore, we can use exactly the same time schedules as the Hermitian positive definite case to perform AQC(p) and AQC(exp) schemes, and properties of AQC(p) and AQC(exp) are stated in the following theorems (see Appendices D and E for the proof). N×N Theorem 3. Let A ∈ C be a Hermitian matrix (not necessarily positive definite) with condition number κ. For any choice of 1 < p < 2, the AQC(p) scheme gives |ψ (s)ψ (s)|−|0,+, x0,+, x| ≤ Cκ/T. (13) T T 2 Therefore, in order to prepare an ϵ−approximation of the solution of QLSP, it suffices to choose the run- time T = O(κ/ϵ). Furthermore, when p = 1, 2, the bound of the runtime becomes T = O(κ log(κ)/ϵ). ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. 5:10 D. An and L. Lin Table 2. Numerical Scaling of the Runtime as a Function of the Condition Number and the Accuracy, Respectively, for the Hermitian Positive Definite Example methods scaling w.r.t. κ scaling w.r.t. 1/ϵ vanilla AQC 2.2022 / RM 1.4912 1.3479 AQC(1) 1.4619 1.0482 AQC(1.25) 1.3289 1.0248 AQC(1.5) 1.2262 1.0008 AQC(1.75) 1.1197 0.9899 AQC(2) 1.1319 0.9904 AQC(exp) 1.3718 0.5377 AQC(exp) / 1.7326 (w.r.t. log(1/ϵ)) QAOA 1.0635 0.4188 QAOA / 1.4927 (w.r.t. log(1/ϵ)) N×N Theorem 4. Let A ∈ C be a Hermitian matrix (not necessarily positive definite) with condition number κ. Then, for large enough T > 0, the final time error |ψ (1)ψ (1)|−|0,+, x0,+, x| T T 2 of the AQC(exp) scheme is bounded by κ log κ C log(κ) exp −C . (14) Therefore, for any κ > e, 0 < ϵ < 1, in order to prepare an ϵ−approximation of the solution of QLSP, log κ 2 4 it suffices to choose the runtime T = O(κ log (κ) log ( )). N×N For a most general square matrix A ∈ C , following [17] we may transform it into the Her- mitian case at the expense of further doubling the dimension of the Hilbert space. Consider an extended QLSP A|x = |b in dimension 2N where 0 A A = σ ⊗ A+ σ ⊗ A = , |b = |1,b. + − A 0 Note that A is a Hermitian matrix of dimension 2N , with condition number κ and A = 1, and |x := |1, x solves the extended QLSP. Therefore, we can directly apply AQC(p) and AQC(exp) for Hermitian matrix A to prepare an ϵ-approximation of x. The total dimension of the Hilbert space becomes 8N for non-Hermitian matrix A (therefore three ancilla qubits are needed). 8 NUMERICAL RESULTS We first report the performance of AQC(p), AQC(exp), and QAOA for a series of Hermitian positive definite dense matrices with varying condition numbers, together with the performance of RM and vanilla AQC. The details of the setup of the numerical experiments are given in Appendix F. Figure 1 shows how the total runtime T depends on the condition number κ and the accuracy ϵ for the Hermitian positive definite case. The numerical scaling is reported in Table 2. For the κ dependence, despite that RM and AQC(1) share the same asymptotic linear complexity with respect to κ, we observe that the preconstant of RM is larger due to its Monte Carlo strategy and the mixed state nature resulting in the same scaling of errors in fidelity and density (see Appendix C for a detailed explanation). The asymptotic scaling of the vanilla AQC is at least O(κ ). When higher ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Computing 5:11 Fig. 1. Simulation results for the Hermitian positive definite example. Top (left/right): the runtime to reach desired fidelity 0.99/0.999 as a function of the condition number. Bottom: a log-log plot of the runtime as a function of the accuracy with κ = 10. fidelity (0.999) is desired, the cost of vanilla AQC becomes too expensive, and we only report the timing of AQC(p), AQC(exp), and QAOA. For the κ dependence tests, the depth of QAOA ranges from 8 to 60. For the ϵ dependence test, the depth of QAOA is fixed to be 20. We find that the runtime for AQC(p), AQC(exp), and QAOA depends approximately linearly on κ, while QAOA has the smallest runtime overall. It is also interesting to observe that although the asymptotic scalings of AQC(1) and AQC(2) are both bounded byO(κ logκ) instead ofO(κ), the numerical performance of AQC(2) is much better than AQC(1); in fact, the scaling is very close to that with the optimal value of p. For the ϵ dependence, the scaling of RM and AQC(p) is O(1/ϵ), which agrees with the error bound. Again the preconstant of RM is slightly larger. Our results also confirm that AQC(exp) only depends poly logarithmically on ϵ. Note that when ϵ is relatively large, AQC(exp) requires a longer runtime than that of AQC(p), and it eventually outperforms AQC(p) when ϵ is small enough. 1.5 The numerical scaling of QAOA with respect to ϵ is found to be only O(log (1/ϵ)) together with the smallest preconstant. Figure 2 and Table 3 demonstrate the simulation results for non-Hermitian matrices. We find that numerical performances of RM, AQC(p), AQC(exp), and QAOA are similar to that of the Hermitian positive definite case. Again QAOA obtains the optimal performance in terms of the runtime. The numerical scaling of the optimal AQC(p) is found to beO(κ/ϵ), while the time complexity of QAOA and AQC(exp) is only O(κ poly(log(κ/ϵ))). ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. 5:12 D. An and L. Lin Fig. 2. Simulation results for the non-Hermitian example. Top: the runtime to reach 0.999 fidelity as a func- tion of the condition number. Bottom: a log-log plot of the runtime as a function of the accuracy with κ = 10. Table 3. Numerical Scaling of the Runtime as a Function of the Condition, Number, and the Accuracy, Respectively, for the Non-Hermitian Example methods scaling w.r.t. κ scaling w.r.t. 1/ϵ vanilla AQC 2.1980 / RM / 1.2259 AQC(1) 1.4937 0.9281 AQC(1.25) 1.3485 0.9274 AQC(1.5) 1.2135 0.9309 AQC(1.75) 1.0790 0.9378 AQC(2) 1.0541 0.9425 AQC(exp) 1.3438 0.4415 AQC(exp) 0.9316 (w.r.t. log(1/ϵ)) QAOA 0.8907 0.3283 QAOA / 0.7410 (w.r.t. log(1/ϵ)) 9 DISCUSSION By reformulating QLSP into an eigenvalue problem, AQC provides an alternative route to solve QLSP other than those based on phase estimation (such as HHL) and those based on approximation of matrix functions (such as LCU and QSP). However, the scaling of the vanilla AQC is at least O(κ /ϵ), which is unfavorable with respect to both κ and ϵ. Thanks to the explicit information of the energy gap along the adiabatic path, we demonstrate that we may reschedule the AQC and dramatically improve the performance of AQC for solving QLSP. When the target accuracy ϵ is relatively large, the runtime complexity of the AQC(p) method (1 < p < 2) is reduced to O(κ/ϵ); for highly accurate calculations with a small ϵ, the AQC(exp) method is more advantageous, and its runtime complexity is O(κ poly(log(κ/ϵ))). To our knowledge, our ACP(exp) method provides the first example that an adiabatic algorithm can simultaneously achieve near linear dependence on the spectral gap, and poly-logarithmic dependence on the precision. Due to the close connection between AQC and QAOA, the runtime complexity of QAOA for solving QLSP is also bounded by O(κ poly(log(κ/ϵ))). Both AQC and QAOA can be implemented on gate-based quantum computers. Our numerical results can be summarized using the following ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Computing 5:13 relation: QAOA AQC(exp) AQC(p) < RM < vanilla AQC. Here, A < B means that the runtime of A is smaller than that of B. The two exceptions are: QAOA AQC(exp) means that the runtime of QAOA is smaller only when the optimizer θ is found, while AQC(exp) AQC(p) holds only when ϵ is sufficiently small. While the runtime complexity of AQC(exp) readily provides an upper bound of the runtime complexity of QAOA, numerical results indicate that the optimizer of QAOA often involves a much smaller depth, and hence the dynamics of QAOA does not necessarily follow the adiabatic path. Therefore, it is of interest to find alternative routes to directly prove the scaling of the QAOA runtime without relying on AQC. In the work [20], our AQC based algorithm has been combined with the eigenvector filtering technique. Reference [ 20] also proposed another AQC inspired quantum linear system solver, which is based on the quantum Zeno effect. Both methods can scale linearly in κ and logarithmically in 1/ϵ. We expect our AQC based QLSP solvers may serve as useful subroutines in other quantum algorithms as well. APPENDICES ATHEGAPOF H (f (s)) FOR HERMITIAN POSITIVE DEFINITE MATRICES The Hamiltonian H (f ) can be written in the block matrix form as 0 ((1− f )I + fA)Q H (f ) = . (15) Q ((1− f )I + fA) 0 Let λ be an eigenvalue of H,then λI −((1− f )I + fA)Q 0 = det −Q ((1− f )I + fA) λI 2 2 = det λ I − ((1− f )I + fA)Q ((1− f )I + fA) , where the second equality holds because the bottom two blocks are commutable. Thus, λ is an 2 2 eigenvalue of ((1−f )I+fA)Q ((1−f )I+fA),and Δ (f ) equals the smallest non-zero eigenvalue of ((1− f )I+fA)Q ((1− f )I+fA). Applying a proposition of matrices that XY andYX have the same 2 2 non-zero eigenvalues, Δ (f ) also equals the smallest non-zero eigenvalue of Q ((1− f )I+ fA) Q . b b Now we focus on the matrix Q ((1− f )I + fA) Q .Notethat |b is the unique eigenstate corre- b b sponding to the eigenvalue 0, and all eigenstates corresponding to non-zero eigenvalues must be orthogonal to with |b. Therefore, 2 2 Δ (f ) = inf φ Q ((1− f )I + fA) Q φ b b b|φ =0, φ|φ =1 = inf φ ((1− f )I + fA) φ b|φ=0,φ|φ=1 ≥ inf φ ((1− f )I + fA) φ φ|φ=1 = (1− f + f/κ) , and Δ(f ) ≥ Δ (f ) = 1− f + f/κ. B RELATIONS AMONG DIFFERENT MEASUREMENTS OF ACCURACY We show two relations that connect the error of density matrix distance and the error of fidelity, which are used in our proof for AQC(p) and AQC(exp). Following the notations in the main text, let |x (s) denote the desired eigenpath of H (f (s)) corresponding to the 0 eigenvalue, and ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. 5:14 D. An and L. Lin Null(H (f (s))) = {|x (s),|b}. P (s) denotes the projector onto the eigenspace corresponding to the 0 eigenvalue. Lemma 5. (i) The following equation holds, 2 2 |1−ψ (s)|P (s)|ψ (s)| = 1− ψ (s)|x (s) = |ψ (s)ψ (s)|−|x (s)x (s)| . (16) T 0 T T T T (ii) Assume that |1−ψ (s)|P (s)|ψ (s)| ≤ η (s). T 0 T Then the fidelity can be bounded from below by 1− η (s), and the 2-norm error of the density matrix can be bounded from above by η(s). Proof. It suffices only to prove part (i). Note that |b is the eigenstate for both H and H cor- 0 1 ¯ ¯ responding the 0 eigenvalue, we have H (f (s))|b = ((1 − f (s))H + f (s)H ) |b = 0, and thus 0 1 ¯ ¯ ¯ b|ψ (s) = 0. Together with the initial condition b|ψ (0) = 0, the overlap of |b and |ψ (s) T T T ds ¯ ¯ ¯ remains to be 0 for the whole time period, i.e., b|ψ (s) = 0. Since P (s) = |x (s)x (s)| + |bb|, T 0 we have P (s)|ψ (s) = |x (s)x (s)|ψ (s)). Therefore, 0 T T |1−ψ (s)|P (s)|ψ (s)| = |1−ψ (s)|x (s)x (s)|ψ (s)| = 1− ψ (s)|x (s) . T 0 T T T T To prove the second equation, let M = |ψ (s)ψ (s)|−|x (s)x (s)|.Notethat M = T T † † λ (M M), we study the eigenvalues of M M by first computing that max M M = |ψ (s)ψ (s)| + |x (s)x (s)|−ψ (s)|x (s)|ψ (s)x (s)|−x (s)|ψ (s)|x (s)ψ (s)|. T T T T T T ⊥ † Since for any |y∈ span{|ψ (s),|x (s)} , M M |y = 0, and M M |ψ (s) = (1− ψ (s)|x (s) ) |ψ (s), T T T M M |x (s) = (1− ψ (s)|x (s) ) |x (s), 2 † we have M = λ (M M) = 1− ψ (s)|x (s) . max T C DIFFERENCE BETWEEN THE SCALINGS OF AQC(P) AND RM WITH RESPECT TO INFIDELITY In our numerical test, we observe that to reach a desired fidelity, RM encounters a much larger pre-constant than AQC(p). This is due to the following reason. Although the runtime of both RM and AQC(p) scales as O(1/ϵ) where ϵ is the 2-norm error of the density matrix, the scalings with respect to the infidelity are different. More specifically, Lemma 5 showsthatforAQC, thesquare of the 2-norm error is exactly equal to the infidelity. Thus, in order to reach infidelity 1 − F using AQC(p), it suffices to choose the runtime to be T = O(κ/ 1− F ). Meanwhile, it has been proved in [6] that the runtime complexity of RM isO(κ/(1− F )). Therefore, to reach a desired fidelity, the runtime of AQC(p) will be smaller than that of RM, as shown in our numerical examples. We further verify the statement above by studying the relation between the 2-norm error of the density matrix and the infidelity for AQC(p), AQC(exp), and RM using the positive definite exam- ple with κ = 10. In AQC(p) and AQC(exp), we change the runtime to obtain approximations with different errors and infidelity. In RM we vary the number of exponential operators to obtain dif- ferent levels of accuracy. All other numerical treatments remain unchanged. As shown in Figure 3, the infidelity is exactly the square of 2-norm error in the case of AQC(p) and AQC(exp), while the infidelity of RM only scales approximately linearly with respect to 2-norm error. This also verifies that although the runtime of both AQC(p) and RM scales linearly with respect to ϵ, the runtime of AQC(p) can be much smaller to reach desired fidelity. ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Computing 5:15 Fig. 3. Relation between 2-norm error and infidelity of AQC and RM. D PROOF OF THEOREMS 1 AND 3 (1) The proof of Theorems 1 and 3 rests on some delicate cancellation of the time derivatives H , (2) H and the gap Δ(f (s)) in the error bound, and can be completed by carefully analyz- ing the κ-dependence of each term in η(s) given in Equation (4). Note that in both cases H (f ) = (1− f )H + fH , and we let Δ (f ) = (1− f + f/κ)/ 2 since such a choice of Δ can serve 0 1 ∗ ∗ as a lower bound of the spectrum gap for both the cases of Theorems 1 and 3. We first compute the derivatives of H (f (s)) by chain rule as d dH (f (s)) df (s) (1) H (s) = H (f (s)) = = (H − H )c Δ (f (s)), 1 0 p ds df ds and d d (2) (1) H (s) = H (s) = (H − H )c Δ (f (s)) 1 0 p ds ds dΔ (f (s)) df (s) p−1 ∗ = (H − H )c pΔ (f (s)) 1 0 p df ds 2p−1 = √ (−1+ 1/κ)(H − H )c pΔ (f (s)). 1 0 p ∗ Then the first two terms of η(s) can be rewritten as (1) (1) (1) (1) H (0) H (s) H (0) H (s) 2 2 2 2 + ≤ + 2 2 2 2 T Δ (0) T Δ (f (s)) T Δ (0) T Δ (f (s)) ∗ ∗ p p (H − H )c Δ (f (0)) (H − H )c Δ (f (s)) 1 0 p ∗ 2 1 0 p ∗ 2 = + 2 2 T Δ (0) T Δ (f (s)) ∗ ∗ p−2 p−2 ≤ c Δ (0) + c Δ (f (s)) p ∗ p ∗ p−2 p−2 ≤ c Δ (0) + c Δ (1) . p ∗ p ∗ ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. 5:16 D. An and L. Lin Here, C stands for a general positive constant independent of s, Δ, and T . To compute the remaining two terms of η(s), we use the following change of variable u = f (s ), du = f (s )ds = c Δ (f (s ))ds , ds and the last two terms of η(s) become ˆ ˆ s s (2) (2) 1 H 1 H 2 2 ds ≤ ds 2 2 T Δ T Δ 0 0 ∗ 2p−1 ˆ 2 s (−1+ 1/κ)(H − H )c pΔ (f (s )) 1 0 ∗ 2 = ds Δ (f (s )) 2p−1 1 2 f (s) (−1+ 1/κ)(H − H )c pΔ (u) 1 0 ∗ 2 1 du 2 p Δ (u) 0 c Δ (u) ∗ p f (s) p−3 ≤ (1− 1/κ)c Δ (u)du p−3 ≤ (1− 1/κ)c Δ (u)du , and similarly ˆ ˆ (1) 2 (1) 2 s s H H 1 1 2 2 ds ≤ ds 3 3 T Δ T 0 0 ˆ p (H − H )c Δ (f (s )) 1 1 0 p = ds Δ (f (s )) f (s) (H − H )c Δ (u) 1 du 1 0 p ∗ T Δ (u) c Δ (u) 0 ∗ p f (s) p−3 ≤ c Δ (u)du p−3 ≤ c Δ (u)du . Summarize all terms above, an upper bound of η(s) is ˆ ˆ 1 1 p−2 p−2 p−3 p−3 η(s) ≤ c Δ (0) + c Δ (1) + (1− 1/κ)c Δ (u)du + c Δ (u)du p p p p ∗ ∗ ∗ ∗ 0 0 ˆ ˆ 1 1 p−3 p−3 −(p−2)/2 2−p = 2 c + c κ + (1− 1/κ)c Δ (u)du + c Δ (u)du . p p p ∗ p ∗ 0 0 Finally, since for 1 < p < 2 p/2 2 κ −p p−1 c = Δ (u)du = (κ − 1), p − 1 κ − 1 and −(p−3)/2 2 κ p−3 2−p Δ (u)du = (κ − 1), 2− p κ − 1 ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Computing 5:17 we have C κ κ p−1 2−p η(s) ≤ (κ − 1) + (κ − κ ) T κ − 1 κ − 1 κ κ p−1 2−p p−1 2−p + (κ − 1)(κ − 1) + (κ − 1)(κ − 1) . κ − 1 κ − 1 The leading term of the bound is O(κ/T ) when 1 < p < 2. Now we consider the limiting case when p = 1, 2. Note that the bound for η(s) can still be written as ˆ ˆ 1 1 p−2 p−2 p−3 p−3 η(s) ≤ c Δ (0) + c Δ (1) + (1− 1/κ)c Δ (u)du + c Δ (u)du p p p p ∗ ∗ ∗ ∗ 0 0 −(p−2)/2 2−p = 2 c + c κ + (1− 1/κ)c c + c c . p p p 3−p p 3−p Straightforward computation shows that −1 c = Δ (u)du = 2 log(κ), κ − 1 and −2 c = Δ (u)du = 2 (κ − 1). κ − 1 Hence, when p = 1, 2, C κ log(κ) −(p−2)/2 2−p η(s) ≤ 2 c + c κ + (1− 1/κ)c c + c c ≤ C . p p 1 2 1 2 T T This completes the proof of Theorems 1 and 3. E PROOF OF THEOREMS 2 AND 4 We provide a rigorous proof of the error bound for the AQC(exp) scheme. We mainly follow the methodology of [25] and a part of technical treatments of [15]. Our main contribution is carefully revealing an explicit constant dependence in the adiabatic theorem, which is the key to obtain the O(κ) scaling. In the AQC(exp) scheme, the Hamiltonian H (s) = (1 − f (s))H + f (s)H with 0 1 H ,H ≤ 1and 0 1 1 1 f (s) = exp − ds . (17) c s (1− s ) The normalization constant c = exp(− )dt ≈ 0.0070. Let U (s) denote the corresponding e T 0 t (1−t ) unitary evolution operator, and P (s) denote the projector onto the eigenspace corresponding to 0. We use Δ (f ) = (1− f + f/κ)/ 2 since this can serve as a lower bound of the spectrum gap for both the cases of Theorems 2 and 4. We first restate the theorems universally with more technical details as following. Theorem 6. Assume the condition number κ > e. Then the final time adiabatic error |1 − ψ (1)|P (1)|ψ (1)| of AQC(exp) can be bounded by η where T 0 T (a) for arbitrary N , κ log κ 2 2 4 η = A D log κ C N , 1 2 where A , D, and C are positive constants which are independent of T , κ and N , 1 2 ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. 5:18 D. An and L. Lin (b) if T is large enough such that 2 2 4π κ log κ −1 16eA D ≤ 1, 3 T then κ log κ 2 2 η = C log κ exp − C , 1 2 where A , D,C , and C are positive constants which are independent of T and κ. 1 1 2 Corollary 7. For any κ > e, 0 < ϵ < 1, to prepare an ϵ-approximation of the solution of QLSP log κ 2 4 using AQC(exp), it is sufficient to choose the runtime T = O(κ log κ log ( )). Proof. We start the proof by considering the projector P (s) onto an invariant space of H,then P (s) satisfies i ∂ P (s) = [H (s), P (s)], P (s) = P (s). (18) We try the ansatz (only formally) −j P (s) = E (s)T . (19) j=0 Substitute it into the Heisenberg equation and match terms with the same orders, we get [H (s), E (s)] = 0, i∂ E (s) = [H (s), E (s)], E (s) = E (s)E (s). (20) 0 s j j+1 j m j−m m=0 It has been proved in [25] that the solution of Equation (20) with initial condition E = P is given 0 0 by −1 −1 E (s) = P (s) = −(2πi) (H (s) − z) dz, (21) 0 0 Γ(s) (1) −1 −1 −1 E (s) = (2π ) (H (s) − z) [E (s), P (s)](H (s) − z) dz + S (s) − 2P (s)S (s)P (s), (22) j 0 j 0 j 0 j−1 Γ(s) where Γ(s) = {z ∈ C : |z| = Δ(s)/2} and j−1 S (s) = E (s)E (s). (23) j m j−m m=1 Furthermore, given E = P , such a solution is unique. 0 0 In general, Equation (19) does not converge, so for arbitrary positive integer N we define a truncated series as −j P (s) = E (s)T . (24) N j j=0 Then N N 1 1 (1) (1) −j −j −(N+1) (1) i P − [H, P ] = i E T − [H, E ]T = iT E . N j N j N T T j=0 j=0 ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Computing 5:19 In Lemma 10,weprove that P (0) = P (0) and P (1) = P (1), then the adiabatic error becomes N 0 N 0 −1 |1−ψ (1)|P (1)|ψ (1)| = |ψ (0)|P (0)|ψ (0)−ψ (0)|U (1) P (1)U (1)|ψ (0)| T 0 T T 0 T T T 0 T T −1 ≤P (1) − U (1) P (0)U (1) 0 T 0 T −1 = P (1) − U (1) P (0)U (1) N T N T −1 = ds U P U . N T ds Straightforward computations show that d d T T −1 −1 −1 −1 −1 −1 (U ) = −U (U )U = −U HU U = − U H, T T T T T T T T ds ds i i d d d d −1 −1 −1 −1 U P U = (U )P U + U (P )U + U P (U ) N T N T N T N T T T T T ds ds ds ds T T T (1) −1 −1 −1 −N −1 = − U HP U + U [H, P ]U + U T E U + U P HU N T N T T N T T T T T i i i (1) −N −1 = T U E U , therefore, (1) (1) −N −1 −N |1−ψ (1)|P (1)|ψ (1)| ≤ T U E U ds ≤ T max E . T 0 T T N N s∈[0,1] In Lemma 15, we prove that (the constant c = 4π /3) [(N + 1)!] (1) E ≤ A A A 1 3 N 2 2 (1+ 1) (N + 1) N 4 A 16 [(N + 1)!] 2 −1 3 2 = D log κ A c D log κ 1 f 4 Δ (N + 1) A N [(N + 1)!] 2 −1 3 2 ≤ D log κ 16A Dc κ log κ 4 (N + 1) 2 −1 3 2 4 ≤ A D log κ 16A Dc κ log κN , 4 2 4N where the last inequality comes from the fact that [(N + 1)!] /(N + 1) ≤ 4N . This completes the proof of part (a). When T is large enough, we now choose ⎢ ⎥ ⎢ 4 ⎥ κ log κ ⎢ ⎥ −1 3 N = ⎢ 16eA Dc ⎥ ≥ 1, 1 f ⎢ ⎥ ⎢ ⎥ ⎣ ⎦ then κ log κ 2 −1 3 4 |1−ψ (1)|P (1)|ψ (1)| ≤ A D log κ 16A Dc N T 0 T 1 1 f κ log κ 2 −1 3 ≤ A D log κ exp − 16eA Dc . This completes the proof of part (b). ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. 5:20 D. An and L. Lin The remaining part is devoted to some preliminary results regarding H, E, and the technical es- timates for the growth of E . It is worth mentioning in advance that in the proof we will encounter many derivatives taken on a contour integral. In fact all such derivatives taken on a contour in- −1 tegral will not involve derivatives on the contour. Specifically, since (H (s) − z) is analytic for any 0 < |z| < Δ(s), for any s ∈ (0, 1), there exists a small enough neighborhood B (s ) such that 0 δ 0 ¸ ¸ −1 −1 ∀s ∈ B (s ), G(s, (H (s)−z) )dz = G(s, (H (s)−z) )dz for any smooth mapping G.This δ 0 Γ(s) Γ(s ) means locally the contour integral does not depend on the smooth change of the contour, and thus the derivatives will not involve derivatives on the contour. In the spirit of this trick, we write the −1 (k ) resolvent R(z, s, s ) = (H (s) − z) for 0 ≤ s ≤ 1, 0 ≤ s ≤ 1, z ∈ C and |z| = Δ(s )/2 and let R 0 0 0 (k ) denote the partial derivative with respect to s, i.e., R(z, s, s ), which means by writing R we ∂s only consider the explicit time derivatives brought by H. ∞ (k ) (k ) Lemma 8. (a) H (s) ∈ C with H (0) = H (1) = 0 for all k ≥ 1. (b) There is a gap Δ(s) ≥ Δ (s) = ((1− f (s))+ f (s)/κ)/ 2, which separates 0 from the rest of the spectrum. The following lemma gives the bound for the derivatives of H. Lemma 9. For every k ≥ 1, 0 < s < 1, (k!) (k ) k H (s)≤ b(s)a(s) , (25) (k + 1) where 2e 1 2 b(s) = exp − [s(1− s)] , a(s) = . c s(1− s) s(1− s) Proof. We first compute the derivatives of f .Let д(s) = −s(1 − s) and h(y) = exp(1/y),then −1 f (s) = c h(д(s)). By the chain rule of high order derivatives (also known as Faà di Bruno’s formula), k! j (k+1) −1 (m +m +···+m ) (j) 1 2 f (s) = c h (д(s)) д (s) , m m m 1 2 m !1! m !2! ···m !k! 1 2 k j=1 where the sum is taken over all k-tuples of non-negative integers (m ,...,m ) satisfying 1 k k (j) jm = k.Notethat д (s) = 0for j ≥ 3, and the sum becomes j=1 k! m m 1 2 (k+1) −1 (m +m ) (1) (2) 1 2 f (s) = c h (д(s)) д (s) д (s) m m 1 2 m !1! m !2! 1 2 m +2m =k 1 2 k! −1 (m +m ) m m 1 2 1 2 = c h (д(s)) (2s − 1) 2 m !m !2 1 2 m +2m =k 1 2 k! −1 (m +m ) m 1 2 1 = c h (д(s)) (2s − 1) . m !m ! 1 2 m +2m =k 1 2 ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Computing 5:21 To compute the derivatives of h, we use the chain rule again to get (the sum is over jn = m) j=1 m n d (1/y) m! (m) h (y) = exp(1/y) n n n j 1 2 m n !1! n !2! ··· n !m! dy 1 2 m j=1 m! j j −j−1 = exp(1/y) (−1) j!y n n n 1 2 m n !1! n !2! ··· n !m! 1 2 m j=1 (−1) m! −m− n = exp(1/y)y . n !n !··· n ! 1 2 m Since 0 ≤ n ≤ m/j, the number of tuples (m ,...,m ) is less than (m + 1)(m/2 + 1)(m/3 + j 1 n 2m 2m 1)... (m/m + 1) = < 2 ,sofor 0 < y < 1and m ≤ k we have (m) 2k −2k |h (y)|≤ 2 k!exp(1/y)y . (k+1) Therefore f can be bounded as 2k k! 1 1 (k+1) −1 2k m |f (s)|≤ c 2 k!exp − |2s − 1| m !m ! s(1− s) s(1− s) 1 2 m +2m =k 1 2 2k 1 2 1 −1 2 ≤ c exp − (k!) s(1− s) s(1− s) m ! m ≤k 2k 1 2 −1 2 ≤ ec exp − (k!) . s(1− s) s(1− s) Substitute k + 1by k and for every k ≥ 1 2(k−1) 1 2 (k ) −1 2 |f (s)|≤ ec exp − ((k − 1)!) s(1− s) s(1− s) 2(k−1) 1 2 (k!) −1 ≤ 4ec exp − . s(1− s) s(1− s) (k + 1) (k ) (k ) Noting that H ≤ 1,H ≤ 1and H = (H − H ) f , we complete the proof of bounds for 0 1 1 0 (k ) H . The following result demonstrates that E ’s for all j ≥ 1 vanish on the boundary. (k ) (k ) (k ) (k ) Lemma 10. (a) For all k ≥ 1, E (0) = P (0) = 0, E (1) = P (1) = 0. 0 0 0 0 (k ) (k ) (b) For all j ≥ 1, k ≥ 0, E (0) = E (1) = 0. j j (k ) (k ) Proof. We will repeatedly use the fact that R (0) = R (1) = 0. This can be proved by taking the k-th order derivative of the equation (H − z)R = I and k k k k (k ) (l ) (k−l ) (l ) (k−l ) R = −R (H − z) R = −R H R . l l l=1 l=1 (k ) (a) This is a straightforward result by the definition of E and the fact that R ’s vanish on the boundary. (b) We prove by induction with respect to j. For j = 1, Equation (22) tells that (1) −1 E = (2π ) R[P , P ]Rdz. 1 0 ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. 5:22 D. An and L. Lin Therefore, each term in the derivatives of E must involve at least one of the derivative of R or the derivative of P , which means the derivatives of E much vanish on the boundary. 0 1 Assume the conclusion holds for < j,thenfor j, first each term of the derivatives of S must involve the derivative of some E with m < j, which means the derivatives of S must vanish on m j the boundary. Furthermore, for the similar reason, Equation (22) tells that the derivatives of E must vanish on the boundary. Before we process, we recall three technical lemmas introduced in [15, 25]. Throughout let c = 4π /3 denote an absolute constant. Lemma 11. Let α > 0 be a positive real number, p,q be non-negative integers, and r = p+ q.Then, 1+α 1+α k [(l + p)!(k − l + q)!] [(k + r )!] ≤ c . 2 2 2 l (l + p + 1) (k − l + q + 1) (k + r + 1) l=0 Lemma 12. Let k be a non-negative integer, then 1 1 ≤ c . 2 2 2 (l + 1) (k + 1− l) (k + 1) l=0 Lemma 13. Let A(s), B(s) be two smooth matrix-valued functions defined on [0, 1] satisfying 1+α 1+α [(k + p)!] [(k + q)!] (k ) k (k ) k A (s)≤ a (s)a (s) , B (s)≤ b (s)b (s) , 1 2 1 2 2 2 (k + 1) (k + 1) for some non-negative functions a , a ,b ,b , non-negative integers p,q, and for all k ≥ 0. Then, for 1 2 1 2 every k ≥ 0, 0 ≤ s ≤ 1, 1+α [(k + r )!] (k ) k (A(s)B(s)) ≤ c a (s)b (s) max{a (s),b (s)} , f 1 1 2 2 (k + 1) where r = p + q. Next we bound the derivatives of the resolvent. This bound provides the most important im- provement of the general adiabatic bound. Lemma 14. For all k ≥ 0, 2 (k!) (k ) 2 R (z, s , s )≤ D log κ , 0 0 Δ(s ) (k + 1) where 2048 2e D = c . Proof. We prove by induction, and for simplicity we will omit explicit dependence on argu- ments z, s, and s . The estimate obviously holds for k = 0. Assume the estimate holds for< k.Take the k-th order derivative of the equation (H − z)R = I and we get k k k k (k ) (l ) (k−l ) (l ) (k−l ) R = −R (H − z) R = −R H R . l l l=1 l=1 Using Lemma 9 and the induction hypothesis, we have 2 4 2 k (l!) 2 k−l [(k − l)!] (k ) l 2 R ≤ ba D log κ . 2 2 Δ l (l + 1) Δ (k − l + 1) l=1 ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Computing 5:23 −1 l To proceed we need to bound the term Δ ba for l ≥ 1. Let us define exp(− ) s (1−s) −1 l F (s) = Δ (s)b(s)a(s) = . 2l−2 2l (1− f (s) + f (s)/κ)[s(1− s)] 2 2 2e Note that F (0) = F (1) = 0, F (s) > 0for s ∈ (0, 1) and F (1/2+ t) > F (1/2− t) for t ∈ (0, 1/2),then there exists a maximizer s ∈ [1/2, 1) such that F (s) ≤ F (s ),∀s ∈ [0, 1]. Furthermore, F (s ) = 0. ∗ ∗ ∗ Now we compute the F as 2l−2 2 [(1− f + f/κ)[s(1− s)] ] F (s) 1 1− 2s 2l−2 = exp − (1− f + f/κ)[s(1− s)] 2 2 s(1− s) s (1− s) 2l−2 2l−3 − exp − (−f + f /κ)[s(1− s)] + (1− f + f/κ)(2l − 2)[s(1− s)] (1− 2s) s(1− s) 2l−4 = exp − [s(1− s)] s(1− s) −1 2 2 × (1− f + f/κ)(1− 2s)[1− (2l − 2)s(1− s)]− exp − c (−1+ 1/κ)s (1− s) s(1− s) 2l−4 = exp − [s(1− s)] G(s), s(1− s) where −1 2 2 G(s) = (1− f + f/κ)(1− 2s)[1− (2l − 2)s(1− s)]+ exp − c (1− 1/κ)s (1− s) . s(1− s) The sign of F (s) for s ∈ (0, 1) is thesameasthe sign of G(s). We now show that s cannot be very close to 1. Precisely, we will prove that for all s ∈ [1 − , 1) with c = c /4 ≈ 0.021, G(s) < 0. For such s,wehave l log κ 1− f + f/κ ≥ f (1/2)/κ > 0, 1− 2s < −1/2, and 2c 1− (2l − 2)s(1− s) ≥ 1− (2l − 2)(1− s) ≥ 1− ≥ 1/2, logκ then f (1/2) 1 (1− f + f/κ)(1− 2s)[1− (2l − 2)s(1− s)]≤− = − . 4κ 8κ On the other hand, −1 1 c l logκ exp − ≤ exp − 1− s(1− s) l logκ c c −1 l −(1− ) l logκ c = κ −l /c ≤ κ −1 ≤ κ , ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. 5:24 D. An and L. Lin then −1 2 2 exp − c (1− 1/κ)s (1− s) s(1− s) 1 1 c κ c l logκ ≤ . 16κ c c Therefore, for all s ∈ [1− , 1] we have G(s) ≤−1/(16κ) < 0, which indicates s ≤ 1− . l log κ l log κ We are now ready to bound F (s). From the equation G(s ) = 0, we get exp − s (1−s ) (1− 2s )[1− (2l − 2)s (1− s )] ∗ ∗ ∗ ∗ ∗ = , −1 2 1− f + f/κ c (−1+ 1/κ)s (1− s ) e ∗ which gives F (s) ≤ F (s ) (1− 2s )[1− (2l − 2)s (1− s )] ∗ ∗ ∗ −1 2l c (−1+ 1/κ)[s (1− s )] e ∗ ∗ 2s − 1 −1 2l c (1− 1/κ)[s (1− s )] e ∗ ∗ 2l −2l ≤ 2c · 2 (1− s ) e ∗ 2l l logκ 2l ≤ 2c · 2 2l 2l = 2c (logκ) l 2c 64e 2l 2 ≤ (logκ) (l!) . e c l l−1 The last inequality comes from the fact l ≤ e l!, which can be derived from the fact that log i ≥ log x dx = n log n − (n − 1). i=1 By definition of F (s) we immediately get √ √ 2 2e 4 2 256e −1 l l 2l 2 Δ ba ≤ 4 F ≤ (logκ) (l!) . c e c e e (k ) Now we go back to the estimate of R . By Lemma 11, 2 4 2 k (l!) 2 k−l [(k − l)!] (k ) l 2 R ≤ ba D log κ 2 2 Δ l (l + 1) Δ (k − l + 1) l=1 k l 2 2 4 2 k 8 2 256e (l!) k−l [(k − l)!] 2l 2 2 ≤ (logκ) (l!) D log κ 2 2 Δ l e c (l + 1) (k − l + 1) l=1 ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Computing 5:25 4 4 2 k (l!) [(k − l)!] 2 k −1 ≤ (D log κ) c f 2 2 Δ l (l + 1) (k − l + 1) l=1 2 (k!) 2 k ≤ (D log κ) . Δ (k + 1) This completes the proof. The next lemma is the main technical result, which gives the bound of derivatives of E defined in Equation (20). Lemma 15. (a) For all k ≥ 0, (k!) (k ) (k ) 2 k |E | = |P |≤ (D log κ) . (26) 0 0 2 (k + 1) (b) For all k ≥ 0, j ≥ 1, [(k + j)!] (k ) j k E ≤ A A A , (27) j 2 2 2 (k + 1) (j + 1) with 1 −1 2 2 A = c 1+ 2c , f f −1 3 2 A = A c D log κ, A = D log κ. Remark 16. Thechoiceof A , A can be rewritten as 1 2 3 2 c D log κ = A A , 1 2 2 2 c 1+ 2c A = . f f Furthermore, using c > 1, we have 16 A 1 c = A ≤ . Δ A 2 These relations will be used in the proof later. Proof. (a) By Lemma 14, (k!) (k ) −1 (k ) 2 k |P (s )| = (2πi) R (z, s , s )dz≤ (D log κ) 0 0 0 0 2 (k + 1) Γ(s ) (b) We prove by induction with respect to j. For j = 1, Equation (22) tells k k d Δ d (k ) (1) (1) −1 E = (2π ) (R[P , P ]R)dz≤ (R[P , P ]R). 0 0 1 0 0 k k ds ds By Lemmas 13 and 14, 2 4 2 [(k + 1)!] (k ) 3 2 2 k E ≤ Δc D log κ(D log κ) 1 f Δ (k + 1) [(k + 1)!] ≤ A A A . 1 2 2 2 (k + 1) (1+ 1) ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. 5:26 D. An and L. Lin Now assume < j the estimate holds, for j, by Lemmas 12, 13, and the induction hypothesis, j−1 [(k + j)!] (k ) j−m m k S ≤ c A A A A A f 1 1 2 3 j 2 2 2 2 (k + 1) (m + 1) (j − m + 1) m=1 j−1 [(k + j)!] 1 2 k = A A A c 1 2 3 2 2 2 (k + 1) (m + 1) (j − m + 1) m=1 [(k + j)!] 2 2 k ≤ c A A A . 1 3 f 2 2 2 (k + 1) (j + 1) Again by Lemmas 13, 14, and the induction hypothesis, k k k d d d (k ) −1 (1) E ≤ (2π ) R[E , P ]Rdz + S + 2P S P 0 j 0 j 0 j j−1 k k k ds ds ds 4 4 2 1 [(k + j)!] [(k + j)!] j−1 j 3 k 2 2 k ≤ Δc A A A A + c A A A 1 3 f 2 3 f 2 3 2 2 2 2 Δ j (k + 1) (k + 1) (j + 1) 1 [(k + j)!] 2 2 2 j k + 2c c A A A 1 3 f f 2 2 2 (j + 1) (k + 1) 4 4 16 [(k + j)!] [(k + j)!] j−1 j 3 k+1 2 2 k ≤ c A A A + c A A A f 2 3 f 1 2 3 2 2 2 2 Δ (k + 1) (j + 1) (k + 1) (j + 1) [(k + j)!] 4 2 k + 2c A A A f 2 3 2 2 (k + 1) (j + 1) 16 A [(k + j)!] 3 2 2 j k = c + c 1+ 2c A × A A A 1 1 f f f 2 2 2 Δ A (k + 1) (j + 1) [(k + j)!] ≤ A A A . 2 3 2 2 (k + 1) (j + 1) F DETAILS OF THE NUMERICAL TREATMENTS AND EXAMPLES For simulation purpose, the AQC schemes are carried out using and the induction hypothesis with a time step size 0.2. We use the gradient descent method to optimize QAOA and record the running time corresponding to the lowest error in each case. In QAOA, we also use the true fidelity to measure the error. RM is a Monte Carlo method, and each RM calculation involves performing (i) 200 independent runs to obtain the density matrix ρ for i-th repetition, then we use the averaged (i) density ρ¯ = 1/n ρ to compute the error. We report the averaged runtime of each single RM rep calculation. We perform calculations for a series of 64-dimensional Hermitian positive definite dense matrices A , and 32-dimensional non-Hermitian dense matrices A with varying condition 1 2 number κ. For concreteness, for the Hermitian positive definite example, we choose A = U ΛU .Here, U is an orthogonal matrix obtained by Gram–Schmidt orthogonalization (implemented via a QR factorization) of the discretized periodic Laplacian operator given by 1 −0.5 −0.5 −0.51 −0.5 −0.51 −0.5 L = . (28) . . . . . . . . . −0.51 −0.5 −0.5 −0.51 ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. Quantum Linear System Solver Based on Time-optimal Adiabatic Quantum Computing 5:27 Λ is chosen to be a diagonal matrix with diagonals uniformly distributed in [1/κ, 1]. More precisely, Λ = diag(λ , λ ,..., λ ) with λ = 1/κ+ (k−1)h,h = (1−1/κ)/(N−1). Such construction ensures 1 2 N k A to be a Hermitian positive definite matrix which satisfies A = 1 and the condition number N N of A is κ.Wechoose |b = u / u where {u } is the set of the column vectors of U . k k 2 k k=1 k=1 Here, N = 64. For the non-Hermitian positive definite example, we choose A = U ΛV .Here, U is the same as those in the Hermitian positive definite case, except that the dimension is reduced to N = 32. Λ = diag(λ , λ ,..., λ ) with λ = (−1) (1/κ+ (k−1)h),h = (1−1/κ)/(N−1). V is an orthogonal 1 2 N k matrix obtained by Gram–Schmidt orthogonalization of the matrix 2 −0.5 −0.5 −0.52 −0.5 −0.52 −0.5 K = . (29) . . . . . . . . . −0.52 −0.5 −0.5 −0.52 Such construction ensures A to be non-Hermitian, satisfying A = 1 and the condition number of A is κ. We choose the same |b as that in the Hermitian positive definite example. ACKNOWLEDGMENTS We thank Rolando Somma, Yu Tong and Nathan Wiebe for helpful discussions. REFERENCES [1] Tameem Albash and Daniel A. Lidar. 2018. Adiabatic quantum computation. Rev. Mod. Phys. 90, 1 (2018), 015002. [2] Andris Ambainis. 2012. Variable time amplitude amplification and quantum algorithms for linear algebra problems. In Proceedings of the STACS’12 (29th Symposium on Theoretical Aspects of Computer Science). Vol. 14. LIPIcs, Paris, France, 636–647. [3] Seraph Bao, Silken Kleer, Ruoyu Wang, and Armin Rahmani. 2018. Optimal control of superconducting gmon qubits using pontryagin’s minimum principle: Preparing a maximally entangled state with singular bang-bang protocols. Phys. Rev. A 97, 6 (2018), 062343. [4] Dominic W. Berry, Andrew M. Childs, Richard Cleve, Robin Kothari, and Rolando D. Somma. 2015. Simulating hamil- tonian dynamics with a truncated taylor series. Phys. Rev. Lett. 114, 9 (2015), 090502. [5] Dominic W. Berry, Andrew M. Childs, and Robin Kothari. 2015. Hamiltonian simulation with nearly optimal depen- dence on all parameters. In Proceedings of the 2015 IEEE 56th Annual Symposium on Foundations of Computer Science. IEEE, Piscataway, NJ, 792–809. [6] Sergio Boixo, Emanuel Knill, and Rolando D. Somma. 2009. Eigenpath traversal by phase randomization. Quantum Info. Comput. 9 (2009), 833–855. [7] Sergio Boixo and Rolando D. Somma. 2010. Necessary condition for the quantum adiabatic approximation. Phys. Rev. A 81, 3 (2010), 032308. [8] Carlos Bravo-Prieto, Ryan LaRose, M. Cerezo, Yigit Subasi, Lukasz Cincio, and Patrick J. Coles. 2020. Variational Quantum Linear Solver. arXiv:1909.05820. Retrieved from https://arxiv.org/abs/1909.05820. [9] Marin Bukov, Alexandre G. R. Day, Dries Sels, Phillip Weinberg, Anatoli Polkovnikov, and Pankaj Mehta. 2018. Re- inforcement learning in different phases of quantum control. Phys. Rev. X 8, 3 (2018), 031086. [10] Yudong Cao, Anargyros Papageorgiou, Iasonas Petras, Joseph Traub, and Sabre Kais. 2013. Quantum algorithm and circuit design solving the poisson equation. New J. Phys. 15, 1 (2013), 013021. [11] Shantanav Chakraborty, András Gilyén, and Stacey Jeffery. 2019. The power of block-encoded matrix powers: Im- proved regression techniques via faster hamiltonian simulation. In Proceedings of the 46th International Colloquium on Automata, Languages, and Programming (ICALP 2019) (Leibniz International Proceedings in Informatics (LIPIcs), Vol. 132). Schloss Dagstuhl–Leibniz-Zentrum fuer Informatik, Dagstuhl, Germany, 33:1–33:14. [12] Andrew M. Childs, Robin Kothari, and Rolando D. Somma. 2017. Quantum algorithm for systems of linear equations with exponentially improved dependence on precision. SIAM J. Comput. 46, 6 (2017), 1920–1950. ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022. 5:28 D. An and L. Lin [13] Andrew M. Childs, Yuan Su, Minh C. Tran, Nathan Wiebe, and Shuchen Zhu. 2021. Theory of trotter error with commutator scaling. Physical Review X 11, 1 (2021), 011020. [14] Edward Farhi, Jeffrey Goldstone, and Sam Gutmann. 2014. A Quantum Approximate Optimization Algorithm. arXiv:1411.4028. Retrieved from https://arxiv.org/abs/1411.4028. [15] Yimin Ge, András Molnár, and J. Ignacio Cirac. 2016. Rapid adiabatic preparation of injective projected entangled pair states and gibbs states. Phys. Rev. Lett. 116 (2016), 080503. [16] András Gilyén, Yuan Su, Guang Hao Low, and Nathan Wiebe. 2019. Quantum singular value transformation and beyond: Exponential improvements for quantum matrix arithmetics. In Proceedings of the 51st Annual ACM SIGACT Symposium on Theory of Computing (STOC 2019). Association for Computing Machinery, New York, NY, 193–204. [17] Aram W. Harrow, Avinatan Hassidim, and Seth Lloyd. 2009. Quantum algorithm for linear systems of equations. Phys.Rev.Lett. 103, 15 (2009), 150502. [18] Itay Hen. 2019. How quantum is the speedup in adiabatic unstructured search? Quant. Inf. Proc. 18, 6 (2019), 162. [19] Sabine Jansen, Mary-Beth Ruskai, and Ruedi Seiler. 2007. Bounds for the adiabatic approximation with applications to quantum computation. J. Math. Phys. 48, 10 (2007), 102111. [20] Lin Lin and Yu Tong. 2020. Optimal polynomial based quantum eigenstate filtering with application to solving quan- tum linear systems. Quantum 4 (2020), 361. [21] Joseph W. H. Liu. 1992. The multifrontal method for sparse matrix solution: Theory and practice. SIAM Rev. 34, 1 (1992), 82–109. [22] Guang Hao Low and Isaac L. Chuang. 2017. Optimal hamiltonian simulation by quantum signal processing. Phys. Rev. Lett. 118, 1 (2017), 010501. [23] Guang Hao Low and Nathan Wiebe. 2019. Hamiltonian Simulation in the Interaction Picture. arXiv:1805.00675. Re- trieved from https://arxiv.org/abs/1805.00675. [24] Yvon Maday and Gabriel Turinici. 2003. New formulations of monotonically convergent quantum control algorithms. J. Chem. Phys. 118, 18 (2003), 8191–8196. [25] Gheorghe Nenciu. 1993. Linear adiabatic theory exponential estimates. Comm. Math. Phys. 152, 3 (1993), 479–496. [26] Murphy Yuezhen Niu, Sergio Boixo, Vadim N. Smelyanskiy, and Hartmut Neven. 2019. Universal quantum control through deep reinforcement learning. npj Quantum Info. 5, 1 (2019), 33. [27] Ali T. Rezakhani, W.-J. Kuo, Alioscia Hamma, Daniel A. Lidar, and Paolo Zanardi. 2009. Quantum adiabatic brachis- tochrone. Phys.Rev.Lett. 103 (2009), 080502. [28] Jérémie Roland and Nicolas J. Cerf. 2002. Quantum search by local adiabatic evolution. Phys. Rev. A 65, 4 (2002), [29] Yousef Saad. 2003. Iterative Methods for Sparse Linear Systems. Vol. 82. SIAM, Philadelphia, PA. [30] Yiğit Subaşı, Rolando D. Somma, and Davide Orsucci. 2019. Quantum algorithms for systems of linear equations inspired by adiabatic quantum computing. Phys.Rev.Lett. 122, 6 (2019), 060504. [31] Wim van Dam, Michele Mosca, and Umesh Vazirani. 2001. How powerful is adiabatic quantum computation? In Proceedings 42nd IEEE Symposium on Foundations of Computer Science. IEEE, Piscataway, NJ, 279–287. [32] Nathan Wiebe and Nathan S. Babcock. 2012. Improved error-scaling for adiabatic quantum evolutions. New J. Phys. 14, 1 (2012), 1–10. [33] Leonard Wossnig, Zhikuan Zhao, and Anupam Prakash. 2018. Quantum linear system algorithm for dense matrices. Phys.Rev.Lett. 120, 5 (2018), 050502. [34] Xiaosi Xu, Jinzhao Sun, Suguru Endo, Ying Li, Simon C. Benjamin, and Xiao Yuan. 2021. Variational algorithms for linear algebra. Science Bulletin in press (2021). [35] Zhi-Cheng Yang, Armin Rahmani, Alireza Shabani, Hartmut Neven, and Claudio Chamon. 2017. Optimizing varia- tional quantum algorithms using pontryagin’s minimum principle. Phys.Rev.X 7, 2 (2017), 021027. [36] Wusheng Zhu and Herschel Rabitz. 1998. A rapid monotonically convergent iteration algorithm for quantum optimal control over the expectation value of a positive definite operator. J. Chem. Phys. 109, 2 (1998), 385–391. Received April 2020; revised September 2021; accepted November 2021 ACM Transactions on Quantum Computing, Vol. 3, No. 2, Article 5. Publication date: February 2022.
Journal
ACM Transactions on Quantum Computing – Association for Computing Machinery
Published: Mar 4, 2022
Keywords: Quantum linear system problem