Abstraqt: Analysis of Quantum Circuits via Abstract Stabilizer Simulation

Stabilizer simulation can efficiently simulate an important class of quantum circuits consisting exclusively of Clifford gates. However, all existing extensions of this simulation to arbitrary quantum circuits including non-Clifford gates suffer from an exponential runtime. To address this challenge, we present a novel approach for efficient stabilizer simulation on arbitrary quantum circuits, at the cost of lost precision. Our key idea is to compress an exponential sum representation of the quantum state into a single abstract summand covering (at least) all occurring summands. This allows us to introduce an abstract stabilizer simulator that efficiently manipulates abstract summands by over-approximating the effect of circuit operations including Clifford gates, non-Clifford gates, and (internal) measurements. We implemented our abstract simulator in a tool called Abstraqt and experimentally demonstrate that Abstraqt can establish circuit properties intractable for existing techniques.


Introduction
Stabilizer simulation [1] is a promising technique for efficient classical simulation of quantum circuits consisting exclusively of Clifford gates.Unfortunately, generalizing stabilizer simulation to arbitrary circuits including non-Clifford gates requires exponential time [2,3,4,5,6,7].Specifically, the first such generalization by Aaronson and Gottesman [2, §VII-C] tracks the quantum state ρ at any point in the quantum circuit as a sum whose number of summands m grows exponentially with the number of non-Clifford gates: Here, while c i , P i , b ij , and Q j can be represented efficiently (see §2), the overall representation is inefficient due to exponentially large m.
Abstraction.The key idea of Abstraqt is to avoid tracking the exact state ρ of a quantum system and instead only track key aspects of ρ.To this end, we rely on the established framework of abstract interpretation [8,9], which is traditionally used to analyze classical programs [10,11] or neural networks [12] by describing sets of possible states without explicitly enumerating all of them.Here, we use abstract interpretation to describe the set of quantum states that could occur at a specific point during execution of a circuit, by over-approximating the summands that could occur in any of those quantum states ρ.
Merging Summands.This allows us to curb the exponential blow-up of stabilizer simulation by merging multiple summands in Eq. ( 1) into an abstract single summand which over-approximates all summands, at the cost of lost precision.The key technical challenge addressed by our work is designing a suitable abstract domain to describe sets of summands, accompanied by the corresponding abstract transformers to over-approximate the actions performed by the original exponential stabilizer simulation on individual summands.

Background
We first introduce the necessary mathematical concepts.
Basic Notation.We use Z n := Z/(nZ), define B := Z 2 , and write 2 S for the power set of the set S.
We represent a pure n-qubit quantum state ψ ∈ C 2 n as a density matrix ρ ∈ C 2 n ×2 n , defined as ρ = ψψ † , where ψ † denotes the conjugate transpose of ψ.For a mixed state, i.e., a distribution over pure states ψ i with probability p i , the corresponding density matrix is ρ = i p i ψ i ψ † i .Because both ψ and ρ store exponentially many values, they cannot be represented explicitly for large n.We denote the embedding of a k-qubit gate U ∈ U(2 k ) as an n-qubit gate by U (i) := I 2 i ⊗ U ⊗ I 2 n−i−k , where I l denotes the l × l identity matrix.

Stabilizer Simulation.
The key idea of stabilizer simulation [1,2] is representing quantum states ρ = ψψ † implicitly, by stabilizers Q which stabilize the state ψ, that is Qψ = ψ.As shown in [2], appropriately selecting n stabilizers Q j then specifies a unique n-qubit state ρ = n j=1 I+Qj 2 .In stabilizer simulation, all Q j are Pauli elements from P n of the form i v • P (0) ⊗ • • • ⊗ P (n−1) , where P (j) ∈ {X, Y, Z, I 2 } and v ∈ Z 4 .This directly implies that all stabilizers Q i for the same state ψ commute, that is Q i Q j = Q j Q i , as elements from the Pauli group P n either commute or anti-commute.These elements can be represented efficiently in memory by storing v and P (0) , . . ., P (n−1) .In App.B, we list states stabilized by Pauli matrices (Tab.7) and the results of multiplying Pauli matrices (Tab.8).Further, in this work we use the functions bare b : P n → P n and prefactor f : P n → Z 4 which extract the Pauli matrices without the prefactor and the prefactor, respectively: b(i v P (0) ⊗ • • • ⊗ P (n−1) ) = P (0) ⊗ • • • ⊗ P (n−1) . ( Applying gate U to state ρ can be reduced to conjugating the stabilizers Q j with U : U † [13, Sec. 10.5]  = n j=1 While Eq. ( 4) holds for any gate U , stabilizer simulation can only exploit it if U Q j U † ∈ P n .
Clifford gates such as S, H, CN OT , I, X, Y, and Z satisfy this for any Q j ∈ P n .
To also support the application of non-Clifford gates such as T gates, we follow [2, §VII.C] and represent ρ more generally as for c i ∈ C, P i ∈ P n , b ij ∈ B, and Q j ∈ P n .Here, applying U to ρ amounts to replacing P i by U P i U † and Q j by U Q j U † , which we can exploit if both U P i U † and U Q j U † lie in P n .Otherwise, we decompose2 U to the sum Here, d * q denotes the complex conjugate of d q , + denotes addition modulo 2, and Q j ⋄ R q is the commutator defined as 0 if Q j and R q commute and 1 otherwise.Note that • ⋄ • : P n × P n → B has the highest precedence.
Overall, the decomposition of a k-qubit non-Clifford gate results in at most K = 4 k summands, thus blowing up the number of summands in our representation by at most 4 k • 4 k = 16 k .In practice, the blow-up is typically smaller, e.g., decomposing a T gate only requires 2 summands, while decomposing a CCN OT gate requires 8 summands.

Measurement.
Measuring in bare Pauli basis P ∈ b(P n ) yields one of two possible quantum states.They can be computed by applying the two projections P + := I+P 2 and P − = I−P 2 , resulting in states ρ + = P + ρP + and ρ − = P − ρP − , respectively.For example, collapsing the i th qubit to |0⟩ or |1⟩ corresponds to measuring in Pauli basis Z (i) .The probability of outcome ρ + is tr (ρ + ), and analogously for ρ − .Note that we avoid renormalization for simplicity.We discuss in §5 how measurements are performed in stabilizer simulation [2, Sec.VII.C].
Abstract Interpretation.Abstract interpretation [8] is a framework for formalizing approximate but sound calculation.An abstraction consists of ordered sets (2 X , ⊆) and (X , ≤), where X and X are called concrete set and abstract set respectively together with a concretization function γ : X → 2 X which indicates which concrete elements x = γ(x) ⊆ X are represented by the abstract element x.Additionally, ⊥ ∈ X refers to for all x ∈ X , where f was lifted to operate on subsets of X .This ensures that f ♯ (over-)approximates f , a property referred to as soundness of f ♯ .Abstract transformers can analogously be defined for functions f : X n → X .Further, we introduce join ⊔ : X ×X → X , satisfying γ(x)∪γ(y) ⊆ γ(x⊔y).Throughout this work, we distinguish abstract objects x ∈ X and concrete objects x ∈ X by stylizing them in bold or non-bold respectively.
As an example, a common abstraction is the interval abstraction with X = R.The abstract set is the set of intervals where x = (l, u) is a tuple.The concretization function γ : X → X maps these tuples to sets: Further, ⊤ = (−∞, ∞) and ⊥ = (l, u) for l > u.Common abstract transformers for the interval abstraction are shown in Tab. 1.
The transformers in Tab. 1 are precise, meaning that for f : R → R, we have that f ♯ ((l, u)) = (min l≤v≤u f (v), max l≤v≤u f (v)) and analogously for f : R n → R.An abstract transformer for a composition of functions f • g is the composition of the abstract transformers.Although this is sound, it is not necessarily precise: let g : R → R 2 with g(x) = ( x x ) and f : R ) whereas a precise transformer would map (−2, 2) to (0, 4).

Notational Convention.
In slight abuse of notation, throughout this work we may write the concretization of abstract elements instead of the abstract element itself.For example, for to indicate that it represents an interval.Where clear from context, we omit ♯ and write f for f ♯ .For example, we write

Overview
In this section, we showcase Abstraqt by applying it to the example circuit in Fig. 1.Overall, Abstraqt proceeds analogously to [2, §VII-C], but operates on abstract summands representing many concrete summands.
Example Circuit.We first discuss the circuit in Fig. 1.Both qubits are initialized to |0⟩.The circuit then applies a succession of gates.The abstract representation of the state after the application of each gate is shown in the gray boxes below the circuit.On the final state, the circuit collapses the upper qubit to |−⟩ by applying the projection M − = I−X (0) 2 . Precise circuit simulation shows that the probability of obtaining |−⟩ is 0, in this case.In the following, we demonstrate how Abstraqt computes an over-approximation of this probability.Initial State.The density matrix for the initial state |0⟩ ⊗ |0⟩ can be represented as (see [2]): .
We now explain how each operation in the circuit modifies this abstract state.
Non Clifford Gate Application.Next, the circuit applies gate T on the upper qubit.To this end, we again follow the simulation described in §2.We first decompose T into Pauli elements: , where d 1 ≈ e −0.1+0.4i and d 2 ≈ e −1.0−1.2i .Replacing T with its decomposition, we can then write ρ T = T ρ B T † , using Eq. ( 6), as: Analogously to §2, we can rewrite this to: , where Merging Summands.Unfortunately, simply applying T gates as shown above may thus increase the number of summands in the abstract density matrix by a factor of 4. To counteract this, our key idea is to merge summands, by allowing a single abstract summand to represent multiple concrete ones, resulting in reduced computation overhead at the cost of lost precision.Our abstract representation allows for a straightforward merge: we take the union of sets and join intervals.Specifically, for complex numbers, we join the intervals in their representation, obtaining: Finally, we introduce the symbol ⋆ to denote how many concrete summands an abstract summand represents.Altogether, merging the summands in ρ T yields: .
Note that for an abstract element x, r ⋆ x is not equivalent to r • x.For example, Measurement.After the T gate, the circuit applies two additional CN OT gates, resulting in the updated density matrix: Finally, the circuit applies the projection M − = I−X (0) 2 . To update the density matrix accordingly, we closely follow [2], which showed that measurement can be reduced to simple state updates through a case distinction on M − and the state ρ.If (i) the measurement Pauli (here −X (0) ) commutes with the product Paulis (here (−1) {0,1} X (0) X (1) and (−1) {0} X (1) ) and (ii) the measurement Pauli cannot be written as a product of the product Paulis, the density matrix after measurement is 0. We will explain in §5.2 how our abstract domain allows both of these checks to be performed efficiently.
Here, both conditions are satisfied, and we hence get the final state ρ M 1 = 0. We can then compute the probability of such an outcome by p = tr (ρ M 1 ) = 0. Thus, our abstract representation was able to provide a fully precise result.
Imprecise Measurement.Suppose now that instead of the measurement in Fig. 1, we had collapsed the lower qubit to |0⟩ by applying projection M 0 = . To derive the resulting state, we again follow [2] closely.We note that the measurement Pauli +Z (1) (i) anticommutes with the first product Pauli (−1) {0,1} X (0) X (1) and commutes with the second one (−1) {0} X (0) and (ii) commutes with the initial Paulis {I, Z (0) }.In this case, we get that the density matrix is unchanged, thus ρ M 2 = ρ D .To compute the trace of this matrix, we follow the procedure outlined in §5.4.We omit intermediate steps here and get: Thus, our abstraction here is highly imprecise and does not yield any information on the measurement result (we already knew that the probability must lie in [0, 1]).

Transformers Domains Definition
Eq. ( 12) P Q ∈ P P P P P P n P, Q ∈ P P P P P P n Eq. ( 13) Eq. ( 14)

Abstract Domains
In the following, we formalize all abstract domains (Tab. 2) underlying our abstract representation of density matrices ρ along with key abstract transformers operating on them (Tab.3).We note that all abstract transformers introduced here naturally also support (partially) concrete arguments.
Example Elements.Tab. 2 provides an example element x of each abstract domain, along with an example of its concretization γ(x), where γ : X → 2 X .While Tab. 2 correctly distinguishes abstract elements from their concretization, in the following, when describing operators we write concretizations instead of abstract elements (as announced in §2).

Booleans and Z
are subsets of B, as exemplified in Tab. 2. The addition of two abstract booleans naturally lifts boolean addition to sets and is clearly sound: We define multiplication of abstract booleans analogously.Further, we define the join of two abstract booleans as their set union.Analogously to booleans, our abstract domain Z Z Z Z Z Z 4 consists of subsets of Z 4 , where addition, subtraction, multiplication, and joins works analogously to abstract booleans.Further, we can straight-forwardly embed abstract booleans into Z Z Z Z Z Z 4 by mapping 0 to 0 and 1 to 1.
Real Numbers.We abstract real numbers by intervals of the form [a, a] ⊆ R ∪ {±∞}, and denote the set of such intervals by R R R R R R. Here, a and a indicate the lower and upper bounds of the interval, respectively.Interval addition, interval multiplication, and the cosine and exponential transformer on intervals are defined in their standard way, see §2.Complex Numbers.We parametrize complex numbers c ∈ C in polar coordinates (with magnitude in log-space), as c = e r+φi for r, φ ∈ R. For example, we parametrize 0 as e −∞+0i .
Based on this parametrization, we abstract complex numbers using two real intervals for r and φ respectively, as exemplified in Tab. 2. Formally, we interpret c ∈ C C C C C C as the set of all possible outcomes when instantiating both intervals: We can compute the multiplication and join of two abstract complex numbers c = e [r,r]+[φ,φ]i and c Again, simple arithmetic shows that Eqs. ( 9)-( 10) are sound.We note that to increase precision, we could map complex numbers to a canonical representation before joining them, by exploiting e r+ϕi = e r+(ϕ+2π)i to ensure that φ lies in [0, 2π].
We compute the real part of an abstract complex number c = e [r,r]+[φ,φ]i as where we rely on interval transformers to evaluate the right-hand side.The soundness of Eq. ( 11) follows from the standard formula to extract the real part from a complex number in polar coordinates.We will later use Eq. ( 11) to compute tr (ρ).To this end, we also need the transformer Pauli Elements.Recall that a Pauli element P ∈ P n has the form 1) , for v in Z 4 and P (k) ∈ {I, X, Y, Z}.We therefore parametrize P as a prefactor v (in log i space) and n bare Paulis P (k) .Accordingly, we parametrize abstract Pauli elements P ∈ P P P P P P n as , where v ∈ Z Z Z Z Z Z 4 is a set of possible prefactors and P (k) ⊆ {X, Y, Z, I 2 } are sets of possible Pauli matrices.Formally, we interpret P as the set of all possible outcomes when instantiating all sets: We define the product of two abstract Pauli elements as: To this end, we evaluate the prefactor induced by multiplying Paulis as where we can evaluate the summands in the right-hand side of Eq. ( 14) by precomputing them for all possible sets of Pauli matrices P (i) and Q (i) .Then, we compute the sum using Eq. ( 8).Analogously, we can evaluate b P (i) Q (i) by precomputation.The soundness of Eq. ( 13) follows from applying the multiplication component-wise, and then separating out prefactors from bare Paulis.
We also define the conjugation of an abstract Pauli element P with k-qubit gate U padded to n qubits as: where P (i:j) denotes P (i) ⊗ • • • ⊗ P (j−1) .Because k is typically small, and all possible gates U are known in advance, we can efficiently precompute f(U P (i:i+k) U † ) and b(U P (i:i+k) U † ).We note that this only works if the result of conjugation is indeed an (abstract) Pauli element-if not, this operation throws an error 5 .The soundness from Eq. ( 15) follows from applying U to qubits i through i + k, and then separating out prefactors from bare Paulis.We define the commutator P ⋄ Q of two abstract Pauli elements P and Q as Here, we evaluate the sum using Eq. ( 8), and efficiently evaluate B B by precomputing: The soundness of Eq. ( 16) can be derived from the corresponding concrete equation, which can be verified using standard linear algebra.
We define the join of abstract Pauli elements as where Clearly, this join is sound.Finally, we define an abstract transformer for modifying the sign of an abstract Pauli element P by: The soundness of Eq. ( 18) follows directly from (−1) v = i 2v .

Abstract
Here, r ∈ N, c ∈ C C C C C C, P ∈ P P P P P P n , b j ∈ B B B B B B, and Q j ∈ P n .Note that Q j are concrete Pauli elements, while P is abstract.Further, both P and Q j can have a prefactor, i.e., are not necessarily bare Paulis.Here, the integer counter r records how many concrete summands were abstracted.Specifically, r ⋆ x is defined as r i=1 x.Overall, we interpret ρ as: relying on the previously discussed interpretations of C, P n , and B.

Abstract Transformers
We now formalize the abstract transformers used by Abstraqt to simulate quantum circuits.The soundness of all transformers is straightforward, except for the trace transformer ( §5.4) which we discuss in App. A.
Initialization.We start from initial state ⊗ n i=1 |0⟩, which corresponds to density matrix as established in [2, Sec.III].We note that we can prepare other starting states by applying appropriate gates to the starting state ⊗ n i=1 |0⟩.

Gate Application
Analogously to the concrete case discussed in §2, applying a unitary gate U to ρ yields: for 21) still holds, but we cannot represent the resulting matrices efficiently.In this case, again analogously to §2, we instead decompose the offending gate as U = p d p R p , with R p ∈ P n and obtain for c pq = d p cd * q , P ′ pq = R p P R q , and b jq = b j + Q j ⋄ R q .Overall, we can evaluate Eqs. ( 21)-( 22) by relying on the abstract transformers from §4.

Compression.
To prevent an exponential blow-up of the number of summands and to adhere to the abstract domain of ρ which does not include a sum, we compress all summands to a single one.Two summands can be joined as follows: , and P = P 1 ⊔ P 2 .The key observation here is that the concrete Q j are independent of the summand, and thus need not be joined.
We note that we could also only merge some summands and leave the others precise-investigating the effect of more flexible merging strategies could be interesting future research.with probability tr(ρ − ).We describe how to compute ρ + .Computing ρ − works analogously by using −R instead of R.

Measurement
In the following, we will consider a concrete state ρ as defined in §2 and an abstract state ρ as defined in Eq. ( 19): Concrete simulation of measurement distinguishes two cases: either (i) R commutes with all Q j or (ii) R anti-commutes with at least one Q j .Note that as the Q j are concrete in an abstract state ρ, those two cases translate directly to the abstract setting.We now describe both cases for concrete and abstract simulation.

Background: Concrete Case (i).
In this case, we assume R commutes with all Q j .Focusing on a single summand ρ i of ρ, measurement maps it to (see [2]): Let us first introduce the notation {(−1) bij Q j } ⇝ R, denoting that R can be written as a product of selected Pauli elements from {(−1) bij Q j }.Symmetrically, we write is null.Further, using that R 2 = I, we get from Eq. ( 24) that if P i commutes with R, ρ i,+ is equal to ρ i , otherwise, P i anti-commutes with R and ρ i,+ is null.Putting it all together, we finally get: Abstract Case (i).Let us first define ⇝ u and ̸ ⇝ u for a concrete R, concrete Q j and abstract b j .We say and ̸ ⇝ u are under-approximations, and there can exist some R and {(−1) bj Q j } such that neither apply.Using those two abstract relations, we get the abstract transformer for ρ + : We can evaluate Eq. ( 26) by relying on the abstract transformers from Tab. 3 and by evaluating ⇝ u as discussed shortly.

Background: Concrete Case (ii).
We now suppose R anti-commutes with at least one Q j .In this case, we can rewrite ρ such that R anti-commutes with Q 1 , and commutes with all other Q j .Specifically, we can select any Q j * which anti-commutes with R, swap b ij * and Q j * with b i1 and Q 1 , and replace all other Q j anti-commuting with R by Q 1 Q j (and analogously b ij by b ij + b i1 ), which leaves ρ invariant (see [2]).Assuming ρ is the result after this rewrite, we have: where Overall, after rewriting ρ as above, Eq. ( 27) replaces c i by 1 2 c i , P i by P ′ i , b i1 by 0, and Q 1 by R.

Abstract Case (ii).
In the abstract case, we first apply the same rewrite as in the concrete case, where we pick j * as the first j for which Q j anti-commutes with R. 6 Then, directly abstracting Eq. ( 27) yields: where Here, we replace c by 1 2 c, P by P ′ , b 1 by {0}, and Q 1 by R. When defining P ′ , we follow the two cases from Eq. ( 27) when our abstraction is precise enough to indicate which case we should choose, or join the results of both cases otherwise.Again, we can evaluate Eq. ( 28) by relying on the abstract transformers from Tab. 3.
Joining Both Measurement Results.For measurements occurring within a quantum circuit, stabilizer simulation generally requires randomly selecting either ρ + or ρ − with probability tr(ρ + ) and tr(ρ − ), respectively, and then continues only with the selected state.In contrast, Abstraqt can join both measurement outcomes into a single abstract state ρ + ⊔ ρ − , as the Q j are the same in both.This allows us to pursue both measurement outcomes simultaneously, as we demonstrate in §7.

Efficiently computing ⇝
To simulate the result of a measurement, we introduced the new operator {(−1) bj Q j } ⇝ R, denoting that some Pauli R can be written as a product of {(−1) bj Q j }.We now show how to compute ⇝ efficiently.
Background: Concrete case.We first note that {(−1) bj Q j } ⇝ R holds if and only if there exist some x ∈ B n such that: Further, this solution x would satisfy: Eq. ( 30) has a solution if and only if R commutes with all the Q j , in which case this solution x is unique (see [2]).Hence, to check if {(−1) bj Q j } ⇝ R, we can first verify whether R ⋄ Q j = 0 for all j, and if so, check if the unique x satisfying Eq. ( 30) also satisfies Eq. ( 29).
Background: Finding x for Eq.(30).To compute this solution x, the stabilizer simulation relies critically on an isomorphism g between Pauli matrices {I, X, Y, Z} and B 2 .Specifically, g maps I to ( 0 0 ), X to ( 1 0 ), Y to ( 1 1 ), and Z to ( 0 1 ).Further, g extends naturally to bare Pauli elements R ∈ b(P n ) and tuples n by: where g(R) ∈ B 2n×1 and g(Q) ∈ B 2n×n .We can naturally extend g to P n , by defining g(R) = g(b(R)).This isomorphism g is designed so that the product of bare Pauli elements ignoring prefactors corresponds to a component-wise addition of encodings: Using Eq. ( 31), we can obtain solution candidates x for Eq. ( 30) by solving a system of linear equations using Gaussian elimination modulo 2: (32) Because in our case, g(Q) is over-determined and has full rank, Eq. ( 32) either has no solution, or a unique solution x.
Background: Checking prefactors.Once we have found the unique x (if it exists) satisfying Eq. (30) as described above, we need to check if it also satisfies Eq. ( 29).It is enough to check if the prefactors match: or equivalently: where the subtraction and sum operations are over Z 4 .
Putting it all together, we can define F : where x is the unique value such that g(R) = g(Q)x and indicates there is no such x.We then have that F for abstract b j .For abstract values b j , we define F : Following the same reasoning as above, we have F for abstract b j and R. To compute the trace of a state (see §5.4), we further extend Eq. ( 33) to abstract b j and abstract R, and define F : Here, evaluating Eq. ( 35) requires evaluating Q b j for an abstract boolean b, which we define naturally as Further, Eq. ( 36) requires over-approximating all x which satisfy the linear equation g(R) = g(Q)x.
Here, we naturally extend g to abstract Paulis by joining their images.For instance, we have that . We then view g(R) = g(Q)x as a system of linear equations b = Ax, where the left-hand side consists of abstract booleans b ∈ B B B B B B 2n .We then drop all equations in this equation system where the left-hand side is {0, 1}, as they do not constrain the solution space.This updated system is fully concrete, hence we can solve it using Gaussian elimination.We get either no solution, or a solution space y + p k=1 λ k u k , where y is a possible solution and u 1 , ..., u p is a possibly empty basis of the null solution space.In the case of no solution, x is not needed in Eq. ( 35).Otherwise, we can compute x j as {y j + m k=1 λ k u k,j | λ k ∈ B}.

Trace
Recall that the probability of obtaining state ρ + when measuring ρ is tr (ρ + ).We now describe how to compute this trace using F defined above.
Background: Concrete Trace.Following [2], we compute the trace of a density matrix ρ by: where we define i := 0. Because the trace of a density matrix is always real, Re(•) is redundant, but will be convenient to avoid complex traces in our abstraction.
Abstract Trace.For an abstract state ρ, we define: where we use F(•) as defined in Eq. ( 35).

Implementation
In the following, we discuss our implementation of the abstract transformers from §4 and §5 in Abstraqt.
Language and Libraries.We implemented Abstraqt in Python 3.8, relying on Qiskit 0.40.0 [14] for handling quantum circuits, and a combination of NumPy 1.20.0[15] and Numba 0.54 [16] to handle matrix operations.

Bit Encodings. An abstract density matrix
To encode the concrete Pauli matrices Q j , we follow concrete stabilizer simulation encodings such as [17] and encode Pauli matrices P using two bits g(P ) (see §5.3).To encode abstract elements of a finite set we use bit patterns.For example, we encode , where the least significant bit (i.e. the right-most bit Further, we encode {Z, Y} as 1100 2 , where the indicator bits correspond to Z, Y, X, and I, respectively, from left to right.Hence the abstact Pauli P = ({0, 3}, {Z, Y}, {X}) would be represented as (1001 2 , 1100 2 , 0010 2 ).are small finite domains, we can implement operations in these domains using lookup tables, which avoids the need for bit manipulation tricks.While such tricks are applicable in our context (e.g., [2] uses bit manipulations to compute H (i) P H † (i) for P ∈ P n ), they are generally hard to come up with [18].In contrast, the efficiency of our lookup tables is comparable to that of bit manipulation tricks, without requiring new insights for new operations.

Implementing Transformers. The abstract transformers on abstract density matrices can be implemented using operations in
For example, to evaluate {} + {0} over B B B B B B using Eq. ( 8), we encode the first argument {} as 00 and the second argument {0} as 01.Looking up entry (00, 01) in a two-dimensional pre-computed table then yields 00, the encoding of the correct result {}.We note that we cannot implement this operation directly using a XOR instruction on encodings, as this would yield incorrect results: 00 XOR 01 = 01 ≃ {0}, which is incorrect.
Gaussian Elimination.To efficiently solve equations modulo two as discussed in §5, we implemented a custom Gaussian elimination relying on bit-packing (i.e., storing 32 boolean values in a single 32-bit integer).In the future, it would be interesting to explore if Gaussian elimination could be avoided altogether, as suggested by previous works [2,17].
Testing.To reduce the likelihood of implementation errors, we have complemented Abstraqt with extensive automated tests.We test that abstract transformers f ♯ are sound with respect to concrete functions f , that is to say that We check this inclusion for multiple selected samples of x i and x i ∈ x i (typically corner cases).This approach is highly effective at catching implementation errors, which we have found in multiple existing tools as shown in §7.

Evaluation
We now present our evaluation of Abstraqt, demonstrating that it can establish circuit properties no existing tool can establish.

Benchmarks
To evaluate Abstraqt, we generated 12 benchmark circuits, summarized and visualized in Tab. 4.

Benchmark Circuit Generation.
Each circuit operates on 62 qubits, partitioned into 31 upper qubits and 31 lower qubits.We picked the limit of 62 qubits because our baseline ESS (discussed shortly) only supports up to 63 qubits; Abstraqt is not subject to such a limitation.
Each circuit operates on initial state |0⟩ and is constructed to ensure that all lower qubits are eventually reverted to state |0⟩.We chose this invariant as it can be expressed for most of the evaluated tools, as we will discuss in §7.2.Further, as some tools can only check this for one qubit at a time, we only check if the very last qubit is reverted to |0⟩, instead of running 31 independent checks (which would artificially slow down some baselines).Note that this check is of equivalent difficulty for all lower qubits.Benchmark Details.Tab. 4 details how each benchmark circuit was generated.Most of the circuits are built from three concatenated subcircuits.First, c 1 modifies the upper qubits, then c 2 modifies the lower qubits (potentially using gates controlled by the upper qubits) and finally c 3 reverts all lower qubits to |0⟩, but in a non-trivial way.Circuit CCX+H;Cliff slightly deviates from this pattern, as it also modifies the upper qubits using gates controller by lower qubits.Further, circuits Cliff+T;H;CZ+RX and Cliff+T;H;CZ+RX' additionally apply two layers of H gates to the lower qubits.Finally, circuit MeasureGHZ applies internal measurements, as discussed below.
The majority of circuits revert the lower qubits to |0⟩ by applying c 3 , the inverse of c 2 but optimized using PyZX [19]-this obfuscates the fact that c 2 and c 3 cancel out.Four circuits, marked with a trailing prime ('), generate c 3 by optimizing the un-inverted c 2 .They still reset all lower qubits to |0⟩, but establishing this requires advanced reasoning.Specifically, RZ 2 +H;CX' flips each lower qubit an even number of times. 7Similarly, Cliff+T;CX+T' and CCX+H;CX+T' additionally modify the phase but still flip each lower qubit an even number of times.Finally, Table 4: Description of benchmark circuits, where upper = {1, . . ., 31} and lower = {32, . . ., 62}.
We note that because YP21 does not support CX(a, b) for a > b, we instead encoded such gates as H(a); H(b); CX(b, a); H(b); H(a).
Statevector.Qiskit [14] further provides a simulator based on state vectors, which we also used for completeness.
Abstraqt.In Abstraqt, we can establish that a qubit is in state |0⟩ by measuring the final abstract state ρ in basis Z (i) and checking if the probability of obtaining |1⟩ is 0.
Experimental Setup.We executed all experiments on a machine with 110 GB RAM and 56 cores at 2.6 GHz, running Ubuntu 22.04.Because some tools consumed excessive amounts of memory, we limited them to 12 GB of RAM.This was not necessary for Abstraqt, which never required more than 600 MB of RAM.We limited each tool to a single thread.

Results
Tab. 5 summarizes the results when using all tools discussed in §7.2 to establish that the last qubit in 10 randomly selected instantiations of each benchmark from Tab. 4 is in state |0⟩.Overall, it demonstrates that while Abstraqt can establish this for all benchmarks within minutes, QuiZX can only establish it for a few instances, and all other tools cannot establish it for any benchmark.Further, we found that for some circuits the established simulation tool ESS yields incorrect results.We now discuss the results of each tool in more details.
MeasureGHZ.Importantly, no baseline tool except Abstraqt can simultaneously simulate both outcomes of a measurement, without incurring an exponential blow-up.Therefore, for MeasureGHZ, we consider internal measurements as an unsupported operation in these tools.We note that we could randomly select one measurement outcome and simulate the remainder of the circuit for it, but then we can only establish that the final state is |0⟩ for a given sequence of measurement outcomes.In contrast, a single run of Abstraqt can establish that the final state is |0⟩ for all possible measurement outcomes (see also §5.2).
QuiZX.As QuiZX is the only baseline tool solving some of our benchmark instances, we provide a detailed comparison to it in Tab. 6.
Overall, QuiZX cannot consistently handle any of the benchmarks from Tab. 4, Instead, it often either times out or runs out of memory.Further, QuiZX consistently runs into an internal error Future Abstractions.We expect that for many real-world circuits, existing approaches work better than the current implementation of Abstraqt.However, as Abstraqt only abstracts the first stabilizer simulation generalized to non-Clifford gates [2, §VII-C], we believe it paves the way to also abstract more recent stabilizer simulators.For example, ESS [5] operates on so-called CH-forms which, like the generalized stabilizer simulation underlying Abstraqt, can be encoded using bits and complex numbers.Hence, it seems plausible that our ideas could be adapted to abstract ESS.QuiZX operates on ZX-diagrams consisting of graphs whose nodes are parametrized by rotation angles α.Again, a promising direction for future research is introducing abstract ZXdiagrams that support abstract rotation angles.This is particularly promising because both ESS and QuiZX scale better in number of T gates than [2, §VII-C]: with 2 n instead of 4 n .
We note however that not all concrete simulation techniques are directly amenable to abstraction.For example, when naively abstracting the Clifford simulation by Aaronson and Gottesmann, applying a measurement requires selecting an entry in an boolean matrix that definitively equals one [2, Case I in §III]-it is unclear how to generalize this to abstract boolean matrices whose entries may be {0, 1}.
Improving Abstraqt.Another promising route towards better abstractions in incrementally improving Abstraqt itself.For example, it would be interesting to consider the effect of keeping more than one abstract summand, abstracting P i or b ij using a custom relational domain (which retains information about the relationship between different values) [22], or a more precise abstraction for complex numbers by taking into account that restricted gate sets such as Clifford+T only induce matrices over finite sets of values.
Summary.Overall, we believe that all tools in Tab. 5 are valuable to analyze quantum circuits.We are hoping that addressing some limitations of the considered baselines (e.g., fixing bugs in QuiZX and ESS) and cross-pollinating ideas (e.g., extending QuiZX by abstract interpretation) will allow the community to benefit from the fundamentally different mathematical foundations of all tools.

Related Work
Here, we discuss works related to the goal and methods of Abstraqt.
Quantum Abstract Interpretation.Some existing works have investigated abstract interpretation for simulating quantum circuits [21,23,24].As [21] is not specialized for Clifford circuits, it is very imprecise on the circuits investigated in §7: it cannot derive that the lower qubits are |0⟩ for any of them.While [23,24] are inspired by stabilizer simulation, they only focus on determining if certain qubits are entangled or not, whereas Abstraqt can extract more precise information about the state.Further, both tools are inherently imprecise on non-Clifford gatesin contrast, a straight-forward extension of Abstraqt can treat some non-Clifford gates precisely at the exponential cost of not merging summands.
Stabilizer Simulation.The Gottesman-Knill theorem [1] established that stabilizers can be used to efficiently simulate Clifford circuits.Stim [17] is a recent implementation of such a simulator, which only supports Clifford gates and Pauli measurements.
Stabilizer simulation was extended to allow for non-Clifford gates at an exponential cost, while still allowing efficient simulation of Clifford gates [2, §VII-C].Various works build upon this insight, handling Clifford gates efficiently but suffering from an exponential blow-up on non-Clifford gates [3,4,5,6,7].In our evaluation, we demonstrate that Abstraqt extends the reach of state-of-the-art stabilizer simulation by comparing to two tools from this category, ESS [5] (chosen because it is implemented in the popular Qiskit library) and QuiZX [4] (chosen because it is a recent tool reporting favorable runtimes).
Verifying Quantum Programs.Another approach to establishing circuit properties is end-toend formal program verification, as developed in [25] for instance.However, this approach often requires new insights for each program it is applied to.Even though recent works have greatly improved verification automation, proving even the simplest programs still requires a significant time investment [26], whereas our approach can analyze it without any human time investment.
The work [27] automatically generates rich invariants, but is exponential in the number of qubits, limiting its use to small circuits.Finally, [20] can automatically verify the equivalence of two given circuits, but times out on the benchmarks considered in §7.

Conclusion
In this work, we have demonstrated that combining abstract interpretation with stabilizer simulation allows to establish circuit properties that are intractable otherwise.
Our key idea was to over-approximate the behavior of non-Clifford gates in the generalized stabilizer simulation of Aaronson and Gottesman [2] by merging summands in the sum representation of the quantum states density matrix.Our carefully chosen abstract domain allows us to define efficient abstract transformers that approximate each of the concrete stabilizer simulation functions, including measurement.

K
p=1 d p R p , where d p ∈ C and R p ∈ b(P n ) are bare Pauli elements, which have a prefactor of i 0 = 1.Then,

Figure 1 :
Figure 1: Overview of Abstraqt, where we define c and c 1 -c 4 in §3.
Density Matrices.The concrete and abstract domains introduced previously allow us to represent an abstract density matrix ρ ∈ D D D D D D as follows: We now describe how to perform Pauli measurements, by extending the (concrete) stabilizer simulation to abstract density matrices.The correctness of the concrete simulation was previously established in [2, Sec.VII.C], while the correctness of the abstraction is immediate.SimulatingMeasurement.Applying a Pauli measurement in basis R ∈ b(P n ) has a probabilistic outcome and transforms ρ to ρ + = I+R 2 ρ I+R 2 with probability tr(ρ + ) or ρ − = I−R 2 ρ I−R 2

Table 1 :
Transformers for the interval abstraction.

Gate Application. First
, the circuit applies one Hadamard gate H to each qubit.

Table 2 :
Example elements on abstract domains.

Table 5 :
Success rates when running simulators on benchmarks from Tab. 4.