Fast and converged classical simulations of evidence for the utility of quantum computing before fault tolerance

Tomislav Begušić; Johnnie Gray; Garnet Kin-Lic Chan

Source-linked AI summary

Fast and converged classical simulations of evidence for the utility of quantum computing before fault tolerance

Tomislav Begušić, Johnnie Gray, Garnet Kin-Lic Chan

arXiv:2308.05077v3Aug 9, 2023quant-ph

TL;DR

The paper addresses whether approximate classical methods can accurately simulate observables from a 127-qubit kicked Ising experiment that exceeded exact classical simulation. It combines sparse Pauli dynamics with tensor-network methods, including mixed PEPS/PEPO evolution and Bethe-free-entropy belief propagation, and obtains converged observables with effective bond dimension above 16 million and absolute error below 0.01 without extrapolation.

Problem
The 127-qubit kicked Ising experiment was argued to exceed exact classical simulation and some approximate methods, while earlier approximations left uncertainty in the observables.
Method
The paper develops sparse Pauli dynamics and tensor-network algorithms using mixed Schrödinger/Heisenberg PEPS/PEPO representations with Bethe-free-entropy belief propagation.
Results
The methods converge the studied 127-qubit observables beyond experimental accuracy, reaching effective bond dimension 16,777,216 and absolute error below 0.01 without extrapolation.
Takeaways & Limitations
For this experiment, the classical data may be considered converged for practical purposes, while the work identifies inaccuracies in the experimental extrapolations.

Abstract

from arXiv · show

A recent quantum simulation of observables of the kicked Ising model on 127 qubits implemented circuits that exceed the capabilities of exact classical simulation. We show that several approximate classical methods, based on sparse Pauli dynamics and tensor network algorithms, can simulate these observables orders of magnitude faster than the quantum experiment, and can also be systematically converged beyond the experimental accuracy. Our most accurate technique combines a mixed Schrödinger and Heisenberg tensor network representation with the Bethe free entropy relation of belief propagation to compute expectation values with an effective wavefunction-operator sandwich bond dimension >16,000,000, achieving an absolute accuracy, without extrapolation, in the observables of <0.01, which is converged for many practical purposes. We thereby identify inaccuracies in the experimental extrapolations and suggest how future experiments can be implemented to increase the classical hardness.

INTRODUCTION

Exact classical simulations are becoming inadequate for large quantum circuits, motivating approximate methods for the 127-qubit kicked Ising experiment. The paper develops sparse Pauli and tensor-network approaches that converge observables beyond the experiment’s accuracy.

The motivating experiment used a 127-qubit Eagle heavy-hex processor to measure Pauli observables in kicked Ising circuits that were argued to exceed exact and some approximate classical simulations.
Earlier sparse Pauli simulations matched experimental extrapolations faster than the experiment, but residual classical and experimental errors left uncertainty in the observables.
Approximate classical methods converge all studied kicked Ising observables, including 20-step dynamics on 127 qubits, with uncertainty substantially below experimental error bars.
A mixed Schrödinger/Heisenberg tensor-network method with Bethe free entropy reaches effective bond dimension 16,777,216 and absolute observable error below 0.01 without extrapolation.For many practical applications, this is considered fully converged.
Sparse Pauli dynamics represents observables as weighted sums of Pauli operators, whose terms can proliferate under non-Clifford rotations but not under Clifford rotations.

2. Tensor network simulations

The tensor-network approach splits circuit evolution between a forward PEPS state and a backward PEPO operator, then contracts their overlap with belief propagation. A mixed split limits entanglement growth in both representations and converges fastest.

2. Tensor network simulations: The MIX method evolves PEPS and PEPO representations for half the circuit, limiting entanglement growth in both state and operator and converging fastest.Pure PEPS and PEPO methods place all evolution in one picture, whereas MIX uses τ = T/2.
2. Tensor network simulations: The method’s three stages are forward PEPS evolution, backward PEPO evolution, and contraction of the resulting PEPS/PEPO/PEPS expectation-value sandwich.
2. Tensor network simulations: L2BP compresses PEPS and PEPO bonds using local singular-value decompositions with belief-propagation gauges, scaling as O(χ^4) on heavy-hex.
2. Tensor network simulations: L1BP estimates the final overlap from the Bethe free entropy of converged belief-propagation messages, avoiding expensive exact contraction and scaling as O(χ^6) for MIX.
2. Tensor network simulations: Compared with earlier PEPS evolution, the lazy implementation evolves and compresses one or more full gate layers simultaneously, which the authors find improves accuracy.

B. Simulation of expectation values of the kicked Ising circuit

Approximate sparse-Pauli and tensor-network methods reproduce the kicked-Ising observables efficiently, converge beyond experimental accuracy, and expose errors in some experimental and classical extrapolations.

SPD and MIX simulations match the zero-noise-extrapolated quantum data while running faster than the quantum processor, with SPD about three orders of magnitude faster for selected observables.SPD averages 10 seconds per point on a laptop core, and MIX completes the reported fast simulations in under three minutes on a consumer GPU.
SPD reproduces shallow-circuit exact benchmarks with maximum error around 10^-3, while the other three observables reach maximum error around 10^-4.The fast SPD results take about 10 seconds per point on one laptop core; numerically exact magnetization requires about 15 minutes per point.
The converged classical methods show a maximum spread below 0.045 for 20-step ⟨Z62⟩, substantially smaller than the experimental extrapolation error bars.For the 9-step circuit, all methods agree within 0.01 with the exact tensor-network benchmark.
The MIX method is the most converged tensor-network approach, while PEPS and MIX converge monotonically with bond dimension and PEPO and SPD show non-monotonic behavior.Supplementary estimates place MIX’s residual error from bond dimension well below 10^-2, potentially as low as 10^-3.
An absolute accuracy better than 0.01 is achieved without extrapolation for the 20-step ⟨Z62⟩ observable using the MIX tensor-network method.The largest MIX calculations required approximately one day on an Nvidia A100 GPU.
The converged MIX estimates identify a significant experimental extrapolation error at 6π/32 and show that MPS and MPO results are not well converged.The reduced 31-qubit model also deviates from the reference by approximately 0.1 at intermediate rotation angles.

DISCUSSION

The study finds that approximate classical algorithms reproduce the 127-qubit kicked Ising experiment faster than the quantum experiment and beyond its accuracy, while identifying conditions for harder future demonstrations.

Classical algorithms simulate the 127-qubit kicked Ising expectation values faster than, and well beyond the accuracy of, current quantum experiments.
Low estimated global fidelity does not necessarily prevent accurate expectation values for some circuit parameters.
Future quantum-utility experiments should seek physically relevant observables that are more sensitive to global fidelity.
The reported classical-simulation advances may support quantum simulation beyond this specific kicked Ising setting.

Sparse Pauli dynamics

Sparse Pauli dynamics approximates Heisenberg evolution by retaining a thresholded sum of Pauli operators, exploiting Clifford structure while controlling truncation and computational cost.

Sparse Pauli dynamics: Sparse Pauli dynamics represents the observable as a coefficient-weighted sum of n-qubit Pauli operators and evolves it in the Heisenberg picture.
Sparse Pauli dynamics: Clifford recompilation separates rotations into Clifford and non-Clifford parts, allowing Clifford transformations of the observable and rotation generators before evolution.
Sparse Pauli dynamics: A non-Clifford rotation expands each anticommuting Pauli term with σP and updates coefficients using cosine and sine factors.
Sparse Pauli dynamics: The Pauli representation can grow exponentially in the worst case, with scaling O(2^N), motivating threshold-based truncation.
Sparse Pauli dynamics: After each gate, terms with coefficient magnitude below δ are removed, while newly generated terms meeting the threshold are inserted.
Sparse Pauli dynamics: Expectation values are obtained from identity-Z Pauli terms, while the operator Frobenius norm is the coefficient-vector 2-norm.

Tensor Network simulations

The tensor-network approach splits evolution between a forward PEPS state and a backward PEPO operator, then contracts their expectation-value sandwich. Lazy belief propagation compresses these networks and evaluates the final contraction while exploiting low entanglement and device geometry.

Tensor Network simulations: The MIX method evolves a PEPS forward and a PEPO backward before contracting their expectation-value sandwich at an intermediate time.The PEPS, PEPO, and MIX methods correspond to τ=T, τ=0, and τ=T/2, respectively.
Tensor Network simulations: Lazy 2-norm belief propagation compresses evolving PEPS and PEPO networks, while lazy 1-norm belief propagation contracts the final expectation value.The latter approximately computes non-local quantities, including high-weight observables and tensor-network norms.
Tensor Network simulations: The method represents the initial state as a PEPS, the target observable as a PEPO, and circuit layers as site-associated tensors before forward and backward evolution.Two-qubit gates are spatially decomposed so each tensor belongs to a single site and layer.
Tensor Network simulations: Belief propagation compresses tensor networks through iterative message passing while avoiding dense formation of effective site tensors.The lazy implementation improves memory and cost scaling, with exact contractions of unstructured networks as the main computational step.
Tensor Network simulations: The MIX overlap sandwich with χ=256 would require an effective flattened 2D bond dimension of approximately 16,000,000.The contraction is evaluated with L1BP rather than explicitly flattening the network.
Tensor Network simulations: Importance sampling from BP-decimated configurations provides unbiased Monte Carlo estimates for checking L1BP expectation values.Approximate sampling probabilities are combined with exact probabilities obtained by contraction.
Tensor Network simulations: The approach is implemented for arbitrary geometries and gates, although larger blocks, MPS messages, and improved environment correlations remain possible fidelity enhancements.The authors report sufficient fidelity across the relevant observables without these improvements.

Detailed descriptions of the lazy TN contractions

The detailed contraction scheme applies lazy belief propagation separately to PEPS evolution, PEPO evolution, and the final expectation-value network. Each case uses geometry-specific message structures and optimized contractions to control computational scaling.

Detailed descriptions of the lazy TN contractions: The tensor construction assigns one tensor per site and layer, with physical and unitary virtual indices of dimension 2 and evolved PEPS or PEPO bonds up to χ.Spatial SVD decomposition of two-qubit gates creates the site-local representation.
Detailed descriptions of the lazy TN contractions: PEPS evolution forms a two-layer 2-norm network whose lazy message update scales as O(χ^4), instead of explicitly forming tensors of size 64χ^6.The construction traces selected future-layer indices and compresses the evolving state with bond dimension χ.
Detailed descriptions of the lazy TN contractions: PEPO evolution compresses the operator together with the next gate layer using 2-norm belief propagation and messages of size 16χ^2.The optimized pairwise contractions scale as O(χ^4).
Detailed descriptions of the lazy TN contractions: The final MIX expectation network has three χ-dimensional bonds between neighboring sites, producing messages of total size χ^3.For this network, optimized message updates scale as O(χ^6), compared with χ^9 for direct contraction of three site tensors.
Detailed descriptions of the lazy TN contractions: Lazy BP updates messages in parallel and stops when the maximum message change falls below 5 × 10^-6 in single precision.The stopping criterion is applied across all messages.

Error analysis

The MIX method’s convergence was assessed with high-bond-dimension variability, extrapolation comparisons, and normalized-versus-unnormalized expectation values. All reported error estimates remained below 0.01.

Error analysis: The largest standard deviation among the three highest-χ results is 1.2 × 10^-3 at θh=10π/32.This standard-deviation metric is evaluated separately for each θh.
Error analysis: All MIX convergence error estimates remain below 0.01, including the larger errors obtained from linear extrapolation.Figure S1 compares results across θh and χ with extrapolation fits and associated error estimates.
Error analysis: For the 20-step ⟨Z62⟩ observable, unnormalized and norm-normalized expectation values are compared with their average.The normalized value divides by the appropriate PEPS–PEPO norm.

1. TN timing analysis

Timing analysis reports total runtimes for the MIX, PEPS, and PEPO tensor-network methods, including all algorithmic components. Automatic contraction-tree optimization is the second-largest time cost after the L2BP and L1BP contractions.

1. TN timing analysis: Total runtimes for MIX, PEPS, and PEPO include all algorithm components, while automatic contraction-tree optimization is the second-largest cost after the BP contractions.The timing data are not fully comparable because some main-text results used different CPUs or GPUs and were not all timed.

Supplementary Figures

The supplementary figures assess convergence, normalization, computational scaling, exact agreement, and timings for SPD and tensor-network simulations across circuit depths and observables.

Convergence: MIX convergence is evaluated by comparing the highest-χ expectation value with a linear extrapolation from the three largest available χ values.The difference Δ and standard deviation σ quantify extrapolation deviation and variation among the final three points.
Convergence: For MIX, unnormalized, normalized, and averaged expectation values are compared, with truncation dominating convergence for θh ≤ 6π/32.The truncation cutoff is κ = 5 × 10^-6.
Circuit-depth comparisons: Supplementary figures extend the main-text comparisons to broader x-axis ranges and to 9-step circuits for the relevant observables.Figures S3–S5 reproduce or extend main-text plots for alternative circuit depths.
Computational scaling: SPD scaling is measured through Pauli counts, Z-type Pauli counts, computational time, and observable Frobenius norm at three θh values.The caption states that total Pauli count and time scale roughly quadratically with inverse threshold, while Z-type count scales inversely.
Exact comparison: Exact and MIX simulations are also compared for the 31-qubit model.
Timing comparisons: Timing figures report costs for MIX, PEPS, and PEPO tensor-network algorithms across multi-qubit observables and ⟨Z62⟩ at depths 9 and 20.The timings include circuit construction, belief-propagation compression or contraction, and contraction-tree optimization.

Supplementary Tables

The supplementary tables document SPD thresholds, simulation timings, Pauli-operator counts, and tensor-network bond dimensions used in the reported experiments.

SPD parameters: Table S1 lists the SPD thresholds used for the simulations presented in Figures 3 and 4.
SPD timings: Table S2 reports wall times in minutes on four CPU cores for SPD simulations of high-weight Pauli observables.
SPD simulations: Table S3 gives SPD thresholds, final Pauli-operator counts, and six-core wall times for the 9-step and 20-step simulations of Figure 5.The table notes that the number of Pauli operators processed during simulation can substantially exceed the final count.
Tensor-network parameters: Table S4 records the maximum bond dimensions χ used for the tensor-network simulations and the common truncation threshold κ = 5 × 10^-6.For 9 steps, sufficient χ values are 256 for PEPS, 64 for PEPO, and 16 for MIX; accessible 20-step values vary with method and angle.