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

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.

💡 Research Summary

The authors present a phase‑space based simulation technique, termed Phase‑Space Approximation (PSA), for the continuous‑time dynamics of large qubit registers. PSA builds on a statistical ensemble of independent mean‑field (MF) trajectories: each qubit is represented by a Bloch vector (x, y, z) that evolves according to the MF equations of motion, but the initial Bloch vectors are sampled from a probability distribution chosen so that the ensemble average reproduces the exact initial quantum state and its moments. Because each trajectory is separable (the full density matrix is a tensor product of single‑qubit density matrices) the computational cost scales as O(N_traj · L²), where L is the number of qubits; the method is trivially parallelizable and therefore can handle thousands of qubits on a conventional workstation.

The paper benchmarks PSA on the k‑local transverse‑field Ising model (TFIM) with Hamiltonian
H = −h ∑₁ᴸ Z_i − J ∑_{i<j,|i−j|≤k} X_i X_j.
The dimensionless parameter η = J·k/(h·L) controls the relative strength of interactions and connectivity. By varying η and the locality parameter k (from nearest‑neighbor to all‑to‑all), the authors explore weak‑coupling (η ≪ 1), strong‑coupling (η ≫ 1), and intermediate regimes, including the quantum phase transition between paramagnetic and ferromagnetic phases.

Key findings:

1. Single‑qubit observables (⟨X_i⟩, ⟨Y_i⟩, ⟨Z_i⟩) are reproduced with high fidelity. PSA consistently outperforms plain MF, reducing average errors by a factor of two to three across a broad range of η and k. The method is especially accurate when (a) η ≪ 1 with large k (weak coupling but many neighbors) and (b) η ≫ 1 (strong coupling where the initial state is close to an eigenstate of H).

2. Multi‑qubit correlators (e.g., strings like X₁X₂…Y_L) are not reliably captured. Because the sampling treats x, y, z as independent real variables, non‑commuting operator relations (e.g., ⟨X_iY_i⟩ = i⟨Z_i⟩) cannot be satisfied, leading to substantial deviations in higher‑order moments. Consequently, PSA’s predictive power diminishes near critical points where entanglement proliferates.

3. Scalability: Simulations up to L = 2000 qubits in one dimension, as well as 2D (32 × 32) and 3D (16 × 16 × 16) lattices, demonstrate that the quadratic scaling remains manageable. Memory usage stays linear in L because the full density matrix never needs to be stored; only the Bloch vectors for each trajectory are kept.

4. Comparison with other classical methods: Tensor‑network approaches (MPS, PEPS) excel for low‑entanglement, low‑dimensional systems but become prohibitive for large L or higher dimensions. PSA, by contrast, is agnostic to geometry and dimensionality, offering a fast “sketch” of the dynamics. However, it cannot replace tensor networks when precise multi‑qubit correlations are required.

5. Methodological insights: The PSA is mathematically equivalent to the discrete truncated Wigner approximation (dTWA) for spin‑½ systems, but the authors emphasize a density‑matrix formulation that aligns better with quantum‑computing terminology. The initial sampling respects all moments of single‑qubit operators, ensuring exact reproduction of the initial state, while the stochastic MF evolution introduces effective quantum fluctuations through ensemble averaging.

In conclusion, the paper establishes PSA as a practical tool for simulating the dynamics of very large qubit arrays with modest computational resources. It delivers accurate single‑qubit observables across a wide parameter space and can be applied to various interaction topologies and spatial dimensions. Its main limitation lies in capturing entanglement‑driven multi‑qubit observables, especially near phase transitions. Future work may focus on improving the sampling scheme to better respect non‑commuting operator algebra or on hybridizing PSA with more exact methods to extend its applicability to strongly correlated regimes.