Quantum PDE Solvers in Practice: Application-Driven Benchmarking of the Heat Equation

Abstract: Quantum PDE solvers are difficult to evaluate in practice because published studies use different discretizations, output models, reconstruction rules, and hardware assumptions. We present a reproducible, application-driven benchmark for the 1-D Dirichlet heat equation that compares eleven kernels under the same problem instances and readout contract. The benchmark covers coherent linear solvers (HHL, QSVT, and QLS-Fourier), VQLS, imaginary-time methods (QITE, var-QITE, and AVQDS), real-time Hamiltonian simulation and unitary dilations (Hamiltonian simulation, Schade-Hamiltonian, and Schr"odingerisation), and the spectral quantum simulation method (QSM). We use three initial conditions, four grid sizes from $n=4$ to $7$ qubits ($N=16$ to $128$), a CFL-like ratio $r\approx0.4$, and final time $T=1$. Statevector, ideal-shot ($10^5$ shots per step), and noisy Aer backends separate algorithmic, sampling, and device-noise errors. On statevector, QSM and Schade-Hamiltonian reproduce the semi-discrete reference to floating-point precision, Schr"odingerisation reaches approximately $10<sup>{-4}$ error, and QITE is the strongest non-transform method for smooth data. Under the fixed-shot setting, HHL degrades to approximately $0.79$ relative $\ell_2$ error, while several low-depth or postselected methods become readout-limited. A norm-mismatch ablation attributes 23--29% of the $n=7$ smooth-initial-condition error of Hamiltonian simulation, AVQDS, and QLS-Fourier to reconstruction normalization. Compact observables, including total thermal energy and individual Fourier-mode weights, require 1--3 orders of magnitude fewer shots than full-field reconstruction. The resulting public benchmark provides a practical guide for selecting quantum PDE solvers.

• The paper establishes a unified benchmark for quantum PDE solvers on the heat equation, eliminating fragmented evaluation protocols by fixing grid, time, and reconstruction parameters.
• It systematically evaluates eleven algorithmic kernels across five paradigms, detailing performance metrics such as circuit depth, relative errors, and scaling challenges.
• The study identifies Pareto-optimal approaches and emphasizes that observable-focused strategies may outperform full-field reconstruction in practical quantum advantage applications.

Cross-Paradigm Benchmarking of Quantum PDE Solvers on the Heat Equation

Introduction and Benchmark Motivation

This paper establishes a unified, application-driven benchmark for quantum algorithms tasked with solving the one-dimensional Dirichlet heat equation. Prior literature has suffered from fragmented evaluation methodologies—algorithms are compared on disparate discretization scales, under varying hardware assumptions, output models, and reconstruction protocols—obfuscating direct performance comparisons and undermining claims of quantum advantage for PDE applications.

To address this, the benchmark assembles eleven algorithmic kernels classifiable into five paradigms: coherent quantum linear solvers (including HHL, QSVT, and QLS-Fourier), a NISQ-friendly variational linear system solver (VQLS), variational imaginary-time evolution methods (QITE, var-QITE, and AVQDS), real-time Hamiltonian simulation/unitary dilation strategies (Hamiltonian simulation, Schade-Hamiltonian, and Schrödingerisation), and a quantum spectral/transform method (QSM). The harness fixes the grid, time horizon, CFL-like ratio, initial conditions, and output reconstruction model, providing an unprecedented controlled environment to interrogate paradigm-dependent strengths and limitations.

Benchmark Construction and Evaluation Protocol

The benchmark solves the 1D heat equation with zero Dirichlet boundaries, encapsulating the canonical parabolic PDE:

$$\frac{\partial u}{\partial t} = \alpha \frac{\partial^2 u}{\partial x^2}$$

Spatial discretization yields a symmetric positive definite Laplacian, while time-stepping employs an implicit Euler integrator—guaranteeing a favorable, scale-stable condition number for the linear system ($\kappa(A)\approx 2.6$ at large $N$ with fixed $r$).

Each kernel is evaluated under three initial conditions designed to stress distinct spectral regimes: a broadband rectangular pulse, a low-pass Gaussian, and a sparse two-mode (bimodal) initial profile. The assessment is conducted across three backend models—statevector simulation, ideal (noiseless) circuit sampling, and a noisy, depolarization-based simulator—enabling isolation of intrinsic algorithmic errors, sampling-induced statistical errors, and decoherence-induced degradation, respectively. Shot budgets, circuit resource counts, and output models (full state tomography versus compact observables) are standardized for all comparisons.

Paradigm and Kernel Taxonomy

Coherent Linear Solvers (Paradigm A):

• HHL employs QPE with ancillary-controlled rotations; its performance is limited by phase estimation resolution and postselection, with high circuit depth.
• QSVT leverages polynomial singular value transformation, efficiently block-encoding the system matrix and offering error control via polynomial degree, yet incurring steep resource growth.
• QLS-Fourier realizes functional inversion through LCU in the Laplacian eigenbasis, made tractable only for small $n$ due to exponential LCU term growth.

NISQ Variational Linear Solver (Paradigm B):

• VQLS operationalizes minimization of a local Hadamard-test objective. Its accuracy is constrained by truncation of the system matrix's non-sparse Pauli expansion and ansatz expressivity; shallow circuits render it suitable for NISQ devices but noncompetitive at larger $n$.

Imaginary-Time Variational Dynamics (Paradigm C):

• QITE approximates the heat semigroup with compressed local Pauli evolutions; errors stem from truncation support and Trotterization.
• var-QITE and AVQDS extend QITE with variational parameterizations, the latter adaptively growing its ansatz in response to residual error, at the cost of increasing classical overhead.

Real-Time Hamiltonian Simulation/Unitary Dilation (Paradigm D):

• Hamiltonian simulation applies Trotterized real-time evolution, with dissipation introduced classically.
• Schade-Hamiltonian and Schrödingerisation embed the dissipative map into an expanded, unitary-conserving circuit, recovering non-unitary dynamics after measurement/postselection; both are benchmarked for their ability to preserve spectral structure in the solution.

Quantum Spectral/Transform Method (Paradigm E):

• QSM performs basis transformation via quantum DST, applies per-mode decay, then inverse transforms, achieving a near-exact mapping of the heat propagator.

Algorithmic Performance and Resource Scaling

Statevector Backend Results:

QSM and Schade-Hamiltonian kernels reproduce the semi-discrete (classical eigendecomposition) baseline to double precision across all $n$ and initial conditions. Schrödingerisation achieves relative errors on the order of $10^{-4}$. QITE is the most accurate non-transform variant for smooth initial data but is impaired by high-frequency content in the pulse case due to local-update truncation.

Notably, HHL and QSVT, while foundational, are strongly limited by phase estimation or polynomial degree, as evidenced by plateauing errors and steep circuit scaling (Figure 1). QLS-Fourier's resource requirements become prohibitive past $n=4$.

Figure 2: Pulse initial condition evolution under different quantum solver paradigms, showing time evolution of $T(x)$ by kernel and grid size.

Ideal Sampling and Noisy Backends:

Introduction of finite sampling (Figure 3) and simulated device noise lifts all error floors without rearranging the qualitative ordering. Methods relying on postselection become severely readout-limited (especially HHL), and noisy transpilation failure emerges for deep/complex circuits (HHL fails beyond $\kappa(A)\approx 2.6$0). Variational and QITE-like algorithms, with shallow circuits and lower 2Q depth, are more robust at small $\kappa(A)\approx 2.6$1 but do not scale to high precision for large grids.

Figure 1: Gate-count scaling (log scale) with register size, highlighting sharp growth for HHL, QSVT, QSM, and unitary dilation methods.

Norm-Mismatch Ablation:

The analysis demonstrates that for some kernels, a substantial fraction of the error for smooth initial conditions is an artifact of reconstruction and norm-tracking, not inherent to the unitary core. In particular, AVQDS, Hamiltonian simulation, and QLS-Fourier share a spurious “error plateau” traceable to output normalization—up to 29% of $\kappa(A)\approx 2.6$2 error at $\kappa(A)\approx 2.6$3 is solely due to norm drift.

Observable-Readout Advantage:

Recovering full-field amplitude vectors requires a number of shots that grows with both system size and the required precision (scaling as $\kappa(A)\approx 2.6$4). In contrast, compact observables, such as total energy or individual Fourier mode weights, are achievable at 1–3 orders of magnitude lower shot budgets—a significant practical implication for quantum advantage claims (Figure 3).

Figure 3: Estimated measurement shot budget (log scale) required to reach $\kappa(A)\approx 2.6$5 relative $\kappa(A)\approx 2.6$6 error for various kernels and output models.

Practical Implications, Theoretical Impacts, and Future Directions

The central practical insight is that quantum PDE solvers, under identical problem setup and measurement protocols, yield vastly different cost-accuracy tradeoffs. QSM and Schade-Hamiltonian—both leveraging structural knowledge of the Laplacian eigenbasis—are Pareto-optimal among tested methods, offering machine-precision solutions at moderate compiled two-qubit resource cost and without the error inflation associated with postselection or deep gate sequences.

No tested algorithm demonstrates credible full-field quantum advantage for the heat equation at grid sizes $\kappa(A)\approx 2.6$7 since tomographic reconstruction dominates total cost compared to efficient classical spectral solvers. However, the path to meaningful, application-level quantum benefit lies in tasks where the end goal is observable estimation rather than full state reconstruction. This includes physical quantities such as energy, boundary fluxes, or mode weights, directly accessible with polynomial resource overhead.

Theoretically, the benchmark’s norm-mismatch ablation exposes a subtle but critical methodological issue: output maps and reconstruction rules can themselves induce algorithmic ceilings or shared error plateaus, potentially misleading comparative assessments if not separated from core kernel performance.

Future extensions should target: (1) real-hardware runs using calibrated device models, (2) generalization to higher-dimensional, nonlinear, or non-self-adjoint PDEs, (3) exploration of better preconditioners and decomposition strategies for coherence-limited algorithms, and (4) integration of quantum-inspired tensor network methods as comparators for low-entanglement evolution.

Conclusion

The benchmark provides the most comprehensive landscape to date for quantum heat-equation solvers, establishing clear regime-dependent performance frontiers and resource scaling laws. QSM and Schade-Hamiltonian are the reference methods for maximal accuracy; QITE is suitable for moderate-depth studies on smooth data. For practical quantum advantage, observable-centric tasks—particularly those insensitive to full state recovery—are the most promising paradigm. The publicly released harness enables extensible, reproducible testing of new algorithms under a rigorously controlled set of conventions.