Numerical methods based on finite differences for solving the wave equation are prone to various errors, a challenge that also applies to the semi-discrete case of a recently proposed quantum algorithm that preserves continuous time while discretizing space into a quantum state vector. The inherent limitations of quantum processing, which affect the generality and extensibility of the quantum algorithm compared to classical numerical solvers, are offset by the logarithmic scaling of the number of qubits with the number of grid points. The higher grid density in combination with the higher-order approximation of the continuous Laplacian operator improves the numerical accuracy. In contrast to the Dirichlet boundary condition, however, the Neumann boundary condition in the quantum algorithm cannot exceed a second-order Laplacian approximation. Moreover, the anisotropy error is inherent to the lattice geometry. In this work, we propose a quantum wave equation solver on a hexagonal grid, which exhibits intrinsic grid dispersion comparable to that of a fourth-order accurate method on a rectangular grid and provides a quadratic improvement in the anisotropy measures. We show that the graph Laplacian on a hexagonal grid embedding the Neumann boundary condition is a bidirected graph. To make it suitable for use as a system Hamiltonian, a symmetrization and factorization based on the asymmetric adjacency matrix is proposed. Dispersion relations for quantum algorithms solving the wave equation on different lattices are derived and compared. The algorithms are analyzed as a limiting case of Yee’s finite-difference time-domain method in the context of electromagnetic wave propagation.
Numerical methods based on finite differences for solving the wave equation are prone to various errors, a challenge that also applies to the semi-discrete case of a recently proposed quantum algorithm that preserves continuous time while discretizing space into a quantum state vector. The inherent limitations of quantum processing, which affect the generality and extensibility of the quantum algorithm compared to classical numerical solvers, are offset by the logarithmic scaling of the number of qubits with the number of grid points. The higher grid density in combination with the higher-order approximation of the continuous Laplacian operator improves the numerical accuracy. In contrast to the Dirichlet boundary condition, however, the Neumann boundary condition in the quantum algorithm cannot exceed a second-order Laplacian approximation. Moreover, the anisotropy error is inherent to the lattice geometry. In this work, we propose a quantum wave equation solver on a hexagonal grid, which exhibits intrinsic grid dispersion comparable to that of a fourth-order accurate method on a rectangular grid and provides a quadratic improvement in the anisotropy measures. We show that the graph Laplacian on a hexagonal grid embedding the Neumann boundary condition is a bidirected graph. To make it suitable for use as a system Hamiltonian, a symmetrization and factorization based on the asymmetric adjacency matrix is proposed. Dispersion relations for quantum algorithms solving the wave equation on different lattices are derived and compared. The algorithms are analyzed as a limiting case of Yee’s finite-difference time-domain method in the context of electromagnetic wave propagation.
We present a quantum algorithm for simulating the wave equation under Dirichlet and Neumann boundary conditions. The algorithm uses Hamiltonian simulation and quantum linear system algorithms as subroutines. It relies on factorizations of discretized Laplacian operators to allow for polynomially improved scaling in truncation errors and improved scaling for state preparation relative to general purpose quantum algorithms for solving linear differential equations. Relative to classical algorithms for simulating the D-dimensional wave equation, our quantum algorithm achieves exponential space savings and achieves a speedup which is polynomial for fixed D and exponential in D. We also consider using Hamiltonian simulation for Klein-Gordon equations and Maxwell's equations.
Quantum methods for solving partial differential equations (PDEs) generally fall into a few general categories, each leveraging different quantum mechanisms to achieve potential speedups. One well-known approach is the quantum linear system solver, epitomized by the Harrow-Hassidim-Lloyd (HHL) algorithm, which uses phase estimation and eigenvalue inversion to solve large, sparse linear systems arising from discretized linear PDEs [1]. This method can yield exponential improvements over classical solvers under certain conditions (mostly when the matrix is well-conditioned and sparse). More recent developments, such as Quantum Singular Value Transformation (QSVT), provide a framework for matrix inversion in a more flexible way, often improving on the cost and error dependence of earlier solvers [2]. Time-dependent PDEs can be tackled via Hamiltonian simulations. The Trotter-Suzuki method decomposes the exponential of a Hamiltonian into products of exponentials of simpler sub-Hamiltonians, allowing a piecewise approximation of continuous dynamics [3]. Continuous-time quantum walks, which treat a PDE operator as a Hamiltonian that governs the walk, similarly simulate processes like diffusion or wave propagation on discretized grids or graphs [4]. Beyond direct simulation, variational quantum algorithms have emerged as a promising route in the noisy intermediate-scale quantum era. These approaches encode the PDE solution as a parameterized quantum state and optimize the parameters via classical feedback, thereby minimizing a cost function tied to the PDE residual or energy functional [5]. Quantum annealing or adiabatic methods, meanwhile, can solve optimization problems derived from PDE constraints by adiabatically evolving a simple initial Hamiltonian to the final problem Hamiltonian and staying in the system's ground state [3]. In our study, we implement a new parameterized quantum model, specifically designed to solve large linear systems on quantum hardware. This model is less deep and computationally faster than the proposed one in [5], reducing the complexity and speeding up convergence times. Such improvements are particularly valuable in cases where classical algorithms may fail due to the problem's size or complexity. Ongoing research is focusing on generalized quantum models that can be adapted to practical PDE problems for solving large electromagnetic problems where classical computers can fail.
The discrete time quantum walk is a quantum cellular automaton whose wavefunction comprises pairs of complex numbers assigned to uniformly spaced points on a line. The wavefunction evolves through the application of an alternating sequence of unitary operators: streaming of wavefunction values to adjacent points, and a Hadamard-type unitary matrix to blend pairs of values at individual points. Each operator generates the exact evolution due to part of the Hamiltonian for the one-dimensional Dirac equation over a finite time step. Composing these operators thus creates a discrete approximation to the Dirac equation. However, the composition of two non-commuting operators creates a global splitting error proportional to the length of the time step. The global error can be reduced from first order to second order in the time step by a unitary pre- and post-processing of the initial conditions and final output. The algorithm then becomes equivalent to a symmetric composition, a Strang splitting, between the two operators. This paper describes a fourth-order accurate composition scheme using nine stages, the fewest possible when the lengths of the time steps employed in the different stages are constrained to be integer multiples of some base time step. Each stage is itself a symmetric composition between two operators. This fourth-order scheme produces quantitatively smaller errors for a typical benchmark problem on spatial lattices with 1024 or more points, and shows the expected fourth-order convergence on sufficiently fine lattices. It has greater accuracy, over sufficiently long times, than three better-known fourth-order composition schemes using fewer stages, but with lengths related by irrational coefficients. The truncation error for plane-wave solutions is due to an operator that separates into a resonant part proportional to the Hamiltonian, and a non-resonant part orthogonal to the Hamiltonian. The resonant part commutes with the exact evolution operator, so its error accumulates to grow linearly with time. The orthogonal part produces oscillations that remain bounded over many time steps. The nine-stage integer scheme has the smallest resonant truncation error of the four schemes, despite being the only scheme that can be implemented using local operations. The other schemes implement streaming by irrational fractions of the lattice spacing using discrete Fourier transforms. • A fourth-order accurate quantum algorithm using local operations on regular lattices. • Higher quantitative accuracy for relevant problem sizes than existing schemes. • Higher long-time accuracy than standard fourth-order schemes with global operations. • Discrete invariants give fourth and higher order approximations to the energy.
Modelling of electromagnetic wave propagation in telecommunications has evolved from empirical models to highly deterministic ray tracing and numerical methods. The extraordinary computational effort of these methods and their wave nature seem ideally suited to be approached by quantum algorithms. We first examine current progress by reviewing recent proposals in the field. We scrutinize potentially advantageous quantum architectures, ranging from mainstream gate-level computers and near-term, mid-scale noisy architectures with limited capabilities, to adiabatic and annealing approaches that are already in commercial use. We analyze the weaknesses and strengths of recent proposals. Beyond the core algorithm, mechanisms to bridge the quantum and classical worlds are of particular interest. Extremely diverse algorithm specifications, from those based on Hamiltonian simulations and emulation of variational optimization to the unconstrained binary formulation, are compared with the use of pure gate-level circuits and known quantum subroutines. We show that the graph Laplacian, given its ability to integrate boundary conditions, is uniquely suited for quantum propagation modelling algorithms rooted in differential numerical methods. Quantum computers could overcome the temporal and spatial limitations of classical methods for larger computational domains and, to some extent, address the problems of dispersion and stability in finite-difference approximations. The ability to express the solution of a problem as an eigenvalue problem turns out to be an advantage in the quantum world, where eigenvalues and eigenvectors are inextricably intertwined with quantum mechanics. In this paper, we identify the most promising techniques and scenarios that hold the greatest potential.
This thesis describes quantum algorithms for Hamiltonian simulation, ordinary differential equations (ODEs), and partial differential equations (PDEs). Product formulas are used to simulate Hamiltonians which can be expressed as a sum of terms which can each be simulated individually. By simulating each of these terms in sequence, the net effect approximately simulates the total Hamiltonian. We find that the error of product formulas can be improved by randomizing over the order in which the Hamiltonian terms are simulated. We prove that this approach is asymptotically better than ordinary product formulas and present numerical comparisons for small numbers of qubits. The ODE algorithm applies to the initial value problem for time-independent first order linear ODEs. We approximate the propagator of the ODE by a truncated Taylor series, and we encode the initial value problem in a large linear system. We solve this linear system with a quantum linear system algorithm (QLSA) whose output we perform a post-selective measurement on. The resulting state encodes the solution to the initial value problem. We prove that our algorithm is asymptotically optimal with respect to several system parameters. The PDE algorithms apply the finite difference method (FDM) to Poisson's equation, the wave equation, and the Klein-Gordon equation. We use high order FDM approximations of the Laplacian operator to develop linear systems for Poisson's equation in cubic volumes under periodic, Neumann, and Dirichlet boundary conditions. Using QLSAs, we output states encoding solutions to Poisson's equation. We prove that our algorithm is exponentially faster with respect to the spatial dimension than analogous classical algorithms. We also consider how high order Laplacian approximations can be used for simulating the wave and Klein-Gordon equations. We consider under what conditions it suffices to use Hamiltonian simulation for time evolution, and we propose an algorithm for these cases that uses QLSAs for state preparation and post-processing.
Discrete Laplacian operators arise ubiquitously in scientific computing and frequently appear in quantum algorithms for tasks such as linear algebra, Hamiltonian simulation, and partial differential equations. Block encoding provides the standard method for accessing matrix data within quantum circuits. Efficient implementations of such algorithms require efficient block encodings of the discretized operator. While several general-purpose techniques exist for block encoding arbitrary matrices, they usually require deep quantum circuits. Moreover, existing efficient constructions that exploit Laplacian structure are limited in scope, typically assuming fixed boundary conditions or uniform grid resolutions. In this work, we present a unified framework for efficiently block encoding finite-difference discretizations of the Laplacian that supports Dirichlet, periodic, and Neumann boundary conditions in arbitrary spatial dimensions. Our construction allows different boundary conditions and grid sizes to be specified independently along each coordinate axis, enabling mixed-boundary and anisotropic discretizations within a single modular circuit architecture. We provide analytical gate-complexity estimates and perform circuit-level benchmarks after transpilation to an IBM hardware gate set. Across one-, two-, and three-dimensional examples, the resulting circuits exhibit substantially lower gate counts and higher success probabilities when compared to certain existing approaches.