We present a quantum algorithm to solve systems of linear equations of the form Ax = b , where A is a tridiagonal Toeplitz matrix and b results from discretizing an analytic function, with a circuit complexity of O (1/√ε, poly (log κ, log N )), where N denotes the number of equations, ε is the accuracy, and κ the condition number. The repeat-until-success algorithm has to be run O (κ/(1-ε)) times to succeed, leveraging amplitude amplification, and needs to be sampled O (1/ε 2 ) times. Thus, the algorithm achieves an exponential improvement with respect to N over classical methods. In particular, we present efficient oracles for state preparation, Hamiltonian simulation, and a set of observables together with the corresponding error and complexity analyses. As the main result of this work, we show how to use Richardson extrapolation to enhance Hamiltonian simulation, resulting in an implementation of Quantum Phase Estimation (QPE) within the algorithm with 1/√ε circuits that can be run in parallel each with circuit complexity 1/√ ε instead of 1/ε. Furthermore, we analyze necessary conditions for the overall algorithm to achieve an exponential speedup compared to classical methods. Our approach is not limited to the considered setting and can be applied to more general problems where Hamiltonian simulation is approximated via product formulae, although our theoretical results would need to be extended accordingly. All the procedures presented are implemented with Qiskit and tested for small systems using classical simulation as well as using real quantum devices available through the IBM Quantum Experience.

We consider a bounded Lipschitz-polyhedron $\Omega\subset\mathbb{R}^3$ of general topology equipped with a tetrahedral triangulation that induces a mesh $\Gamma_h$ of the surface $\partial\Omega$. We seek a maximal set of surface edge cycles that are independent in $H_1(\Gamma_h,\mathbb{Z})$ and bounding with respect to the exterior of $\Omega$. We present an algorithm for constructing suitable 1-cycles in $\Gamma_h$: First, representatives of a basis of the homology group $H_1(\Gamma_h,\mathbb{Z})$ are constructed, merely using the combinatorial description of the surface mesh $\Gamma_h$. Then, a duality pairing based on linking numbers is used to determine those combinations that are bounding with respect to $\mathbb{R}^3\setminus\Omega$. This is the key to circumventing a triangulation of the exterior region $\mathbb{R}^3\setminus\Omega$ in the computations. For shape-regular, quasi-uniform families of meshes, the asymptotic complexity of the algorithm is shown to be O(N 2 ), where N is the number of edges of $\Gamma_h$. The scheme provides an essential preprocessing step for all boundary element methods for eddy current simulation, which rely on discrete divergence-free vectorfields and their description through stream functions.