Abstract We introduce a family of identities that express general linear non-unitary evolution operators as a linear combination of unitary evolution operators, each solving a Hamiltonian simulation problem. This formulation can exponentially enhance the accuracy of the recently introduced linear combination of Hamiltonian simulation (LCHS) method [An, Liu, and Lin, Physical Review Letters, 2023]. For the first time, this approach enables quantum algorithms to solve linear differential equations with both optimal state preparation cost and near-optimal scaling in matrix queries on all parameters.
Abstract We introduce a family of identities that express general linear non-unitary evolution operators as a linear combination of unitary evolution operators, each solving a Hamiltonian simulation problem. This formulation can exponentially enhance the accuracy of the recently introduced linear combination of Hamiltonian simulation (LCHS) method [An, Liu, and Lin, Physical Review Letters, 2023]. For the first time, this approach enables quantum algorithms to solve linear differential equations with both optimal state preparation cost and near-optimal scaling in matrix queries on all parameters.
Quantum linear systems algorithms (QLSAs) aim to solve linear systems $A\bb{x}=\bb{b}$ exponentially faster than classical methods under certain conditions. In this work, we develop quantum algorithms for solving linear algebraic equations from an ODE-based perspective. Inspired by the linear combination of Hamiltonian simulation (LCHS) representation in the Fourier approach \cite{Childs2017QLSA}, we express the solution $\bb{x}$ as a linear combination of solutions to a system of linear convection equations, which become Schrödinger-type equations with unitary evolutions in the Fourier domain. We refer to this representation as LC-Schrödingerization. Based on this result, we construct an LCHS-based quantum algorithm with two LCHS instances: one for time-marching and one for numerical integration. The key construction uses the derivative of a Gaussian-smoothed hat function and recovers the solution over a fixed auxiliary interval. This permits a truncation time independent of the target accuracy and avoids the loss in success probability from selecting a single grid point. Periodization and explicit Fourier coefficients provide the corresponding projection error bounds. Under the stated oracle assumptions and given a constant-factor estimate of the solution norm, direct simulation of the select operators and block preconditioning achieve the optimal matrix-query complexity O(κAlogε1) without using variable-time amplitude amplification (VTAA). The same upper bound holds for queries to the right-hand-side preparation oracle.
We discover two complementary linear-combination-of-Hermitian-matrices (LCHM) formulations to achieve a general non-normal matrix eigenvalue transformation g(A). Firstly, for A=L+iH with Hermitian L and H, the vanilla LCHM formula represents g(A) as a kernel integral of g(i(H+kL)), and it contains linear-combination-of-Hamiltonian-simulation (LCHS) [An, Liu, Lin, Phys. Rev. Lett. 2023] as the special case for matrix exponentials. Secondly, for the angular Hermitian Xθ=cosθL+sinθH, the Weyl LCHM formula expresses g(A) via integrating g(eiθ(Xθ±i(I−Xθ2)1/2)). For the matrix power g(A)=Am, the Fourier projection of Weyl LCHM gives Am=π2∫0πeimθTm(Xθ)dθ=N2∑j=0N−1eimθjTm(Xθj),θj=Nπj,for every N>m with Chebyshev polynomial of Hermitian Tm(Xθ) and N samples. The discrete formula is exact, introduces no truncation and angular quadrature error, and offers O(1) post-selection weights.
LCHM formulas lead to new quantum eigenvalue transformation (QET) algorithms. For a degree-d polynomial pd(A) on ∣ψ⟩, our QET algorithm can achieve optimal Θ(d) circuit depth and optimal O(∣∣pd∣∣∞/∣∣pd(A)∣ψ⟩∣∣) post-selection repetitions. LCHM-based QETs unify various quantum linear algebraic problems with near-optimal O(dlog(d/ε)) Clifford+T gates, including driven ODEs (reduced to standard LCHS), iterative methods, resolvents, log(I+A), (λI+A)ν, Sign and ReLU transforms, and Faber approximation on noncircular domains.
Linear combination of Hamiltonian simulation (LCHS) provides an efficient method for implementing matrix exponentials e−tA on quantum computers. In this paper, we develop LCHS formulas for computing general matrix functions f(A) when f is analytic on the numerical range of A, with A possibly non-normal. The essential technical tool is Weyl calculus, which reduces the construction of LCHS formulas for noncommuting operators to scalar Fourier approximation problems. Our construction yields a quantum eigenvalue transformation algorithm with optimal O(logε1) query complexity scaling. Furthermore, our Weyl-calculus-based theory gives rise to an ansatz-free convex optimization framework that directly produces discrete LCHS formulas. This circumvents the inefficiencies of traditional quadrature rules and yields formulas highly optimized for coherent implementation on quantum computers. In addition, both our theory and optimization framework apply to the simulation of time-dependent dissipative ODE dtdψ(t)=−A(t)ψ(t), for which we achieve a 2.1× cost reduction over prior art.
Eigenvalue transformations, which include solving time-dependent differential equations as a special case, have a wide range of applications in scientific and engineering computation. While quantum algorithms for singular value transformations are well studied, eigenvalue transformations are distinct, especially for non-normal matrices. We propose an efficient quantum algorithm for performing a class of eigenvalue transformations that can be expressed as a certain type of matrix Laplace transformation. This allows us to significantly extend the recently developed linear combination of Hamiltonian simulation (LCHS) method [An, Liu, Lin, Phys. Rev. Lett. 131, 150603, 2023; An, Childs, Lin, arXiv:2312.03916] to represent a wider class of eigenvalue transformations, such as powers of the matrix inverse, A−k, and the exponential of the matrix inverse, e−A−1. The latter can be interpreted as the solution of a mass-matrix differential equation of the form Au′(t)=−u(t). We demonstrate that our eigenvalue transformation approach can solve this problem without explicitly inverting A, reducing the computational complexity.
We present an efficient, nearly optimal quantum algorithm for solving linear matrix differential equations, with applications to the simulation of open quantum systems and beyond. For unitary or dissipative dynamics, the algorithm computes an entry of the solution matrix with query complexity O(νLt/ε), where the constant ν depends on the problem parameters, L involves a time integral of upper bounds on the norms of evolution operators, and ε is the error. In particular, νL is linear in t for unitary dynamics and can be a constant for dissipative dynamics. Our result contrasts prior quantum approaches for differential equations that typically require exponential time for this problem due to the encoding in a quantum state, which can lead to exponentially small amplitudes. We demonstrate the utility of the algorithm through an end-to-end application, namely the simulation of dissipative dynamics for non-interacting fermions, which can be extended to other quantum and classical systems. We compare with classical algorithms and give evidence of polynomial quantum speedups for systems in a lattice, which become more pronounced for systems with long-range interactions and can be shown to be exponential in general. We also provide a lower bound of Ω(νLt/ε) for unitary or dissipative dynamics that proves our algorithm is optimal up to logarithmic factors.
Linear combination of Hamiltonian simulation (LCHS) provides an asymptotically optimal framework for the quantum simulation of linear non-unitary dynamics. The LCHS approach approximately implements the Peano--Baker propagator through an integral representation determined by a kernel function. We establish an LCHS approximation theorem for a broadened class of admissible kernel functions, defined by natural complex-analytic conditions. In the standard time-independent block-encoding access model, we explicitly construct, for each target error ε, the unique kernel function in this class minimizing the query-cost functional and show that the minimum of this functional has the asymptotic expansion π2e(logε1−loglogε1+o(1)) as ε↓0. The provably optimal leading coefficient 2e/π improves by a factor of approximately 2.1 on the coefficient obtained empirically by Low and Somma through numerical minimization over their family of optimally scaling kernel functions. We construct a spectrally convergent quadrature with explicit bounds on operator error and linear combination of unitaries (LCU) subnormalization, proving that this asymptotic improvement persists after discretization into a finite linear combination of unitaries suitable for implementation. Our results are based on Hardy space theory and convex analysis.