PaperlayerView exact revision on arXiv

Source-linked AI summary

Improved classical simulation of quantum circuits dominated by Clifford gates

Sergey Bravyi, David Gosset

arXiv:1601.07601v3quant-ph

TL;DR

The paper asks how far classical simulation can extend beyond efficiently simulable Clifford circuits for Clifford+T computations. It develops probability-estimation and approximate-sampling algorithms based on stabilizer decompositions, obtaining mild exponential dependence on T-count and demonstrating a 40-qubit hidden-shift simulation.

  • Problem

    Classical simulation becomes impractical as generic state-vector costs grow exponentially, motivating better verification methods for quantum circuits dominated by Clifford gates.

  • Method

    The paper combines Clifford simulation with linear combinations of stabilizer states to estimate output probabilities and sample approximately from Clifford+T circuit distributions.

  • Results

    The algorithms have polynomial dependence on circuit size and mild exponential dependence on T-count, with approximate stabilizer decompositions achieving γ ≈ 0.228; a 40-qubit hidden-shift circuit with t ≤48 was simulated.

  • Takeaways & Limitations

    The methods may serve as verification tools for small quantum computations dominated by Clifford gates, including fault-tolerant circuits with costly logical T-gates.

Abstract

from arXiv · show

The Gottesman-Knill theorem asserts that a quantum circuit composed of Clifford gates can be efficiently simulated on a classical computer. Here we revisit this theorem and extend it to quantum circuits composed of Clifford and T gates, where T is the single-qubit 45-degree phase shift. We assume that the circuit outputs a bit string x obtained by measuring some subset of w qubits. Two simulation tasks are considered: (1) computing the probability of a given output x, and (2) sampling x from the output probability distribution. It is shown that these tasks can be solved on a classical computer in time $poly(n,m)+2^{0.5 t} t^3$ and $poly(n,m)+2^{0.23 t} t^3 w^3$ respectively, where t is the number of T-gates, m is the total number of gates, and n is the number of qubits. The proposed simulation algorithms may serve as a verification tool for medium-size quantum computations that are dominated by Clifford gates. The main ingredient of both algorithms is a subroutine for approximating the norm of an n-qubit state which is given as a linear combination of $χ$ stabilizer states. The subroutine runs in time $χn^3 ε^{-2}$, where $ε$ is the relative error. We also develop techniques for approximating tensor products of "magic states" by linear combinations of stabilizer states. To demonstrate the power of the new simulation methods, we performed a classical simulation of a hidden shift quantum algorithm with 40 qubits, a few hundred Clifford gates, and nearly 50 T-gates.

I. INTRODUCTION

The paper develops classical simulation algorithms for Clifford+T circuits whose cost is polynomial in qubits and Clifford gates but exponential in T-count, enabling verification of small Clifford-dominated computations. It reports improved stabilizer-decomposition bounds and a 40-qubit hidden-shift simulation with up to 48 T-gates.

  • Motivation: Classical state-vector simulation is limited to roughly 30 qubits, whereas Clifford circuits remain efficiently simulable through the Gottesman–Knill theorem.The Clifford+T gate set adds the non-Clifford T gate and is relevant to fault-tolerant quantum computing.
  • Main results: The paper presents two Clifford+T simulation algorithms whose runtimes are polynomial in n and Clifford-gate count but exponential in T-count.The algorithms target computing an output probability and sampling from an approximately correct output distribution.
  • Main results: The probability algorithm approximates Pout(x) for a specified output string with relative error ε and failure probability pf.Its runtime uses an exponential factor 2^(βt), with β ≤ 1/2 under the stated implementation-dependent bound.
  • Main results: The sampling algorithm produces a distribution ε-close to Pout in L1-norm, with practical scaling expected when the output width w is small and ε is not too small.For single-bit outputs and fixed precision, the runtime benefits from the favorable w = 1 setting.
  • Techniques: The algorithms represent magic-state tensor products as linear combinations of stabilizer states and simulate each Clifford-only term using Gottesman–Knill methods.The stabilizer-rank bound uses β ≈ 0.47 for exact decompositions, while approximate decompositions achieve γ ≈ 0.228 with χt(δ) = O(2^(γt)δ^-1).
  • Demonstration: A MATLAB implementation simulated hidden-shift circuits on 40 qubits with a few hundred Clifford gates and T-counts t = 40 and t = 48.The output probabilities concentrated on the corresponding hidden-shift bits, with the first 20 bits perfectly recovered for the chosen bent functions.

II. SKETCH OF TECHNIQUES

The simulation replaces T gates with Clifford-based gadgets, represents magic states using stabilizer decompositions, and estimates post-measurement norms. This decomposition supports sampling even when individual probabilities cannot be computed accurately.

  • Each T gate is replaced by a Clifford gadget consuming one magic state, with measurement outcomes handled through controlled Clifford operations.
  • The magic state A^⊗t is represented exactly or approximately as a linear combination of χ stabilizer states.Each resulting Clifford circuit branch remains efficiently simulable using the Gottesman-Knill theorem.
  • Measurement and postselection preserve stabilizer structure, leaving a final state expressed as a combination of at most χ stabilizer states.
  • The most time-consuming step computes the norm of the post-measurement state to obtain the required simulation quantity.
  • Stabilizer decompositions cannot accurately compute individual output probabilities but can support sampling with small statistical error.

III. STABILIZER FORMALISM

The stabilizer formalism represents Clifford-generated states through Pauli stabilizer groups and provides efficient operations on their standard representations. Stabilizer states also form a structured basis used by the simulation methods.

  • A stabilizer group is an abelian Pauli subgroup that excludes −I and is generated by pairwise commuting self-adjoint Pauli operators.
  • A pure state is a stabilizer state exactly when it equals U|0^n⟩ for a Clifford unitary U.
  • Each stabilizer state is uniquely specified up to global phase by an n-dimensional stabilizer group whose only stabilized state is that state.
  • The stabilizer-state set forms a 2-design under Haar-measure comparison with all normalized states.
  • In standard form, stabilizer-state inner products can be computed in O(n^3), alongside efficient projections onto stabilizer codespaces.

IV. CLASSICAL SIMULATION ALGORITHMS

The algorithms simulate Clifford+T circuits by replacing T gates with Clifford gadgets, decomposing magic states into stabilizer states, and estimating resulting norms or output distributions. Their costs are polynomial in circuit size but depend exponentially on T-count through the stabilizer-decomposition rank.

  • Probability estimation: The probability algorithm computes a stabilizer group and integer u, then approximates Pout(x)=∥ψ∥2 using randomized inner products with stabilizer states.The norm-estimation procedure uses random stabilizer states and median amplification to reduce failure probability.
  • Algorithmic setup: The T-gate gadget replaces each T gate with Clifford operations, a computational-basis measurement, and one consumed magic state |A⟩.Measurement outcomes occur uniformly, and the S gate is classically controlled by the outcome.
  • Algorithmic setup: Postselecting all gadget measurements at 0 removes classically controlled S gates and produces a Clifford circuit with c+t gates on n+t qubits.The ancillary magic-state qubits are appended to the computational qubits.
  • Algorithmic setup: The initial magic state A⊗t is represented by χ stabilizer states, allowing each Clifford-evolved term to be simulated efficiently with the Gottesman–Knill theorem.The resulting state is a linear combination of stabilizer states, with χ controlling the dominant runtime.
  • Probability estimation: O(χt^3ϵ^-2 log(p_f^-1)) time suffices to approximate the output probability with relative error ϵ and failure probability p_f.Each stabilizer-state inner product is computable in O(t^3), and the probability algorithm combines this with the preliminary stabilizer-group computation.

V. APPROXIMATING MAGIC STATES

This section develops approximate stabilizer decompositions for tensor powers of magic states, using random linear subspaces to control approximation error and stabilizer count. It establishes exponential scaling bounds and compares the chosen basis with the computational basis.

  • Approximate decompositions: χ = O(2^γt ε^-2) stabilizer states suffice for a decomposition with approximation parameter ε, where γ is specified by Eq. (4).The construction is presented as a systematic method for finding approximate decompositions.
  • Approximate decompositions: The magic state |A⟩ is Clifford-equivalent to |H⟩ up to a global phase, so the construction can target tensor powers of |H⟩.The paper uses |A⟩ = e^iπ/8 HS†|H⟩ to reduce the decomposition problem.
  • Random-subspace construction: The approximation is chosen as a uniform superposition over a k-dimensional linear subspace, using χ = 2^k stabilizer states.The relevant subspaces belong to L(t,k), the set of k-dimensional linear subspaces over the binary vector space.
  • Random-subspace construction: 2^k = O(δ^-1ν^-2t) can be achieved for constant δ, matching the best possible scaling of k with t for this construction.A lower bound follows from Z(L) ≥ 1, while randomized subspace selection produces a suitable L⋆.
  • Basis comparison: χ ∼ 2^tH2(ν2) ≈ 2^0.6t in the computational basis, which is significantly worse than the alternative stabilizer basis.The computational-basis expansion is concentrated near Hamming weight (1−ν2)t ± O(t^1/2).
  • Lower bounds: The paper proves a lower bound χ ≥ ν^-2t f^2∥z∥^-2 for normalized approximate decompositions and conjectures an Ω(ν^-2t) bound at constant error.The conjecture concerns any approximate stabilizer decomposition of H⊗t with constant approximation error.

APPENDIX A: QUADRATIC FORMS

Appendix A defines Z8-valued quadratic forms over affine binary spaces and develops canonical-basis techniques for evaluating their associated exponential sums. These sums provide an O(k^3) subroutine used in stabilizer-state inner products.

  • Definitions: A Z8-valued quadratic form q on an affine space is represented through a symmetric bilinear function J and basis-dependent coefficients.Arithmetic is performed modulo eight, while binary-vector addition is modulo two.
  • Quadratic-form structure: J(x,y) is symmetric and bilinear, with values 0 or 4 modulo 8; consequently, q(x⊕z)−q(z) lies in {0,2,4,6}.These properties follow from the defining identity for quadratic forms.
  • Stabilizer representation: Any stabilizer state has a unique affine-space and quadratic-form representation, linking stabilizer states to the coefficient data used by the appendix algorithms.The representation uses an affine space K and a quadratic form q on K.
  • Exponential sums: O(k^3) time suffices to transform a basis into canonical form and evaluate the associated exponential sum W(q).Canonicalization uses sparse basis changes, after which the sum factorizes into O(k) constant-time terms.

APPENDIX B: STANDARD FORM OF STABILIZER STATES

Appendix B gives a standard data representation for stabilizer states and an algorithm for restricting a state to an affine subspace. The SHRINK procedure updates the representation or reports when the restriction is empty or unchanged.

  • Restriction algorithm: The SHRINK operation restricts the affine space using a vector ξ and target bit α, producing an empty space, the original space, or a space of dimension k−1.The quadratic form is restricted to the resulting affine space.
  • Restriction algorithm: The algorithm updates basis vectors, dual vectors, the shift vector, and quadratic-form coefficients before removing the constrained coordinate.Its reported runtime is O(kn).
  • Restriction algorithm: The lazy SHRINK* variant changes the affine-space data without updating the quadratic-form coefficients.The appendix uses this variant when coefficient updates can be deferred.

APPENDIX C: THE INNER PRODUCT

Appendix C develops an algorithm for computing inner products between stabilizer states by intersecting their affine supports and evaluating the difference of their quadratic forms. The generic runtime is O(n^3), with a linear-in-n special case when both support dimensions are constant.

  • Inner-product construction: The inner-product algorithm first intersects the affine supports of the two stabilizer states and expresses both quadratic forms in a common basis.The resulting form is q ≡ q1 − q2.
  • Complexity: O(n^3) is the overall runtime of the inner-product algorithm in the general case.The runtime combines support-intersection, coefficient-update, and exponential-sum costs.
  • Inner-product construction: The intersection and coefficient updates are combined with an exponential-sum calculation to obtain the stabilizer-state inner product.The final expression uses W(Q,D,J).
  • Complexity: The algorithm is linear in n when both stabilizer-state support dimensions k1 and k2 are O(1), but most stabilizer states have kα ≈ n.Thus the constant-dimension shortcut provides no speedup in the generic case.
  • Timing analysis: The timing analysis evaluates inner products ⟨x̃|φ⟩ for uniformly random strings x and uniformly distributed stabilizer states φ.This input family is chosen because the simulation algorithm requires inner products of this form.

APPENDIX D: RANDOM STABILIZER STATES

The appendix constructs uniformly random n-qubit stabilizer states by sampling affine spaces and quadratic forms, with an average runtime of O(n^2) and worst-case runtime of O(n^3).

  • State construction: The algorithm samples stabilizer states by first choosing an affine space K and then a uniformly random quadratic form q on K.K is generated from a random kernel, shifted by a uniformly random vector, and represented in standard form.
  • State construction: The dimension parameter d is sampled from a distribution computed from counts of stabilizer states with different support dimensions.The appendix derives the relevant counts from linear subspaces, affine shifts, and quadratic-form choices.
  • Implementation: Rank testing and repeated kernel sampling produce a uniformly distributed subspace, while Gaussian elimination computes rank in O(nd^2) time.Full rank occurs after O(1) attempts with high probability when d=O(1).
  • Implementation: The construction is summarized as RandomStabilizerState(n), which samples d, generates K, chooses quadratic-form coefficients, and returns the standard-form state description.The returned description includes the affine-space and quadratic-form parameters needed to represent the stabilizer state.
  • Algorithm and runtime: O(n^2) average and O(n^3) worst-case runtime are achieved by the random stabilizer-state generation algorithm.The runtime depends on sampling dimension d, with d=O(1) occurring with high probability.

APPENDIX E: PAULI MEASUREMENTS

The appendix develops an O(n^2)-runtime procedure for measuring Pauli operators on stabilizer states, returning both the projected-state norm and, when nonzero, its standard form.

  • Projected-state computation: O(n^2) runtime suffices to compute the normalization and standard form of P+|K,q⟩ for a Pauli projector.The procedure also applies to products of up to n projectors describing stabilizer-code spaces.
  • Case analysis: The algorithm separates the cases ξ∈L(K) and ξ∉L(K), checking membership through inner products and coordinate reconstruction in O(kn) time.The second case extends the affine space by one dimension when the Pauli shift lies outside L(K).
  • Case analysis: When ξ∈L(K), the projected state is represented using a modified quadratic form obtained from the Pauli phase and the original form coefficients.The coefficients are computed in the chosen basis using the bilinear form J and the vector ζ.
  • Case analysis: When ξ∉L(K), the states |K,q⟩ and P|K,q⟩ have disjoint computational-basis supports, and the projected state is represented on an extended affine space M.The resulting normalized state can be computed in O(n^2) time.
  • Projected-state computation: MeasurePauli takes a stabilizer state and Pauli operator as input, returns the projected-state norm Γ, and standardizes the projected state when Γ≠0.The algorithm modifies the data describing the input state during computation.

APPENDIX F: SIMULATION OF THE HIDDEN SHIFT ALGORITHM

The appendix implements the hidden-shift benchmark using Clifford+T decompositions and verifies simulation on 40-qubit circuits with controlled T-counts. The reported output distributions recover the hidden shift bits with high concentration.

  • Circuit construction: The benchmark circuits used five and six CCZ gates, alternating with random sequences of 200 Clifford gates, producing the left and right plots respectively.The total T-count of the full circuit U is twice the T-count of Og.
  • Circuit construction: Each Toffoli decomposition used a Jones gadget with four T-gates, two |0⟩ ancillas, Clifford gates, and a classically controlled final Clifford operation.Ancilla measurements were simulated by postselection on randomly chosen output bits.
Loading 1601.07601v3…