Pushing the Classical Frontier of 1D Fermi-Hubbard Quench Dynamics Beyond Current Quantum Simulations

Abstract: Establishing quantum advantage requires comparison against the best achievable classical simulation. The Q-CTRL team recently simulated quench dynamics of the one-dimensional Fermi-Hubbard model on an IBM processor, completing a $L=60$ evolution to time $t=6$ in under three minutes and claiming a $3000\times$ speedup over classical Time-Dependent Variational Principle (TDVP) simulation at bond dimension $χ=4096$. Their classical benchmark required over 160 hours on a CPU cluster, failed to converge in the high-entanglement regime $t\in[5.2,6]$, and left the most challenging window of the experiment unverified. Here, we push the boundaries of classical simulation by exploiting the full $\mathrm{U}(1)\times\mathrm{SU}(2)$ symmetry of the Fermi-Hubbard Hamiltonian combined with GPU-accelerated tensor contractions. Reaching bond dimensions up to $χ\approx62{,}000$ on four NVIDIA H200 GPUs -- among the largest ever achieved in TDVP simulations and fifteen times larger than Q-CTRL's classical baseline -- we achieve fully converged results across the entire simulation window, including rigorous certification of the previously unresolved high-entanglement regime $t\in[5.2,6]$. We further advance the classical frontier to $t=7$, which lies beyond the quantum hardware experiment and any previously verified classical evolution of the full wavefunction. At the bond dimension comparable to Q-CTRL's best classical run, our GPU implementation completes in $\sim!100$ minutes, directly reducing the claimed $3000\times$ quantum advantage to $\sim!36\times$. These results substantially narrow the quantum-classical performance gap and establish a new standard for tensor-network benchmarking of large-scale quantum simulations.

• The paper demonstrates that exploiting U(1)xSU(2) symmetry in TDVP-MPS simulations drastically reduces computational complexity, enabling convergence in high-entanglement regimes.
• It shows that GPU acceleration and adaptive 1-site/2-site switching cut simulation times from over 160 hours to about 100 minutes, reducing the claimed quantum advantage from 3000× to 36×.
• The fully converged, site-resolved density profiles validate the simulation results across the critical time window, establishing a rigorous classical benchmark for quantum devices.

Pushing the Classical Boundary in 1D Fermi-Hubbard Quench Dynamics

Introduction

The delineation of quantum advantage hinges critically on the robustness of classical simulation baselines. The Fermi-Hubbard Model (FHM) is emblematic for benchmarking quantum devices, specifically due to the rapid entanglement growth upon quenching, Luttinger liquid phenomenology, and its tractable mapping onto qubit architectures. Prior demonstrations, notably by the Q-CTRL team, leveraged a 60-site ($L=60$) chain, with quantum hardware showing apparent speedup over Time-Dependent Variational Principle (TDVP)-based matrix product state (MPS) approaches at large but ultimately insufficient bond dimension ($\chi=4096$), leaving the high-entanglement regime $t\in[5.2,6]$ unverified. This paper delivers an explicit, rigorous closure of that certification gap, demonstrating that advanced symmetric, large-scale, and GPU-accelerated TDVP simulations can now comprehensively benchmark quantum hardware in these settings.

Model and Classical Simulation Methodology

The study uses the 1D Fermi-Hubbard Hamiltonian with open boundaries,

$$H = -t_h \sum_{\langle i,j\rangle,\sigma}\! \left(c^\dagger_{i\sigma}c_{j\sigma}+\mathrm{h.c.}\right) + U\sum_i n_{i\uparrow}n_{i\downarrow}$$

focusing on attractive $U=-2$ and half-filling, matching the quantum experiment protocol. The simulation proceeds from a half-filled Néel state using 30 Trotter steps with $\Delta t=0.2$.

The principal technical innovation lies in exploiting the larger $\mathrm{U}(1)\times\mathrm{SU}(2)$ symmetry of the FHM, as opposed to the $\mathrm{U}(1)\times\mathrm{U}(1)$ symmetry restriction applied in prior work. This symmetry, enforced via the Wigner–Eckart theorem, allows grouping MPS tensor indices into multiplets, entailing a drastic reduction in computational complexity while retaining full expressivity.

GPU-Accelerated Symmetric TDVP and Algorithmic Advances

Three key developments underpin the improved simulation:

• Full Non-Abelian Symmetry Exploitation: By encoding $\mathrm{U}(1)\times\mathrm{SU}(2)$ symmetry in MPS tensors, the number of degrees of freedom is reduced $\mathcal{O}(\chi^{1/3})$–$\chi=4096$0 compared to a non-symmetric representation, yielding a substantial reduction in effective bond dimension at fixed physical accuracy.
• High-Performance GPU Execution: The PyTorch-based implementation dispatches block-diagonal tensor contractions as batched GEMM operations, maximally utilizing GPU throughput, with multi-GPU parallelism allowing bond dimensions up to $\chi=4096$1. The simulation at $\chi=4096$2 (the expressivity-matched value to $\chi=4096$3 with symmetry) completes in roughly 100 minutes, a dramatic improvement over the >160 hours required by the previous ITensor+CPU baseline.
• Adaptive 1-site/2-site TDVP Switching: The simulation automatically transitions from 2-site to 1-site updates when bond-dimension growth saturates, eliminating unnecessary SVDs and yielding additional wall-clock speedups (by factors of $\chi=4096$4–$\chi=4096$5).

Main Results

The central numerical outcome is the achievement of fully converged TDVP-MPS simulation of the full quantum hardware window, including the previously inaccessible high-entanglement interval $\chi=4096$6, as well as extension beyond the quantum computation window to $\chi=4096$7. Absolute deviations at $\chi=4096$8 lie below $\chi=4096$9 for all $t\in[5.2,6]$0, and RMSE remains small for $t\in[5.2,6]$1 on key observables, securing the convergence claim. The well-known exponential $t\in[5.2,6]$2 growth is quantitatively resolved by controlling the SVD truncation threshold for various error tolerances.

Figure 1: (a) Time evolution of $t\in[5.2,6]$3 for $t\in[5.2,6]$4 up to $t\in[5.2,6]$5, comparing $t\in[5.2,6]$6 TDVP with previous baselines; (b) Absolute deviation from largest-$t\in[5.2,6]$7 reference; (c) RMSE versus Q-CTRL hardware data; (d) Wall-clock timings—CPU/ITensor, GPU/TDVP, QPU.

In terms of performance, the previous Q-CTRL claim of a $t\in[5.2,6]$8 quantum speedup is quantitatively challenged: at equivalent accuracy and expressivity, the classical simulation now completes in $t\in[5.2,6]$9 minutes, lowering the quantum advantage to $$H = -t_h \sum_{\langle i,j\rangle,\sigma}\! \left(c^\dagger_{i\sigma}c_{j\sigma}+\mathrm{h.c.}\right) + U\sum_i n_{i\uparrow}n_{i\downarrow}$$0.

Full Spatiotemporal Density Benchmark

An additional contribution is the complete, fully converged site-resolved density profile across all time windows validated by the Q-CTRL quantum experiment and its extension.

Figure 2: (a) Site and time-resolved $$H = -t_h \sum_{\langle i,j\rangle,\sigma}\! \left(c^\dagger_{i\sigma}c_{j\sigma}+\mathrm{h.c.}\right) + U\sum_i n_{i\uparrow}n_{i\downarrow}$$1 computed with $$H = -t_h \sum_{\langle i,j\rangle,\sigma}\! \left(c^\dagger_{i\sigma}c_{j\sigma}+\mathrm{h.c.}\right) + U\sum_i n_{i\uparrow}n_{i\downarrow}$$2 showing convergence up to $$H = -t_h \sum_{\langle i,j\rangle,\sigma}\! \left(c^\dagger_{i\sigma}c_{j\sigma}+\mathrm{h.c.}\right) + U\sum_i n_{i\uparrow}n_{i\downarrow}$$3; (b)-(d) Direct comparison at $$H = -t_h \sum_{\langle i,j\rangle,\sigma}\! \left(c^\dagger_{i\sigma}c_{j\sigma}+\mathrm{h.c.}\right) + U\sum_i n_{i\uparrow}n_{i\downarrow}$$4 against quantum hardware data.

The paper also demonstrates (see supplementary results) that $$H = -t_h \sum_{\langle i,j\rangle,\sigma}\! \left(c^\dagger_{i\sigma}c_{j\sigma}+\mathrm{h.c.}\right) + U\sum_i n_{i\uparrow}n_{i\downarrow}$$5 is insufficient for profile convergence at late times Figure 3, with substantial deviations visible. Converged results are achieved for $$H = -t_h \sum_{\langle i,j\rangle,\sigma}\! \left(c^\dagger_{i\sigma}c_{j\sigma}+\mathrm{h.c.}\right) + U\sum_i n_{i\uparrow}n_{i\downarrow}$$6 at $$H = -t_h \sum_{\langle i,j\rangle,\sigma}\! \left(c^\dagger_{i\sigma}c_{j\sigma}+\mathrm{h.c.}\right) + U\sum_i n_{i\uparrow}n_{i\downarrow}$$7 Figure 4.

Validation: Precision and Tolerance

The adequacy of float32 arithmetic is validated by direct comparison against float64, showing absolute differences in observables below $$H = -t_h \sum_{\langle i,j\rangle,\sigma}\! \left(c^\dagger_{i\sigma}c_{j\sigma}+\mathrm{h.c.}\right) + U\sum_i n_{i\uparrow}n_{i\downarrow}$$8; the error floor associated with quantum hardware noise renders double precision unnecessary for benchmarking. SVD truncation error analyses confirm that relative deviations remain within $$H = -t_h \sum_{\langle i,j\rangle,\sigma}\! \left(c^\dagger_{i\sigma}c_{j\sigma}+\mathrm{h.c.}\right) + U\sum_i n_{i\uparrow}n_{i\downarrow}$$9 for recommended tolerances and unrestricted bond growth.

Implications and Outlook

This work delivers unambiguous certification of quantum hardware for 1D Fermi-Hubbard quench dynamics across windows that were previously left unconstrained. Due to symmetry exploitation and algorithm–hardware co-design, the classical baseline is substantially advanced, with bond dimensions and wall-time performance unattainable in earlier approaches.

The claim of a $U=-2$0 quantum advantage for the Q-CTRL experiment is substantively reduced to $U=-2$1 under fair and modern classical baselining. This outcome forces reconsideration of quantum advantage claims for quench dynamics in low dimensions and highlights the need for ongoing development of classical tensor network frameworks, especially those harnessing symmetries and hardware accelerators.

While these advances further compress the quantum-classical frontier, future work may see additional narrowing via multi-GPU clusters, operator-specific classical algorithms (e.g., Majorana propagation), and further improvements in contraction and compression protocols. However, the exponential growth of entanglement entropy sets a hard asymptotic barrier for all classical tensor network algorithms, and simulations at still later times or in higher dimensions will ultimately saturate classical resources.

Conclusion

By combining symmetry-enhanced TDVP, GPU acceleration, and adaptive algorithmic strategies, this work establishes a new state-of-the-art for classical certification of quantum hardware output in the Fermi-Hubbard model. The approach closes all outstanding numerical and certification gaps left by previous benchmarks, reduces quantum advantage claims, and provides the first converged site- and time-resolved density map for $U=-2$2 Fermi-Hubbard dynamics up to $U=-2$3. The methodology sets a rigorous baseline for future quantum simulation experiments and offers insights guiding the next steps in the interplay between quantum hardware progress and classical simulation techniques.