qsimlab.analysis
Structural analysis of circuits: magic, stabilizer rank, simulability, monitored dynamics.
Provisional (phase 2, see python/API.md §1 and §8).
These functions measure why a circuit is hard or easy to simulate, with the exact invariants the engines use, mostly without simulating it:
>>> import qsimlab as qs, qsimlab.analysis as an
>>> c = qs.Circuit(3).h(0).t(0).cx(0, 1).h(2)
>>> p = an.magic_profile(c)
>>> p.t_count, p.d, p.f
(1, 1, 1)
>>> an.branching_rank(c).rank
2
>>> round(an.state_magic(c).m2, 6) # log2(4/3): one T state's worth
0.415037
>>> an.branching_rank(c.copy().t(1)).rank # T on the partner makes it a stabilizer state again
1
magic_profile: the magic atlas (research/simulability/magic-atlas.md): active dimensiond_kof the rotation frame (the exact register size of the compressed-state engine), factored dimensionf_k, stabilizer entanglement of the Clifford skeleton across a cut and the boundE + don the true entanglement, and an affine bound on the support size.state_magic: stabilizer nullity and stabilizer 2-Rényi entropyM2of a state (all4^nPauli expectations;n ≤ 13).branching_rank: the number of stabilizer terms of the exact low-rank simulator (research/theory/theory-rank.md) after every gate, an upper bound on the stabilizer rank.simulability: the planner's features (per-engine log2 work estimates) and its explanation (ranked predicted costs per engine).gaussian: the free-fermion (matchgate) detector: is the circuit a fermionic Gaussian circuit under some Jordan–Wigner order, after undoing SWAP networks, and how far from it (docs/ENGINE_GAUSSIAN.md);gaussian_expectationsruns the Gaussian engine.monitored: exact simulation of Clifford+T circuits with mid-circuit measurements in the rotation frame:d(t), Born probabilities, and cut entropies.
dataclass
¶
MagicProfile(
num_qubits: int,
gates: int,
lowered_gates: int,
two_qubit_gates: int,
toffolis: int,
rotations: int,
t_count: int,
d: int,
f: int,
log2_work: float,
log2_work_factored: float,
cut: int,
e_stab_max: int,
e_bound_max: int,
support_log2: Optional[int],
seconds: float,
d_profile: ndarray,
f_profile: ndarray,
rotation_gate: ndarray,
checkpoints: List[Tuple[int, int, int, int, int, int]],
)
Result of magic_profile (the magic atlas of a unitary circuit on |0^n⟩).
t_count: rotations by odd multiples of π/4;rotations: every non-Clifford rotation after lowering to Clifford + Z rotations (Toffoli = 7 T) and merging half-π multiples into S;d: final active dimension (dimension of the x-span of the rotation axes in the Heisenberg frame);2^dis the compressed-state register;d_profile[j]isdafter rotationj(rotation_gate[j]is its original gate index);f: largest factored component;f_profile[j]the component size touched by rotationj;log2_work=log2 Σ_j 2^{d_j}(amplitude updates of the compressed state; −1 for Clifford circuits),log2_work_factoredthe same forf;e_stab_max: max stabilizer entanglement (bits) of the Clifford skeleton acrosscut;e_bound_max: max ofmin(E + d, cut, n − cut), an upper bound on the entanglement of the true state;support_log2: affine upper bound onlog2 |supp U|0^n⟩|;checkpoints:(gate, rotations, t_count, d, f, e_stab)tuples at evenly spaced gates.
dataclass
¶
nullity: stabilizer nullity ν = n − log2 |{P : |⟨P⟩| = 1}| (Beverland et al.
2020; 0 iff stabilizer state). m2: stabilizer 2-Rényi entropy
M2 = −log2(Σ_P ⟨P⟩⁴ / 2^n) (Leone, Oliviero, Hamma 2022), in bits.
dataclass
¶
BranchingRank(
rank: int,
max_rank: int,
overflow: bool,
branch_events: int,
clifford_events: int,
diagonal_events: int,
merges: int,
pair_merges: int,
cancellations: int,
trace: ndarray,
state: Optional[ndarray] = None,
)
Result of branching_rank.
trace[k] is the number of stabilizer terms after gate k (gates of the
unitary part, in order), rank the final value, max_rank its maximum;
overflow is True if max_terms stopped the run. Events: per term and
non-Clifford gate, branch_events (the term split in two),
clifford_events (updated in place by a Clifford, Theorem R3),
diagonal_events; merges/pair_merges/cancellations: terms
combined. state: the dense state Σ_j c_j|φ_j⟩ (state=True, n ≤ 20).
dataclass
¶
Simulability(
features: Dict[str, Any],
log2_costs: List[Tuple[str, float]],
explanation: Explanation,
)
Result of simulability.
features: the raw feature dict (research/simulability/simulability.md):n, gates, g2, g3, depth2, t_count, rotations, d(active dimension),chi_bits(bound on log2 of the MPS bond),hsf_k(HSF cut bits),sup(affine bound on log2 of the support), and one*_llog2 work estimate per engine;log2_costs:{engine: log2 work}from those estimates, cheapest first;explanation: the planner'sExplanationforrequest(ranked predicted seconds per engine, from fitted cost models).
dataclass
¶
GaussianReport(
num_qubits: int,
exact: bool,
free: bool,
blocks: int,
blocks_2q: int,
gaussian_blocks: int,
gaussian_fraction: float,
max_residual: float,
non_gaussian: int,
nonadjacent: int,
interaction_total: float,
interaction_max: float,
swaps_relabelled: int,
ordering: str,
number_conserving: bool,
seconds: float,
interactions: List[
Tuple[int, Tuple[int, int], Tuple[int, int], float]
],
order: List[int],
paths: List[List[int]],
mode_of_qubit: List[int],
)
Result of gaussian (docs/ENGINE_GAUSSIAN.md §2).
exact: every fused block is Gaussian in the chosen Jordan–Wigner order (totol): the Gaussian engine applies;free: Gaussian up to diagonal interaction phases;gaussian_fraction: Gaussian blocks / all blocks;max_residual: largest distance of a non-interaction block from the Gaussian set (1 for a matchgate on non-adjacent modes or a three-qubit gate);interactions:(block, (wire_a, wire_b), (mode_a, mode_b), g)for every diagonal blockexp(i g n_a n_b)withg != 0;interaction_total=Σ|g|,interaction_max;ordering("identity","paths","greedy_cover"),order[k]= wire (initial qubit) on modek,paths(chains of wires),mode_of_qubit[q]= mode held by qubitqat the end;swaps_relabelled: SWAP gates and SWAP-equivalent blocks turned into wire renamings;number_conserving: every Gaussian block conserves the particle number.
dataclass
¶
GaussianExpectations(
z: ndarray, zz: ndarray, report: GaussianReport
)
Result of gaussian_expectations: z[q] = ⟨Z_q⟩,
zz[i] = ⟨Z_a Z_b⟩ for pairs[i] = (a, b), and the detector's
GaussianReport.
dataclass
¶
MonitoredResult(
engine: str,
components: List[Tuple[int, int, str]],
seed: int,
precision: str,
wall_time: float,
explanation: Optional[Explanation] = None,
d: Optional[ndarray] = None,
final_d: int = 0,
max_d: int = 0,
outcomes: Optional[ndarray] = None,
qubits: Optional[ndarray] = None,
probabilities: Optional[ndarray] = None,
kinds: Optional[ndarray] = None,
entropies: List[
Tuple[
int, List[Tuple[float, float, Optional[float]]]
]
] = list(),
stats: Dict[str, int] = dict(),
state: Optional[ndarray] = None,
)
Bases: Result
Result of monitored (a qsimlab.sim.Result).
d: active dimension after every op (uint32[len(circuit)]); the state isC (|φ⟩ ⊗ |0⟩)with|φ⟩ondvirtual qubits, so2^damplitudes are stored;final_d,max_d;outcomes/qubits/probabilities: one entry per measurement (program order; resets are not recorded): the outcome, its qubit and the Born probability of the observed outcome;kinds: 0 = random in the frame (probability 1/2, no amplitude work), 1 = measured on the register (ddrops by one), 2 = determined;entropies:[(op_index, [(lower, upper, s2), ...per cut])]: bounds on every Rényi entropy of each cut (bits) and the exact Rényi-2 entropys2when the register is tracked and small enough;stats: T gates (activating / in-register), measurement kinds,element_ops(amplitude updates);state: the final state (state=True), up to a global phase.
magic_profile(
circuit: Circuit,
*,
checkpoints: int = 64,
cut: Optional[int] = None,
entanglement: bool = True,
support: bool = True,
threads: Optional[int] = None,
) -> MagicProfile
The magic atlas of circuit in one O(gates · n) pass (no simulation).
Terminal measurements are ignored; mid-circuit measurements, resets and
noise raise UnsupportedOperationError (use monitored).
cut: the skeleton entanglement is across qubits [0, cut) vs the rest
(default n // 2).
state_magic(
x: Union[Circuit, ndarray],
*,
threads: Optional[int] = None,
) -> StateMagic
Stabilizer nullity and stabilizer Rényi entropy of a state (n ≤ 13).
x is a Circuit (its unitary part is simulated from
|0^n⟩) or a normalised state vector of length 2^n. Cost O(4^n n).
stabilizer_nullity(x: Union[Circuit, ndarray]) -> float
StateMagic.nullity of state_magic.
stabilizer_renyi_entropy(
x: Union[Circuit, ndarray],
) -> float
StateMagic.m2 of state_magic.
branching_rank(
circuit: Circuit,
*,
max_terms: int = 1 << 16,
pair_merge: int = 6,
state: bool = False,
threads: Optional[int] = None,
) -> BranchingRank
Run the exact branching-rank (low stabilizer rank) simulator on circuit.
Every non-Clifford gate is a projector gate; a term branches only when the
gate cannot be absorbed as a Clifford on it, and equal rays are merged, so
the term count is an upper bound on the stabilizer rank (exact for many
structured circuits). pair_merge: merge pairs of terms differing in ≤ s
stabilizer generators (0 disables it).
simulability(
circuit: Circuit,
request: Optional[Request] = None,
*,
hsf: bool = True,
budget: Union[Budget, int, str, None] = None,
threads: Optional[int] = None,
) -> Simulability
Simulability features of circuit plus the planner's explanation for request.
The features are computed on the unitary part (terminal measurements
dropped; other non-unitary ops raise UnsupportedOperationError); request defaults to samples(1024). hsf=False skips the
Kernighan–Lin partition, the slowest feature on large circuits.
gaussian(
circuit: Circuit,
*,
tol: float = 1e-10,
relabel_swaps: bool = True,
reorder: bool = True,
threads: Optional[int] = None,
) -> GaussianReport
Free-fermion detector: fuses the gates into blocks and tests each for the matchgate property in a Jordan–Wigner order it searches for.
Takes the unitary part (terminal measurements dropped). relabel_swaps
treats SWAP gates as wire renamings; reorder allows an order other
than the qubit order.
gaussian_expectations(
circuit: Circuit,
pairs: Sequence[Tuple[int, int]] = (),
*,
drop_interactions: bool = False,
tol: float = 1e-10,
threads: Optional[int] = None,
) -> GaussianExpectations
⟨Z_q⟩ of every qubit and ⟨Z_a Z_b⟩ (Wick's theorem) on the
Gaussian engine, in O(gates · n) plus O(1) per pair.
The circuit must be exactly Gaussian (UnsupportedOperationError
otherwise). drop_interactions=True instead sets every diagonal
interaction phase exp(i g n_a n_b) to zero, which gives the free-fermion
part of the circuit. That is an approximation unless report.exact.
monitored(
circuit: Circuit,
*,
seed: Optional[int] = None,
exact: bool = True,
max_d: int = 24,
cuts: Optional[Sequence[Sequence[int]]] = None,
entropy_every: int = 0,
max_cost_log2: int = 24,
state: bool = False,
threads: Optional[int] = None,
) -> MonitoredResult
Simulate a Clifford+T circuit with mid-circuit measurements in the rotation frame.
Clifford gates update a tableau; a T/Rz gate either activates one
virtual qubit (d → d + 1) or rotates the active register in place; a
Z measurement is random in the frame (1/2), determined, or Born-sampled on
the register (d → d − 1). Any gate is accepted (lowered to Clifford +
Z rotations); measurements, resets, c_if and Pauli noise channels are
simulated per shot. One call is one trajectory.
exact=Falsetracks only the tableau:d(t)is exact for every outcome sequence (it does not depend on the outcomes), register outcomes are drawn 50/50 and only entropy bounds are reported; works for anyn.max_d: refuse (ResourceLimitError) beyond2^max_damplitudes (≤ 34).cuts: regions (lists of qubits) whose entropy is reported, at the end and everyentropy_everyops (default: the firstn // 2qubits).