PaperlayerView exact revision on arXiv

Source-linked AI summary

A flexible high-performance simulator for verifying and benchmarking quantum circuits implemented on real hardware

Benjamin Villalonga, Sergio Boixo, Bron Nelson, Christopher Henze, Eleanor Rieffel, Rupak Biswas, Salvatore Mandrà

arXiv:1811.09599v3quant-ph

TL;DR

The paper addresses the need for classical simulation that can both verify NISQ hardware and benchmark quantum-supremacy sampling. It introduces qFlex, a flexible tensor-network simulator for exact and target-fidelity RQCs, and demonstrates large-scale Bristlecone and square-grid simulations, including a 20 PFLOPS run at 64% of peak capability.

  • Problem

    Quantum-supremacy experiments need optimized classical simulation for hardware verification and fair comparisons with noisy quantum-device sampling.

  • Method

    qFlex uses tensor-network contraction with path-based fidelity control and a fast sampling technique to simulate exact and low-fidelity RQCs on regular grids and Bristlecone.

  • Results

    20 PFLOPS at 64% of maximum performance was achieved for a 60-qubit Bristlecone sampling simulation, while qFlex computed up to 10^6 amplitudes for hard Bristlecone sub-lattices.

  • Takeaways & Limitations

    qFlex raises the classical simulation bar for supremacy experiments while providing exact amplitudes for verification and target-fidelity outputs for noisy-device benchmarking.

Abstract

from arXiv · show

Here we present qFlex, a flexible tensor network based quantum circuit simulator. qFlex can compute both exact amplitudes, essential for the verification of the quantum hardware, as well as low fidelity amplitudes, in order to mimic sampling from Noisy Intermediate-Scale Quantum (NISQ) devices. In this work, we focus on random quantum circuits (RQCs) in the range of sizes expected for supremacy experiments. Fidelity $f$ simulations are performed at a cost that is $1/f$ lower than perfect fidelity ones. We also present a technique to eliminate the overhead introduced by rejection sampling in most tensor network approaches. We benchmark the simulation of square lattices and Google's Bristlecone QPU. Our analysis is supported by extensive simulations on NASA HPC clusters Pleiades and Electra. For our most computationally demanding simulation, the two clusters combined reached a peak of 20 PFLOPS (single precision), i.e., $64\%$ of their maximum achievable performance, which represents the largest numerical computation in terms of sustained FLOPs and number of nodes utilized ever run on NASA HPC clusters. Finally, we introduce a novel multithreaded, cache-efficient tensor index permutation algorithm of general application.

I. INTRODUCTION

Quantum supremacy experiments require both optimized classical simulation to establish a computational bar and exact simulation to verify NISQ hardware. qFlex addresses these needs with a flexible tensor-network simulator supporting exact and low-fidelity amplitudes, robust RQC performance, and large-scale Bristlecone benchmarks.

  • NISQ devices may perform tasks beyond classical computers, making optimized classical simulation necessary for credible supremacy comparisons.Classical simulation also supports verification of quantum hardware up to the limits of classical computation.
  • Exact output amplitudes support hardware verification, while equal-fidelity sampling compares quantum devices with classical competitors.XEB requires exact amplitudes for sampled device bit-strings, whereas supremacy sampling compares outputs at the same fidelity.
  • qFlex uses tensor-network contraction to simulate exact and low-fidelity amplitudes for RQCs, including Google’s Bristlecone QPU.The simulator targets circuit sizes relevant to supremacy experiments and is designed for regular tensor-network topologies.
  • 1/f amplitudes at target fidelity f cost the same as one perfect-fidelity amplitude, enabling low-fidelity NISQ sampling without the usual rejection-sampling overhead.Fine-grained cuts balance memory requirements against parallel computations, and an alternative sampling technique provides the same 1/f speedup.
  • Bristlecone sub-lattices provide classically difficult topologies while minimizing qubits and gates, helping preserve quantum-device fidelity.The identified family includes Bristlecone-24, -30, -40, -48, -60, and -70.
  • 20 PFLOPS at 64% of maximum performance was reached for the most demanding 60-qubit Bristlecone simulation using NASA’s Pleiades and Electra clusters.The run used about 90% of the nodes and was described as NASA’s largest computation by peak PFLOPS and nodes used.

II. RESULTS

qFlex was benchmarked on NASA’s Pleiades and Electra clusters for verification and sampling workloads across square lattices and Bristlecone circuits. The experiments report large-scale performance, improved contraction procedures, and trade-offs against competing simulators.

  • Large-scale performance: 20 PFLOPS in single precision, or 64% of the combined Pleiades–Electra maximum, was sustained during the largest simulation.The run used over 3.2 million core-hours and was reported as NASA Ames’s largest simulation by node count and FLOPS rate.
  • Large-scale performance: Skylake nodes generally provided the best performance and belonged to Electra, which was substantially more time- and energy-efficient than Pleiades.The runtime distributions were measured across Broadwell, Haswell, Ivy Bridge, Sandy Bridge, and Skylake architectures.
  • Contraction performance: A revised contraction procedure for Bristlecone-48 and Bristlecone-70 was about twice as fast as the earlier procedure.The improved contractions were identified after the initial simulations and benchmarked separately.
  • Circuit difficulty: For target fidelity close to 0.5%, Bristlecone-60 was almost 10× harder to simulate than Bristlecone-64, while Bristlecone-64 was only 2× harder than the comparison case.These estimates concerned computing 10^6 batches of amplitudes on a single core for different node types.
  • Comparisons: Against Alibaba’s simulator, qFlex was 3.6× to 100× slower depending on the case, although the compared Alibaba circuits were estimated to be about 1000× easier to simulate.The comparison used revised random-circuit prescriptions whose difficulty differs substantially from those used by Alibaba.
  • Comparisons: Against MFIB, qFlex was 7× less efficient for 10^6 amplitudes at 0.51% fidelity on 7×7 grids, but scaled more favorably beyond 8×8 grids and for large-diameter Bristlecone circuits.MFIB had an advantage for many amplitudes but was limited by memory usage; the methods used comparable resources for about 10^5 amplitudes.
  • Comparisons: For the 7×7 circuit instance, qFlex was estimated to be 9.26× more efficient than GPQS, or 4.63× more efficient under an estimated single-precision adaptation of GPQS.The comparison highlights the cost of distributed inter-node tensor contractions relative to qFlex’s cut-based approach.

III. DISCUSSION

qFlex targets hard random quantum circuits with a flexible tensor-contraction approach that supports exact and finite-fidelity amplitudes. Its benchmarks show strong scaling and high utilization on Bristlecone and large regular grids, while performance comparisons depend substantially on circuit design and benchmark fairness.

  • Scope and capabilities: qFlex computes exact and noisy amplitudes for regular-grid RQCs, including Google Bristlecone circuits and sub-lattices.The simulator is optimized for rectangular grids and Bristlecone, and targets both perfect and specified finite fidelity.
  • Robustness and comparisons: The simulator is more robust to circuit modifications than methods exploiting weaknesses in particular RQC ensembles.Its runtimes are directly determined by the number of full lattices of two-qubit gates at a given depth.
  • Large-scale performance: 20 PFLOPS at 64% of maximum performance was reached on Pleiades and Electra for the most demanding Bristlecone simulation.The run used about 90% of the nodes and was reported as NASA’s largest computation by peak PFLOPS and node count at the time.
  • Robustness and comparisons: Compared with Ref. [28], qFlex is 7× less efficient on 7 × 7 grids but scales better beyond 8 × 8 and for large-diameter circuits.The favorable scaling includes Bristlecone, Bristlecone-60, and Bristlecone-70.
  • Robustness and comparisons: qFlex is 37× more efficient than Ref. [40] for Bristlecone-70 and more than 9× more efficient than Ref. [41] for a 7 × 7 circuit.These comparisons use depth 1+32+1 for Bristlecone-70 and depth 1 + 40 + 1 for the 7 × 7 circuit.
  • Large-scale performance: qFlex computed 10^6 amplitudes for 60-qubit Bristlecone sub-lattices at about 0.50% fidelity in effectively well below half a day.The result used Pleiades and Electra combined and targeted depth (1+32+1).
  • Circuit design implications: The paper proposes Bristlecone sub-lattices that increase simulation difficulty while keeping qubit and gate counts small to improve overall fidelity.The theoretical and numerical analyses find Bristlecone more demanding than rectangular grids with the same number of qubits.
  • Scope and future work: The approach is presented as well optimized for regular circuit topologies, with broader circuit classes reserved for future exploration.The authors specifically mention applications involving optimization, machine learning, many-body systems, materials science, and chemistry.

A. Revised set of Random Quantum Circuits

The revised RQC prescription generates hardware-constrained random circuits on square-grid topologies using a fixed cycle of two-qubit layers and conditional single-qubit gates. The construction avoids a simulation shortcut involving T gates after CZ gates and can represent Bristlecone as a tiled-grid subset.

  • Circuit definition: The RQC ensemble acts on n qubits in a Hilbert space of dimension N = 2^n and samples output bit-strings.The circuits are defined by a chosen depth and topology.
  • Circuit construction: The prescription applies an initial layer of H gates, followed by cyclic two-qubit layers and a final H layer.Circuit depth is written 1 + t + 1, explicitly representing the initial and final Hadamard layers.
  • Circuit construction: The two-qubit-gate schedule contains 8 layers cycled in a consistent order and can tile arbitrary 2D square grids.The Bristlecone architecture is a diamond-shaped subset of such a grid, and the simulations use CZ gates.
  • Circuit construction: When a qubit is idle after participating in the previous layer, the prescription randomly applies X1/2 or Y1/2 with equal probability.This conditional single-qubit operation is part of the revised layer-generation rules.
  • Revision rationale: The revised prescription avoids placing T gates after CZ gates because that structure can reduce the computational cost of simulation.CZ is the two-qubit gate used in the simulations.
  • Revision rationale: Replacing CZ gates with iSWAP gates makes equivalent-depth circuits harder to simulate because iSWAP has Schmidt rank 4 rather than 2.The paper states that a CZ circuit of depth 1 + t + 1 is equivalent in simulation cost to an iSWAP circuit of depth 1 + t/2 + 1.
  • Tensor-network representation: Figure 5 represents eight consecutive CZ layers, including single-qubit gates, as a 3D tensor grid with one CZ interaction per neighboring qubit per block.This block structure underlies the tensor-network contraction approach.

B. Overview of the simulator

qFlex represents random quantum circuits as regular tensor networks and contracts them to compute amplitudes. Its key simplification is contracting time blocks first, reducing the network to a two-dimensional grid whose complexity is largely independent of gate randomness.

  • qFlex represents one-qubit gates as rank-2 tensors, two-qubit gates as rank-4 tensors, and larger gates as rank-2n tensors.
  • Every eight circuit layers are contracted in the time direction onto an I × J two-dimensional grid of tensors.The intermediate network is a three-dimensional I × J × K grid, with K = t/8.
  • The resulting tensor-network layout is regular, so contraction complexity does not depend on the particular random-circuit instance.
  • Contracting the time direction first reduces memory requirements and leaves only the resulting two-dimensional grid for high-complexity contractions.

1. Contraction of the 3D tensor network

The simulator contracts the regularized three-dimensional network by first collapsing time columns and then using cuts to keep tensors within node memory limits. For Bristlecone-70, four cuts trade a large tensor for many independently parallel paths.

  • Contraction strategy: Contracting each K-direction column produces a tensor with at most four indexes of dimension 2^K, yielding an I × J two-dimensional grid.
  • Architecture dependence: Bristlecone-60 is a factor of 2^4 harder than Bristlecone-64 despite having four fewer qubits because it requires three rather than two cuts.
  • Architecture dependence: The method applies to Bristlecone sub-lattices because the topology is a sub-lattice of a square grid, with Bristlecone-70 formed by turning off two qubits.
  • Cutting strategy: Each cut fixes an index value and decomposes one contraction into independent lower-complexity contractions, making the resulting paths embarrassingly parallelizable.
  • Cutting strategy: Four cuts for Bristlecone-70 at depth (1+32+1) reduce the largest tensor from 2^(11×4) to 2^(7×4) entries while creating 2^4 contractions.
  • Path contraction: The Bristlecone-70 path contracts regions A and B separately, combines them into AB, and then contracts the path-independent region C to a scalar.
  • Cost considerations: Cut selection balances memory complexity against contraction time and can preserve a small number of high-arithmetic-intensity contractions.

2. Implementation of the simulator

The implementation targets CPU supercomputers by combining matrix multiplication with cache-efficient tensor-index permutations. Decomposing arbitrary permutations into left and right moves improves contraction performance, especially for multithreaded reordering-dominated workloads.

  • Implementation: The engine uses C++ on CPU-based supercomputers and identifies matrix multiplication and tensor reordering as the main implementation bottlenecks.
  • Cache-efficient permutations: Arbitrary tensor-index permutations are decomposed into cache-efficient left and right moves, with L5−R10−L5 working well for binary indexes.
  • Cache-efficient permutations: Figure 8 compares optimized naive, left-move, and right-move reorderings across γ values and thread counts on Broadwell nodes.
  • Performance: Cache-efficient index permutations improve runtimes over naive reordering for both Bristlecone-60 and 7 × 7 grid amplitude contractions.
  • Implementation: Memoization reduces repeated index-map generation because regular tensor networks reuse the same maps.
  • Performance: The resulting speedups range from under 5% for single-threaded matrix-multiplication-dominated contractions to well over 50% for multithreaded reordering-dominated contractions.

C. Fast sampling of bit-strings from low fidelity RQCs

qFlex supports low-fidelity random-circuit sampling to model noisy quantum processors, alongside exact amplitudes for verification. Its fidelity-matching methods provide a 1/f simulation speedup over exact-amplitude computation.

  • Exact amplitudes support quantum-device verification, while low-fidelity amplitudes support classical benchmarking of noisy-device sampling.
  • The simulator matches a target fidelity f using methods that provide a 1/f speedup relative to exact-amplitude computation.
  • Both low-fidelity methods can be adjusted to match the noisy quantum computer’s fidelity and cross entropy, producing equivalent random-circuit sampling.

1. Simulating low fidelity RQCs

The simulator produces low-fidelity amplitudes by summing only a fraction of tensor-network paths or by mixing exact-wave-function samples with random bit-strings. Path-norm variation can make achieved fidelity differ from the target, and few-path simulations may deviate from Porter–Thomas statistics.

  • Path-based simulation: Summing a fraction f of paths produces amplitudes at a computational cost that is fraction f of perfect-fidelity simulation.The path contributions are approximately orthogonal and have similar norms, yielding a wave-function with norm and fidelity f.
  • Path-based simulation: The path-based method’s fidelity can differ from (#paths)^-1 because path contributions have non-negligible norm variation.In the Bristlecone-60 simulation, pairwise mutual fidelity among 4096 paths was about 10^-6, but their norms still varied.
  • Path-based simulation: With only 1–21 paths, randomly selected paths can yield fidelity below or above the target.Summing an extensive subset suppresses norm variation, but the simulations here minimize cuts and therefore use few paths.
  • Distributional behavior: Low-fidelity amplitudes can have a larger tail than the Porter–Thomas distribution when only a small number of paths is used.The authors attribute this to circuit cuts acting as removed gates, increasing effective diameter and requiring greater depth for thermalization.
  • Mixture-based simulation: The alternative method samples from the exact wave-function with probability f and emits a random bit-string with probability 1−f.It avoids summing a fraction of paths while using perfect-fidelity amplitudes for the exact-wave-function samples.
  • Mixture-based simulation: The alternative method uses the same computational cost as the path-fraction method but is more robust in achieving the target fidelity.For Schrödinger–Feynman simulators, however, summing a fraction of paths is preferable because arbitrary-amplitude computation has a small overhead.

2. Fast sampling technique

qFlex removes most rejection-sampling overhead by computing batches of amplitudes at nearly the cost of one amplitude and accepting at most one sample per batch. The resulting procedure approaches one batch per desired sample under suitable batch sizes, while longer-depth circuits show correlations closer to zero.

  • Motivation: 10× rejection-sampling overhead is targeted by computing batches of amplitudes instead of isolated amplitudes.The conventional procedure requires 10×10^6 amplitudes to obtain 10^6 samples when accepting 1/M amplitudes on average with M=10.
  • Batch computation: A batch of 30 amplitudes costs about 10% more than one amplitude for Bristlecone-64 and -60, while 256 amplitudes cost about 15% more for Bristlecone-48 and -70.These measured costs replace theoretical 30× and 256× overheads with modest incremental computation.
  • Sampling procedure: The procedure chooses slightly over 10^6 random AB strings, samples NC C strings per batch, computes their amplitudes, shuffles each batch, and accepts at most one bit-string.Acceptance uses min[1,p(sABC)N/M], and limiting acceptance to one string per batch avoids spurious correlations.
  • Sampling efficiency: For M=10 and NC=30, batch acceptance is 95.76% and 1.045×10^6 batches are needed for 10^6 samples.With NC=60, acceptance rises to 99.82% and 1.002×10^6 batches are needed; NC=256 is virtually 100% acceptance with 1.00×10^6 batches.
  • Sampling condition: The method requires uncorrelated amplitudes for fixed sAB, and Bristlecone-24 tests show correlations approach zero and become Hamming-distance independent at depth (1+32+1).Pearson-coefficient distributions also tend toward Porter–Thomas behavior at longer depth.

3. Sampling from a non fully-thermalized Porter-Thomas distribution

Sampling error increases when the amplitude distribution has not converged to Porter–Thomas statistics and instead has a larger tail. The paper estimates this error from amplitudes exceeding the rejection threshold.

  • Error estimate: Non-Porter–Thomas tails are expected to increase frugal rejection-sampling error beyond the approximately 10^-4 error estimated for M=10 under Porter–Thomas statistics.The numerical estimate sums probability amplitudes larger than M/N, multiplies by N, and divides by the number of amplitudes considered.

D. Simulation of Bristlecone compared to rectangular grids

Bristlecone sub-lattices are harder to simulate than comparable rectangular grids because their topology inherits a large diameter and requires more cuts. qFlex manages this through tensor reuse, alternative contractions, and controlled memory–parallelism trade-offs.

  • Topology and cuts: Bristlecone’s diamond-shaped sub-lattices are particularly hard to simulate because they inherit the diameter of larger rectangular grids.Bristlecone-70 is a sub-lattice of a 10×11 grid, and its depth (1+32+1) simulation requires four cuts to keep tensors manageable.
  • Topology and cuts: Computational cost scales exponentially with the number of cuts, while splitting Bristlecone-60 in two halves cuts 40 CZ gates versus 32 for an 8×8 grid.The comparison uses circuits with the same depth.
  • Contraction strategy: The rectangular-grid contraction uses two cuts to divide a 7×7×(1+40+1) network into four tensors, each of dimension 2^30.The procedure iterates over right and bottom paths, reusing path-independent tensors and contracting AB with CD to obtain each path contribution.
  • Contraction strategy: Tensor reuse reduces repeated contraction work but requires substantial memory, creating a trade-off between reuse and memory usage.The authors report that reuse was profitable in practice despite this trade-off.
  • Fast sampling contraction: An alternative contraction leaves six qubits free for batches of up to 2^6=64 amplitudes while reusing tensor pC.The method contracts D with AB first for a path, then computes amplitudes across the remaining batch strings.
  • Bristlecone contractions: A faster contraction for Bristlecone-48 and -70 runs about twice as fast as the earlier scheme and uses less memory, but cannot be adapted to Bristlecone-60.The earlier contractions remain relevant for the released simulation data.

Appendix C: Qubit complexity of square grids and Bristlecone sub-lattices for Schro¨dinger-Feynman-type simulators

The appendix quantifies qubit complexity for Schrödinger-Feynman-type simulators across different circuit partitions and compares square grids with Bristlecone sub-lattices. Bristlecone instances require higher complexity than comparable square grids, while partitioning choices and runtime contention affect practical performance.

  • Qubit-complexity definition: Qubit complexity measures the effective simulation cost from partition size, cut CZ gates, and sub-circuit wave-function dimensions.For two sub-circuits, each partition path requires simulating both sub-circuits, yielding log2[2^αAB(2^nA + 2^nB)].
  • Partition strategies: For three- and four-region partitions, reusing sub-circuit simulations across path combinations lowers qubit complexity relative to naive enumeration.The appendix gives reduced expressions for both three-region and four-region partitions by avoiding repeated simulation of some sub-circuits.
  • Square grids versus Bristlecone: 65 is the qubit complexity for an 8 × 8 square grid at depth (1+32+1), versus 71 for the best Bristlecone-60 partition.The square grid has four more qubits yet lower complexity, indicating that Bristlecone sub-lattices are harder to simulate than same- or smaller-sized rectangular grids.
  • Square-grid benchmarks: Slightly over 63 and 64 are the qubit complexities for 7 × 7 and 7 × 8 square grids at depth (1+40+1), respectively.Both configurations use two-sub-circuit splits.
  • Runtime considerations: Runtime differences across nodes were attributed to differing contention between the two sockets, despite thread-count equalization and thread pinning.The experimental setup could still leave some cores unused.
Loading 1811.09599v3…