Emulation of large-scale qubit registers with a phase space approach

TL;DR

This work introduces Phase-Space Approximation (PSA) for qubits, a phase-space method that replaces full quantum dynamics with an ensemble of mean-field trajectories to enable classical simulation of large qubit registers. Each trajectory evolves MF-like equations, and the quantum density matrix is approximated by averaging over trajectories, with initial conditions sampled to reproduce exact initial observables and fluctuations; the method scales as $O(N_{\text{traj}} L^2)$ and is highly parallelizable. The authors benchmark PSA on the $k$-local transverse-field Ising model for a wide range of locality and dimensions (1D, 2D, 3D), including up to $L=2000$ qubits, and compare against exact solutions and MPS results. Results show PSA reliably captures one-qubit observables and equilibration trends, improves over MF for many cases, but has limitations for multi-qubit observables and strong all-to-all couplings. Overall, PSA provides a practical, scalable reference for validating quantum simulations on large-scale qubit devices and offers a complementary tool alongside tensor-network approaches.

Abstract

A phase-space approach is used and benchmarked for the simulation of the continuous-time evolution of large registers of qubits. It is based on a statistical ensemble of independent mean-field trajectories, where mean-field is introduced at the level of the qubits, substituting quantum fluctuations/correlations with classical ones. The approach only involves at worse a quadratic cost in the system size, allowing to simulate up to several thousands of qubits on a classical computer. It provides qualitatively accurate description of one-qubit observables evolutions, making it a useful reference in comparison to techniques limited to small qubit numbers. The predictive power is however less robust for multi-qubits observables. We benchmark the method on the $k$-local transverse-field Ising model (TFIM), considering a large variety of systems ranging from local to all-to-all interactions, and from weak to strong coupling regimes, with up to 2000 qubits. To showcase the versatility of the approach, simulations on 2D and 3D Ising models are also made.

Figures (14)

• Figure 1: Mean-field eigenenergies of the TFIM model given by Eq. (\ref{['eq:mf-ansatz-energy']}) for different values of $\eta$, in a chain of length $L = 12$. Thicker points indicate energy minima.
• Figure 2: Comparison of the exact (solid line) and MF (dotted line) evolutions of Pauli matrices expectation values, averaged over the $L = 20$ qubits. Specifically, the evolution of the observables $\hat{S}_O = {1\over L}\sum_i \hat{O}_i$, with $\hat{O} = \hat{X}, \hat{Y}, \hat{Z}$, is shown as a function of time for four different values of $\eta$ and three different $k$, with the initial state $\ket{+}^L$.
• Figure 3: Schematic representation of the PSA approach as the average of multiple stochastic MF trajectories, enabling better reproduction of quantum-state properties and subsequent evolution. Each point corresponds to a trajectory, while different colors guide the eyes. Using the mean-field theory for the potential energy surface shown would lead to no evolution (in this specific example where the initial state is $\ket{0}$). Indeed, in an MF approach (black filled circle), one cannot spontaneously break the symmetry of the problem, which would mean here that the system decides to move to the left or right side.
• Figure 4: Comparison of the exact (solid lines) and PSA (dashed lines) dynamics, obtained by averaging $N_\textrm{traj} = 10^{4}$ trajectories. Same observables and same parameters as in Fig. \ref{['fig:exa-mf-comparison']}. A comparison to Fig. \ref{['fig:exa-mf-comparison']} shows drastic improvements compared to mean field.
• Figure 5: Deviation from the exact dynamics $D_r(T)$, for different $(k,\eta)$ values. From left to right, the initial state is $|0\rangle^L$, $\vert + \rangle^L$ and $\vert 01 \rangle^{L/2}$, and MF (top row) is compared to PSA (bottom row). Results have been obtained for $L=10$ and final time $T=10 \; (h^{-1})$.