Nonabelian lattice-gauge simulation

One layer of a lattice gauge theory simulation: multiply group elements around hundreds of small loops and add up their traces. The discovered circuit cuts the cost of each multiplication from 146 Toffolis to 35 and lets neighboring loops share their work.

In brief

The baseline (Osorio Perez et al., 2025) treats each loop on its own. It multiplies the loop's links, looks up the trace in a table, and runs the loop backward.

The discovered circuit chooses coordinates for the one operation the loops need: multiplying a running product by a link. Neighboring loops share runs of links, so the circuit edits the running product from one loop into the next. A small logic circuit reads the trace directly, which costs less than the table.

TQ = 2·10⁸TQ = 5·10⁸TQ = 10⁹TQ = 2·10⁹1,5002,0002,5002·10⁵5·10⁵10⁶peak qubits Q →Toffoli count T →Discovered circuitOsorio Perez et al. (2025)Live qubits as Toffolis execute1.3M170k1,927
Links and accumulator hold most qubits in both circuits, so the difference lies along the Toffoli axis.
1

Coordinates for one operation

XO_CANDS

Each of the 192 links holds an element of Σ(72×3) as a nine-qubit code for g = ωpCqErV2s+tXu: three trits and three bits. The 192 plaquettes and 384 rectangles need one operation, a right multiplication of a running product W by a link.

First it replaces the generator X by another element Xo of the same coset. The search keeps a group element only if it squares to V2, conjugates V to V−1, and differs from X by an element of ⟨ω, C, E⟩. The build asserts that exactly one element passes. The V and Xo factors then close among themselves, and recoding a link from X to Xo changes only p, q, and r.

Links keep their factor order. The running product uses the opposite placement, W = ωPVMXoυCQER, so the accumulated and the incoming C, E factors meet at the multiplication boundary. W stores the sheared phase P̃ = P + QR mod 3, which folds a recurring correction into the coordinate itself.

# C. The internal X representative and the mirror algebra. The chain works
# with an internal generator Xo (in place of the contract Xc) chosen so
# that X^2 = V^2, X V X^-1 = V^-1, and
# Xc*Xo^-1 lies in Delta(27) = <w,C,E> -- then the
# U-flag composes linearly with the M-counter and the encode transform
# (contract code -> internal code) is cheap. No tie-break is needed: the
# algebraic filter admits EXACTLY ONE candidate (asserted; a group-model
# change that widens the set fails the build).
# ---------------------------------------------------------------------------
XO_CANDS = []
for g in sorted(G):
if g in sig36:
continue
if mpow(g, 2) != V2m:
continue
if mmul(mmul(g, Vg), mdag(g)) != Vinvm:
continue
if mmul(Xc, mdag(g)) in delta:
XO_CANDS.append(g)
assert len(XO_CANDS) == 1, "internal X representative is unique"
Xo = XO_CANDS[0]
2

A walk through all 576 loops

KIND, _decode_prog

The circuit keeps W from one loop and edits it into the next. The real trace does not change under cyclic rotation or inversion of a loop, so each loop has several equivalent representatives. Only the right end of W is editable: multiplication by inverse links cancels an old suffix, and new links append. Inverting W in place brings its other end to the right.

The schedule is a token sequence (KIND). The tokens f and i multiply by a link or by its inverse, l and u open or close a segment at the identity, v inverts W, and P and R fire a plaquette or a rectangle. An exact dynamic program in gen_constants.py chooses the loop order, the representatives, and the restart points. Its optimum within that move set uses 3,136 multiplications, 480 inversions, and 64 anchors.

The optimum is row-periodic. constants.json therefore stores an opening row, a row program, a plane-crossing program, and a closing tail of 260, 266, 266, and 6 tokens. The two repeated programs address links relative to the row they process.

# F. The loop schedule. Conventions pinned by PARAMS-DERIVATION.md
# sections 3 and 7 / reference.py: site s = x + 4y + 16z, link l = 3s + mu,
# inverses on return links.
# ---------------------------------------------------------------------------
# The global chain-walk schedule. Ops: 'f'/'i' push the link / its
# inverse (35-Toffoli multiplies); 'li'/'ul' open/close at the identity
# with a 1-Toffoli inverse-load / its reverse; 'v' inverts the register in
# place (21 Toffoli) so subsequent right-multiplies edit the word's other
# end; 'XP'/'XR' fire a plaquette/rectangle trace-add. Exact DP optimum
# of the net-flip-relaxed morph model (odd inversion counts per morph
# allowed): 3,136 multiplies + 480 inversions + 64 anchors.
#
# The DP optimum is row-periodic: every y-junction row of the x-fastest
# raster is one lattice translate of a single row program, every
# plane-crossing row another, so constants.json stores THREE row
# programs (open_row, then y_row x12 and z_row x3 by translation) plus
# the closing tail, and this build repeats them through a slab conveyor
# (section L) instead of unrolling 4,256 tokens. y_row/z_row link ids
# are FRAME-LOCAL: link 3*(x + 4*dy + 16*dz) + mu addresses the link dy
# rows / dz planes above the row being processed (mod 4); open_row and
# tail are absolute. gen_constants.py runs the search and writes the
# file.
_CONSTS = json.loads((pathlib.Path(__file__).resolve().parent / "constants.json").read_bytes())
def _decode_prog(blobs):
KIND = {"f": "f", "i": "i", "l": "li", "u": "ul", "v": "inv", "P": "XP", "R": "XR"}
out = []
for tok in blobs.split():
kd = KIND[tok[0]]
out.append((kd, int(tok[1:])) if tok[0] in "filu" else (kd,))
return out
OPEN_ROW = _decode_prog(" ".join(_CONSTS["open_row"]))
Y_ROW = _decode_prog(" ".join(_CONSTS["y_row"]))
Z_ROW = _decode_prog(" ".join(_CONSTS["z_row"]))
TAIL = _decode_prog(" ".join(_CONSTS["tail"]))
3

Links, accumulator, and 80 working wires

links, kreg, BANK, CARRY

Only the 1,779 I/O wires persist: 192 links of nine wires and the 51-wire accumulator k. W takes 12 wires. Its trits P̃, Q, and R are one-hot, three wires with exactly one set, and three bits hold M = 2S + T and υ.

On one-hot rails, adding one modulo three permutes three wires, and negation swaps two. Only one running product exists, so only W pays for the wider form. The same form for every link would add 576 qubits, so the links stay compact.

A 38-wire bank holds scratch wires, the 3×3 grid J with cells [Q = q][R = r], row rails for (M, υ), and six trace-class wires. Thirty further wires serve only as adder carries. All 80 working wires are live at each plaquette addition, which sets the peak of 1,859; rectangle additions reach 1,858.

# H. Wire map. Only the 1,779 I/O wires (192 links x 9, then k) are
# persistent. The 80-wire working set lives and dies with its data: the
# nine hot utility wires span the whole program and the 12-wire chain
# state spans the loops stage (explicit stage-boundary lifecycles); the
# other 59 wires -- the J grid, the extraction scaffolding and class
# basis, and the 30 dedicated adder carries -- are allocated when their gadget
# builds them and freed when it (or its measurement mirror) tears them
# down. At every plaquette trace-add ripple all 80 are simultaneously
# live, so the skyline peaks at exactly 1,859 (rectangle adds top out at
# 1,858: LO_R = 1 drops one carry) while breathing everywhere else.
# ---------------------------------------------------------------------------
NIO = 1779
links = [[9 * l + i for i in range(9)] for l in range(192)] # wires 0..1727
kreg = list(range(1728, 1779)) # wires 1728..1778
accP = [1779, 1780, 1781] # one-hot rails of the primed phase trit P'
Qm = [1782, 1783, 1784] # one-hot rails of the C-exponent trit Q
Rm = [1785, 1786, 1787] # one-hot rails of the E-exponent trit R
accS, accT, accU = 1788, 1789, 1790 # V-counter bits M = 2S+T; X flag U
BANK = list(range(1791, 1829)) # the 38-wire pooled ancilla bank
candB = BANK[0:4]
xorB0 = BANK[4]
ctlB = BANK[5:7]
J = BANK[7:16] # 3x3 one-hot product grid J[3q+r] = [Q=q][R=r]
a_st = BANK[16]
r10 = BANK[17]
r01 = BANK[18]
r00 = BANK[19] # M-decoder
rMU1 = BANK[20:24]
rMU0 = BANK[24:28] # (M,U) one-hot row rails
nuq = BANK[28]
sgq = BANK[29] # unit-class / negative-sign flags
scr0 = BANK[30]
t1q = BANK[31]
f2q = BANK[32]
t2q = BANK[33]
fm1q = BANK[34]
cen2q = BANK[35]
f6q = BANK[36] # trace-class basis wires
addd = BANK[37] # the full-width adder's addend wire
ecq = ctlB[0] # e-copy flag lives on ctlB[0]: free during extraction
CARRY = list(range(1829, 1859)) # 30 dedicated full-width adder carries
4

One fired loop

morph_ops, row_chunks

morph_ops lowers the edit between two fires token by token. A multiplication by a link, forward or by its inverse, costs 35 Toffolis, and neither inverts the link register. An anchor opens a segment for one Toffoli and its mirror closes one. An inversion of W costs 21.

row_chunks cuts a row program at each fire. A chunk is the edit, the trace readout, the weighted addition into k, and the mirrored readout that erases the flags. Each row lowers to 36 chunks, one per fired loop.

XF_UNIT = xform_ops(links[0]) + qr_corr_ops(links[0])
XF_FWD = expand(XF_UNIT)
XF_REV = expand(XF_UNIT, rev=True)
def morph_ops(toks):
"""Expand a run of walk tokens (no fires) to flat low-level ops."""
L = []
for op in toks:
k = op[0]
if k == "f":
L += expand(mul_fwd_ops(links[op[1]]))
elif k == "i":
L += expand(mul_inv_ops(links[op[1]]))
elif k == "li":
L += expand(invload_ops(links[op[1]]))
elif k == "ul":
L += expand(invload_ops(links[op[1]]), rev=True)
elif k == "inv":
L += expand(INVERT_M)
else:
raise AssertionError(op)
return L
def row_chunks(form):
"""Lower one row program to its 36 (morph + trace-add) chunks."""
chunks, cur = [], []
for op in form:
if op[0] not in ("XP", "XR"):
cur.append(op)
continue
mats, lo = (MATS_P, LO_P) if op[0] == "XP" else (MATS_R, LO_R)
ops_ = morph_ops(cur)
ops_ += expand(EXT_OPS)
ops_ += expand(kadd_ops(mats, lo=lo, hi=WBITS, pool=KPOOL))
ops_ += expand(EXT_OPS, rev=True)
chunks.append(ops_)
cur = []
assert not cur
return chunks
OPEN_CHUNKS = row_chunks(OPEN_ROW)
Y_CHUNKS = row_chunks(Y_ROW)
Z_CHUNKS = row_chunks(Z_ROW)
TAIL_OPS = morph_ops(TAIL)
assert max(len(c) for c in OPEN_CHUNKS + Y_CHUNKS + Z_CHUNKS) < 0xFFFF
5

One multiplication, 35 Toffolis

mul_fwd_ops

A multiplication W ← Wg reads the nine bits of the stored link code. A fused rotation cluster adds p and the cross terms Qr + 2Rq to P̃. Each cross term is a trit product of three Karatsuba product cells, and push merges adjacent controlled swaps of the same rail pair across terms. The cluster costs 13 Toffolis.

Two rail rotations add q to Q and r to R, three Toffolis each. Fixed stages then move the sum through V and Xo. Under t, a swap of the Q and R rails costs two Toffolis, because rail 0 follows from the other two. Negations of Q under s and of R under s ⊕ t cost one each. Under u, ustage_ops builds J and applies the Xo permutation as six controlled cell swaps, ten Toffolis in all. The phase correction that this stage would need vanishes for Xo; the build asserts it.

The counter update for M and υ adds two Toffolis, for 35 in all; a general normal-form multiplication costs 146. mul_inv_ops multiplies by the inverse link at the same cost.

def core3_ops(rails, c1, c2, xw):
"""One-hot trit rotate by (c1 + 2*c2) mod 3 (binary control wires),
3 Toffoli, using xw as the XOR pivot.
"""
return [
("cswap", c1, rails[1], rails[2]),
("cx", c1, xw),
("cx", c2, xw),
("cswap", xw, rails[0], rails[1]),
("cx", c1, xw),
("cx", c2, xw),
("cswap", c2, rails[1], rails[2]),
]
⋯ 19 lines
def prot_cluster(fams):
"""Fused one-hot accP rotation: total amount = sum of family amounts.
fam ('d', c1w, c2w): amount = c1w + 2*c2w (binary wires).
fam ('k', acc, bits, k): amount = k*A*j (A the one-hot trit on acc
rails, j the binary-coded link trit on bits) -- the CO_Q/CO_R cross
terms. Each 'k' family costs 3 Karatsuba product cells; adjacent
same-pair controlled swaps are XOR-merged across families (the fusion
that makes the cluster cheaper than separate rotations).
"""
pool = [candB[0], candB[1], candB[2], candB[3], ctlB[0], ctlB[1]]
L = []
U = []
toks = []
pi = [0]
def push(kind, ws):
if toks and toks[-1][0] == kind:
toks[-1] = (kind, toks[-1][1] ^ set(ws))
else:
toks.append((kind, set(ws)))
for fam in fams:
if fam[0] == "d":
push("a", {fam[1]})
push("b", {fam[1], fam[2]})
push("a", {fam[2]})
else:
acc, bits, k = fam[1], fam[2], fam[3] % 3
assert k in (1, 2)
w1 = pool[pi[0]]
w2 = pool[pi[0] + 1]
w3 = pool[pi[0] + 2]
pi[0] += 3
pr3 = (((acc[1], acc[2]), (bits[0], bits[1]), 0),)
L += [
("andc", w1, acc[1], bits[0]),
("andc", w2, acc[2], bits[1]),
("bx", w3, pr3, (), 0),
]
U = [
("mx", w3, pr3, (), 0),
("andr", w2, acc[2], bits[1]),
("andr", w1, acc[1], bits[0]),
] + U
d = {w1, w2}
t = {w3}
if k == 1:
push("a", d)
push("b", t)
push("a", d ^ t)
else:
push("a", d ^ t)
push("b", t)
push("a", d)
for kind, ws in toks:
wl = sorted(ws)
if not wl:
continue
pr = (accP[1], accP[2]) if kind == "a" else (accP[0], accP[1])
c0 = wl[0]
for w in wl[1:]:
L.append(("cx", w, c0))
L.append(("cswap", c0, pr[0], pr[1]))
for w in wl[1:]:
L.append(("cx", w, c0))
return L + U
⋯ 75 lines
def qr_swap_ops(tctl):
"""Swap the Q<->R rails under control in 2 Toffoli via the one-hot invariant
rail0 = 1 ^ rail1 ^ rail2 (rail0 derived, fixed up with free CX).
"""
pre = [("cx", Qm[1], Qm[0]), ("cx", Qm[2], Qm[0]), ("cx", Rm[1], Rm[0]), ("cx", Rm[2], Rm[0])]
return pre + [("cswap", tctl, Qm[1], Rm[1]), ("cswap", tctl, Qm[2], Rm[2])] + list(pre)
⋯ 2 lines
def ustage_ops(uctl, swaps):
"""U-conditioned (Q,R) permutation: build J, fold J onto the rails,
apply the 6 controlled J-cell swaps of the X-conjugation permutation,
unfold and unbuild. 10 Toffoli: the 4-Toffoli J build + the 6 swaps
(the measurement-based J teardown is free). The xi-phase correction
that would sit between fold and swaps VANISHES for the selected Xo
(both phase groups empty -- asserted at derivation), which is what
keeps this stage free of an accP rotation.
"""
L = J_build_ops()
L += J_fold_ops()
for i, jj in swaps:
L.append(("cswap", uctl, J[i], J[jj]))
L += J_fold_ops()
L += J_unbuild_ops()
return L
⋯ 2 lines
def stages_fwd_ops(tctl, sctl, uctl):
"""Emit the V/X part of a forward multiply: t-stage (derived-rail Q<->R
swap), s-stage (double negation), u-stage (X permutation).
"""
L = qr_swap_ops(tctl)
L.append(("cswap", sctl, Qm[1], Qm[2]))
L += [
("cx", tctl, xorB0),
("cx", sctl, xorB0),
("cswap", xorB0, Rm[1], Rm[2]),
("cx", tctl, xorB0),
("cx", sctl, xorB0),
]
L += ustage_ops(uctl, SWAPS_Xi)
return L
⋯ 16 lines
def mul_fwd_ops(g):
"""Chain push: state *= link (stored st2 code on wires g[0..8]).
35 Toffoli: fused phase cluster (13) + two rail rotations (3 each) +
V/X stages (14) + the M/U counter update (2).
"""
p0, p1, q0, q1, r0, r1, sB, tB, uB = g
fams = [("d", p0, p1)]
if CO_Q:
fams.append(("k", Qm, (r0, r1), CO_Q))
if CO_R:
fams.append(("k", Rm, (q0, q1), CO_R))
L = prot_cluster(fams)
L += core3_ops(Qm, q0, q1, xorB0)
L += core3_ops(Rm, r0, r1, xorB0)
L += stages_fwd_ops(tB, sB, uB)
L += [
("cx", accU, accT),
("ccx", accT, tB, accS),
("cx", accU, accT),
("cx", tB, accT),
("cx", sB, accS),
("ccx", accU, uB, accS),
("cx", uB, accU),
]
return L
⋯ 2 lines
def mul_inv_ops(g):
"""Chain push of the link INVERSE (not the pop of mul_fwd_ops -- pops
are exact macro mirrors of whatever was pushed).
"""
p0, p1, q0, q1, r0, r1, sB, tB, uB = g
L = stages_inv_ops(tB, sB, uB)
L += [
("cx", accU, accT),
("cx", uB, accT),
("ccx", accT, tB, accS),
("cx", uB, accT),
("cx", accU, accT),
("ccx", accU, uB, accS),
("cx", tB, accT),
("cx", sB, accS),
("cx", uB, accS),
("cx", tB, accS),
("cx", uB, accU),
]
fams = [("d", p1, p0)]
if CO_Q:
fams.append(("k", Qm, (r0, r1), 3 - CO_Q))
if CO_R:
fams.append(("k", Rm, (q0, q1), 3 - CO_R))
L += prot_cluster(fams)
L += core3_ops(Qm, q1, q0, xorB0)
L += core3_ops(Rm, r1, r0, xorB0)
return L
6

Restarts and inversions

invload_ops, invert_macros

An anchor opens a segment at the identity. invload_ops writes the state of g−1 directly from the stored code: P̃ = −p, Q = −q, R = −r, υ = u, T = t, and S = s ⊕ u ⊕ t ⊕ tu. A trit negation is a CNOT onto swapped rails. Only the final AND costs a Toffoli, so each of the 64 anchors costs one Toffoli instead of a 35-Toffoli multiplication.

invert_macros replaces W by W−1 in place for 21 Toffolis. It negates the three trits for free. A υ-controlled phase correction reads its cells from J, and six υ-controlled cell swaps apply the Xo permutation. Stages under t and s follow, and one Toffoli fixes M. The table cells also absorb the correction for the stored shear P̃.

def invload_ops(g):
"""Chain anchor: load the primed state of g^-1 from the identity rails
in 1 Toffoli (a 35-Toffoli multiply would do the same from the
identity; the 64 scheduled anchors save 34 Toffoli each, 2,176
total). From the stored code
(p2,q,r,s,t,u) the inverse state is exactly P' = -p2, Q = -q, R = -r,
U = u, T = t, S = s^u^t^(t&u); trit negations are rail-swapped CX.
"""
p0, p1, q0, q1, r0, r1, sB, tB, uB = g
return [
("cx", p1, accP[1]),
("cx", p0, accP[2]),
("cx", p0, accP[0]),
("cx", p1, accP[0]),
("cx", q1, Qm[1]),
("cx", q0, Qm[2]),
("cx", q0, Qm[0]),
("cx", q1, Qm[0]),
("cx", r1, Rm[1]),
("cx", r0, Rm[2]),
("cx", r0, Rm[0]),
("cx", r1, Rm[0]),
("cx", uB, accU),
("cx", tB, accT),
("cx", sB, accS),
("cx", uB, accS),
("cx", tB, accS),
("ccx", tB, uB, accS),
]
⋯ 317 lines
def invert_macros():
def _sw(a, bq):
return [("cx", a, bq), ("cx", bq, a), ("cx", a, bq)]
L = []
L += _sw(accP[1], accP[2]) + _sw(Qm[1], Qm[2]) + _sw(Rm[1], Rm[2])
L += J_build_ops()
L += J_fold_ops()
# U-conditioned P-correction table stage on the J cells (one-hot: the
# cell group's OR is its XOR, collected with free CX)
for c in _IG1:
L.append(("cx", J[3 * c[0] + c[1]], scr0))
for c in _IG2:
L.append(("cx", J[3 * c[0] + c[1]], xorB0))
L += [("andc", ctlB[0], accU, scr0), ("andc", ctlB[1], accU, xorB0)]
L += core3_ops(accP, ctlB[0], ctlB[1], addd)
L += [("andr", ctlB[1], accU, xorB0), ("andr", ctlB[0], accU, scr0)]
for c in _IG2:
L.append(("cx", J[3 * c[0] + c[1]], xorB0))
for c in _IG1:
L.append(("cx", J[3 * c[0] + c[1]], scr0))
# B swaps on J conditioned on U
for i, jj in SWAPS_X:
L.append(("cswap", accU, J[i], J[jj]))
# re-derive Qm/Rm from the permuted grid; tear J down
L += J_fold_ops()
L += J_unbuild_ops()
# W1 = neg_R then derived-rail Q<->R swap, both conditioned on t
L += [("cswap", accT, Rm[1], Rm[2])]
L += qr_swap_ops(accT)
# W2 = -I conditioned on s
L += [("cswap", accS, Qm[1], Qm[2]), ("cswap", accS, Rm[1], Rm[2])]
# M-fix: accS ^= (accT & ~U) ^ U
L += [("x", accU), ("ccx", accU, accT, accS), ("x", accU), ("cx", accU, accS)]
return L
7

Read the trace, add, erase

build_ext_ops, kadd_ops

The trace class f(W) takes seven values, and six of them are nonzero. build_ext_ops sets six flag wires. Each nonzero value sets one or two of them, and zero sets none. It decodes M into four row wires, splits them by υ into (M, υ) rails, and builds J for (Q, R).

Each class indicator is a table over rails and cells. An exhaustive scan at build time writes it over F2 as rail and cell parities plus as few paid products as possible, one Toffoli per product. The indicators combine into the six flags. A Clifford-only teardown then frees the linear scaffolding, so the adder can borrow those zeroed wires as carries. The readout costs 20 Toffolis.

kadd_consts assigns each flag a fixed constant, so that the XOR of the constants of the set flags equals the weighted addend. It records, for each bit of k, which flags contribute. kadd_ops folds that parity onto one addend wire per bit and ripples with one Toffoli per interior bit. The return sweep releases each carry by measurement. The plaquette ripple costs 50 Toffolis and the rectangle ripple 49, since aR is even. With six Toffolis to erase the flags, a plaquette costs 76 and a rectangle 75, against 478 for a lookup.

def build_ext_ops():
"""Map the chain state to the 6-wire trace-class basis
(cen2, f6, t1, f2, t2, fm1): decode M into row rails, build J, then
synthesize the e-trit / nu / sg indicators from the F2 scans and
combine. Ends with a Clifford-only teardown of the linear scaffolding
so those wires can serve as the full-width adder's carry pool (rebuilt
for free by the mirror on the way back).
"""
kind = EBEST[1] # ty: ignore[not-subscriptable] # EBEST set by the module scan
pa, fa, oa = bx_args(EBEST[2]) # ty: ignore[not-subscriptable] # EBEST set by the module scan
pb, fb, ob = bx_args(EBEST[3]) # ty: ignore[not-subscriptable] # EBEST set by the module scan
pn, fn, on = bx_args(DNU)
ps, fs, osg = bx_args(DSG)
USED = set()
for pp, ff in ((pa, fa), (pb, fb), (pn, fn), (ps, fs)):
for Ax, _Bx, _cf in pp:
USED.update(Ax)
USED.update(ff)
USED &= set(ROWQ)
E: list = [
("andc", a_st, accS, accT),
("al", r10),
("cx", accS, r10),
("cx", a_st, r10),
("al", r01),
("cx", accT, r01),
("cx", a_st, r01),
("al", r00),
("x", r00),
("cx", accS, r00),
("cx", accT, r00),
("cx", a_st, r00),
]
built0 = []
for m, rq in ((0, r00), (1, r01), (2, r10), (3, a_st)):
n0 = (rMU0[m] in USED) or (m == 0)
n1 = (rMU1[m] in USED) or n0
if n1:
E.append(("andc", rMU1[m], rq, accU))
if n0:
E += [("al", rMU0[m]), ("cx", rq, rMU0[m]), ("cx", rMU1[m], rMU0[m])]
built0.append((m, rq))
E += J_build_ops()
w1, w2 = candB[0], candB[1]
if kind == "12":
pre = [("bx", w1, pa, fa, oa), ("bx", w2, pb, fb, ob)]
post = [("mx", w2, pb, fb, ob), ("mx", w1, pa, fa, oa)]
elif kind == "1n":
pre = [("bx", w1, pa, fa, oa), ("bx", w2, pb, fb, ob), ("cx", w1, w2)]
post = [("cx", w1, w2), ("mx", w2, pb, fb, ob), ("mx", w1, pa, fa, oa)]
else:
pre = [("bx", w2, pa, fa, oa), ("bx", w1, pb, fb, ob), ("cx", w2, w1)]
post = [("cx", w2, w1), ("mx", w1, pb, fb, ob), ("mx", w2, pa, fa, oa)]
i1 = (EIDX - 1) % 3
i2 = (EIDX - 2) % 3
ecp = (((w1,), (accP[i1], accP[EIDX]), 0), ((w2,), (accP[i2], accP[EIDX]), 0))
ecf = (accP[EIDX],)
E += pre
E += [("bx", ecq, ecp, ecf, 0)]
E += [("bx", nuq, pn, fn, on), ("bx", sgq, ps, fs, osg)]
E += [
("andc", t1q, nuq, ecq),
("al", t2q),
("cx", nuq, t2q),
("cx", t1q, t2q),
("x", sgq),
("andc", f2q, t1q, sgq),
("andc", fm1q, t2q, sgq),
("x", sgq),
("andc", cen2q, J[0], rMU0[0]),
("andc", f6q, cen2q, ecq),
]
E += [("mx", ecq, ecp, ecf, 0)]
E += post
# Clifford-only teardown: reclaim the linear scaffolding for the
# full-width adder carry pool (rebuilt for free on the way back).
# Deliberately partial -- the four AND cells stay live -- and
# deliberately ordered cols (J[1], J[2]) before rows (J[3], J[6]),
# the reverse of the J_unbuild_ops order. The two
# blocks commute (each writes only cells the other does not read),
# but the free order feeds the Builder LIFO free-list, so
# harmonizing it would permute later qubit ids and move the stream.
for m, rq in built0:
E += [("cx", rq, rMU0[m]), ("cx", rMU1[m], rMU0[m]), ("fr", rMU0[m])]
E += [
("cx", accS, r10),
("cx", a_st, r10),
("fr", r10),
("cx", accT, r01),
("cx", a_st, r01),
("fr", r01),
("x", r00),
("cx", accS, r00),
("cx", accT, r00),
("cx", a_st, r00),
("fr", r00),
]
E += [("cx", J[i], J[0]) for i in range(1, 9)] + [("x", J[0]), ("fr", J[0])]
for bb in (1, 2):
E += [
("cx", Rm[bb], J[bb]),
("cx", J[3 + bb], J[bb]),
("cx", J[6 + bb], J[bb]),
("fr", J[bb]),
]
for a in (1, 2):
E += [
("cx", Qm[a], J[3 * a]),
("cx", J[3 * a + 1], J[3 * a]),
("cx", J[3 * a + 2], J[3 * a]),
("fr", J[3 * a]),
]
return E
⋯ 9 lines
def kadd_consts(mult):
"""Per-bit XOR-parity matrices for adding mult*f through the basis
wires: pats[j] is the constant added when basis wire j is set.
"""
Wm = 1 << WBITS
cm = {v: (v * mult) % Wm for v in OUTV}
pats = [cm[-3], cm[6] ^ cm[-3], cm[-2], cm[2] ^ cm[-2], cm[1], cm[-1] ^ cm[1]]
mats = [[BASISQS[jj] for jj in range(6) if (pats[jj] >> i) & 1] for i in range(51)]
lo = min((cv & -cv).bit_length() - 1 for cv in cm.values() if cv)
return mats, lo
⋯ 19 lines
def kadd_ops(mats, xconst=0, lo=0, hi=51, *, pool):
"""In-place k += addend over bits [lo, hi): per bit the addend parity
is folded onto addd (free CX), one AND makes the carry on a pooled
wire, and the return sweep releases every carry by measurement
(1 Toffoli per interior bit; Gidney ripple).
local: gadgets.add/gidney_add is the same Gidney (2018) 1-T/interior-bit
ripple but requires a register operand wire per addend position
(shared= restores one repeated wire; it cannot refold a per-bit XOR
parity of the 6 basis wires onto the single addd wire), which is
the inexpressible part; its pool= parameter already covers pooled
carries.
Carries come only from the fixed pool, and pool sufficiency is
asserted fail-closed: no dynamic allocation exists, so a wider add
cannot silently push past the 1,859-wire peak.
"""
L = []
def put(i):
for fq in mats[i]:
L.append(("cx", fq, addd))
if (xconst >> i) & 1:
L.append(("x", addd))
carries = {}
c = None
for i in range(lo, hi):
put(i)
if i < hi - 1:
idx = i - lo
assert idx < len(pool), "carry pool exhausted"
nc = pool[idx]
if nc in BREATHING:
L.append(("al", nc))
if c is None:
L.append(("ccx", kreg[i], addd, nc))
L.append(("cx", addd, kreg[i]))
else:
L += [
("cx", c, kreg[i]),
("cx", c, addd),
("ccx", kreg[i], addd, nc),
("cx", c, nc),
("cx", c, addd),
("cx", addd, kreg[i]),
]
carries[i] = nc
c = nc
else:
if c is not None:
L.append(("cx", c, kreg[i]))
L.append(("cx", addd, kreg[i]))
put(i)
for i in range(hi - 2, lo - 1, -1):
cc = carries[i]
cprev = carries.get(i - 1)
put(i)
L.append(("cx", addd, kreg[i]))
if cprev is not None:
L += [("cx", cprev, addd), ("cx", cprev, cc)]
L.append(("andr", cc, kreg[i], addd))
if cprev is not None:
L.append(("cx", cprev, addd))
L.append(("cx", addd, kreg[i]))
put(i)
return L
8

Recode the links

encode

The circuit allocates nine utility wires that live for the whole run. One loop over the 192 links then recodes each link in place. Its body is one recorded unit, XF_FWD, whose link wires are affine in the loop index.

The unit is xform_ops followed by qr_corr_ops, 23 Toffolis per link. xform_ops rewrites the code from X to Xo, which changes only p, q, and r. qr_corr_ops adds the link's own qr term to its stored phase, so no multiplication computes that term again. The recode is paid once per link, and all 3,136 multiplications use the recoded form.

Recode links 4,416 Toffolis · peak 1,788 qubits
with b.phase("encode"):
play([("alloc", _w) for _w in UTIL])
with b.repeat(192) as _i:
play(XF_FWD, wire=lambda w: LNK[9 * _i + w] if w < 9 else None)
9

The walk across rows and planes

loops

The walk allocates the 12 wires of W and sets the zero rail of each trit, which is the identity. The opening row plays its 36 chunks with absolute link ids.

Every later row reuses one recorded body. Between rows, y_shift rotates the link array by one row and z_shift by one plane, with 1,296 Swaps each and no Toffolis or allocations. After a shift, the next row's links sit on the same frame wires, so emit_row replays the same chunks.

Plane 0 runs three translated rows. Each of the planes 1 to 3 runs one plane-crossing row and three translated rows. The 20 shifts return every link to its own wires. The tail then closes W at the identity.

# ---- 2. loops: the conveyor walk + 576 trace-adds ---------------------------
Open the walk 10,030 Toffolis · peak 1,859 qubits
with b.phase("loops/id-set"):
play([("alloc", _w) for _w in CHAINW] + [("x", accP[0]), ("x", Qm[0]), ("x", Rm[0])])
with b.phase("loops/row-open"):
for _ops in OPEN_CHUNKS:
play(_ops)
with b.repeat(3): # plane 0, rows y = 1..3
Shift the links 0 Toffolis · peak 1,800 qubits
with b.phase("loops/shift"):
y_shift()
Row program 122,472 Toffolis · peak 1,859 qubits
with b.phase("loops/row-y"):
emit_row(Y_CHUNKS)
with b.phase("loops/shift"):
y_shift()
with b.repeat(3): # planes z = 1..3
with b.phase("loops/shift"):
z_shift()
Plane program 30,618 Toffolis · peak 1,859 qubits
with b.phase("loops/row-z"):
emit_row(Z_CHUNKS)
with b.repeat(3):
with b.phase("loops/shift"):
y_shift()
with b.phase("loops/row-y"):
emit_row(Y_CHUNKS)
with b.phase("loops/shift"):
y_shift()
with b.phase("loops/shift"):
z_shift()
Close the walk 176 Toffolis · peak 1,809 qubits
with b.phase("loops/close"):
play(TAIL_OPS)
10

Clear W and restore the links

loops/id-clear, decode

With W back at the identity, the circuit clears the zero rails and frees the 12 wires of W.

Decode mirrors encode exactly. The reversed unit runs over the links in reverse order, again 23 Toffolis per link, and the circuit frees the utility wires. Only the links and k remain live, and k holds k + aPSP + aRSR mod 251.

with b.phase("loops/id-clear"):
play([("x", accP[0]), ("x", Qm[0]), ("x", Rm[0])] + [("free", _w) for _w in CHAINW])
# ---- 3/4. decode: exact mirror of encode; release the utilities ------------
Restore links 4,416 Toffolis · peak 1,788 qubits
with b.phase("decode"):
with b.repeat(192) as _i:
play(XF_REV, wire=lambda w: LNK[1719 - 9 * _i + w] if w < 9 else None)
play([("free", _w) for _w in UTIL])
Start (init): 0 Toffolis, peak 1,779 qubitsRecode links (encode): 4,416 Toffolis, peak 1,788 qubitsOpen the walk (open): 10,030 Toffolis, peak 1,859 qubitsShift the links (shift): 0 Toffolis, peak 1,800 qubitsShift the links (shift): 0 Toffolis, peak 1,800 qubitsShift the links (shift): 0 Toffolis, peak 1,800 qubitsShift the links (shift): 0 Toffolis, peak 1,800 qubitsRow program (row-y): 122,472 Toffolis, peak 1,859 qubitsRow program (row-y): 122,472 Toffolis, peak 1,859 qubitsPlane program (row-z): 30,618 Toffolis, peak 1,859 qubitsClose the walk (close): 176 Toffolis, peak 1,809 qubitsRestore links (decode): 4,416 Toffolis, peak 1,788 qubits1,8591,779peak live qubits per phase, above the 1,779 I/O qubits1.7 × 10⁵ Toffolis →
# C. The internal X representative and the mirror algebra. The chain works
# with an internal generator Xo (in place of the contract Xc) chosen so
# that X^2 = V^2, X V X^-1 = V^-1, and
# Xc*Xo^-1 lies in Delta(27) = <w,C,E> -- then the
# U-flag composes linearly with the M-counter and the encode transform
# (contract code -> internal code) is cheap. No tie-break is needed: the
# algebraic filter admits EXACTLY ONE candidate (asserted; a group-model
# change that widens the set fails the build).
# ---------------------------------------------------------------------------
XO_CANDS = []
for g in sorted(G):
if g in sig36:
continue
if mpow(g, 2) != V2m:
continue
if mmul(mmul(g, Vg), mdag(g)) != Vinvm:
continue
if mmul(Xc, mdag(g)) in delta:
XO_CANDS.append(g)
assert len(XO_CANDS) == 1, "internal X representative is unique"
Xo = XO_CANDS[0]
⋯ 208 lines
# F. The loop schedule. Conventions pinned by PARAMS-DERIVATION.md
# sections 3 and 7 / reference.py: site s = x + 4y + 16z, link l = 3s + mu,
# inverses on return links.
# ---------------------------------------------------------------------------
# The global chain-walk schedule. Ops: 'f'/'i' push the link / its
# inverse (35-Toffoli multiplies); 'li'/'ul' open/close at the identity
# with a 1-Toffoli inverse-load / its reverse; 'v' inverts the register in
# place (21 Toffoli) so subsequent right-multiplies edit the word's other
# end; 'XP'/'XR' fire a plaquette/rectangle trace-add. Exact DP optimum
# of the net-flip-relaxed morph model (odd inversion counts per morph
# allowed): 3,136 multiplies + 480 inversions + 64 anchors.
#
# The DP optimum is row-periodic: every y-junction row of the x-fastest
# raster is one lattice translate of a single row program, every
# plane-crossing row another, so constants.json stores THREE row
# programs (open_row, then y_row x12 and z_row x3 by translation) plus
# the closing tail, and this build repeats them through a slab conveyor
# (section L) instead of unrolling 4,256 tokens. y_row/z_row link ids
# are FRAME-LOCAL: link 3*(x + 4*dy + 16*dz) + mu addresses the link dy
# rows / dz planes above the row being processed (mod 4); open_row and
# tail are absolute. gen_constants.py runs the search and writes the
# file.
_CONSTS = json.loads((pathlib.Path(__file__).resolve().parent / "constants.json").read_bytes())
def _decode_prog(blobs):
KIND = {"f": "f", "i": "i", "l": "li", "u": "ul", "v": "inv", "P": "XP", "R": "XR"}
out = []
for tok in blobs.split():
kd = KIND[tok[0]]
out.append((kd, int(tok[1:])) if tok[0] in "filu" else (kd,))
return out
OPEN_ROW = _decode_prog(" ".join(_CONSTS["open_row"]))
Y_ROW = _decode_prog(" ".join(_CONSTS["y_row"]))
Z_ROW = _decode_prog(" ".join(_CONSTS["z_row"]))
TAIL = _decode_prog(" ".join(_CONSTS["tail"]))
⋯ 13 lines
# H. Wire map. Only the 1,779 I/O wires (192 links x 9, then k) are
# persistent. The 80-wire working set lives and dies with its data: the
# nine hot utility wires span the whole program and the 12-wire chain
# state spans the loops stage (explicit stage-boundary lifecycles); the
# other 59 wires -- the J grid, the extraction scaffolding and class
# basis, and the 30 dedicated adder carries -- are allocated when their gadget
# builds them and freed when it (or its measurement mirror) tears them
# down. At every plaquette trace-add ripple all 80 are simultaneously
# live, so the skyline peaks at exactly 1,859 (rectangle adds top out at
# 1,858: LO_R = 1 drops one carry) while breathing everywhere else.
# ---------------------------------------------------------------------------
NIO = 1779
links = [[9 * l + i for i in range(9)] for l in range(192)] # wires 0..1727
kreg = list(range(1728, 1779)) # wires 1728..1778
accP = [1779, 1780, 1781] # one-hot rails of the primed phase trit P'
Qm = [1782, 1783, 1784] # one-hot rails of the C-exponent trit Q
Rm = [1785, 1786, 1787] # one-hot rails of the E-exponent trit R
accS, accT, accU = 1788, 1789, 1790 # V-counter bits M = 2S+T; X flag U
BANK = list(range(1791, 1829)) # the 38-wire pooled ancilla bank
candB = BANK[0:4]
xorB0 = BANK[4]
ctlB = BANK[5:7]
J = BANK[7:16] # 3x3 one-hot product grid J[3q+r] = [Q=q][R=r]
a_st = BANK[16]
r10 = BANK[17]
r01 = BANK[18]
r00 = BANK[19] # M-decoder
rMU1 = BANK[20:24]
rMU0 = BANK[24:28] # (M,U) one-hot row rails
nuq = BANK[28]
sgq = BANK[29] # unit-class / negative-sign flags
scr0 = BANK[30]
t1q = BANK[31]
f2q = BANK[32]
t2q = BANK[33]
fm1q = BANK[34]
cen2q = BANK[35]
f6q = BANK[36] # trace-class basis wires
addd = BANK[37] # the full-width adder's addend wire
ecq = ctlB[0] # e-copy flag lives on ctlB[0]: free during extraction
CARRY = list(range(1829, 1859)) # 30 dedicated full-width adder carries
⋯ 949 lines
XF_UNIT = xform_ops(links[0]) + qr_corr_ops(links[0])
XF_FWD = expand(XF_UNIT)
XF_REV = expand(XF_UNIT, rev=True)
def morph_ops(toks):
"""Expand a run of walk tokens (no fires) to flat low-level ops."""
L = []
for op in toks:
k = op[0]
if k == "f":
L += expand(mul_fwd_ops(links[op[1]]))
elif k == "i":
L += expand(mul_inv_ops(links[op[1]]))
elif k == "li":
L += expand(invload_ops(links[op[1]]))
elif k == "ul":
L += expand(invload_ops(links[op[1]]), rev=True)
elif k == "inv":
L += expand(INVERT_M)
else:
raise AssertionError(op)
return L
def row_chunks(form):
"""Lower one row program to its 36 (morph + trace-add) chunks."""
chunks, cur = [], []
for op in form:
if op[0] not in ("XP", "XR"):
cur.append(op)
continue
mats, lo = (MATS_P, LO_P) if op[0] == "XP" else (MATS_R, LO_R)
ops_ = morph_ops(cur)
ops_ += expand(EXT_OPS)
ops_ += expand(kadd_ops(mats, lo=lo, hi=WBITS, pool=KPOOL))
ops_ += expand(EXT_OPS, rev=True)
chunks.append(ops_)
cur = []
assert not cur
return chunks
OPEN_CHUNKS = row_chunks(OPEN_ROW)
Y_CHUNKS = row_chunks(Y_ROW)
Z_CHUNKS = row_chunks(Z_ROW)
TAIL_OPS = morph_ops(TAIL)
assert max(len(c) for c in OPEN_CHUNKS + Y_CHUNKS + Z_CHUNKS) < 0xFFFF
⋯ 79 lines
def y_shift():
"""Rotate the link array one y-slab down (next row -> frame row):
1,296 Swaps.
"""
with b.repeat(4) as _zb, b.repeat(108) as _w:
_e = 432 * _zb + _w
b.swap(LNK[_e], LNK[_e + 108])
b.swap(LNK[_e + 108], LNK[_e + 216])
b.swap(LNK[_e + 216], LNK[_e + 324])
def z_shift():
"""Rotate the link array one z-slab down (next plane -> frame
plane): 1,296 Swaps.
"""
with b.repeat(432) as _w:
b.swap(LNK[_w], LNK[_w + 432])
b.swap(LNK[_w + 432], LNK[_w + 864])
b.swap(LNK[_w + 864], LNK[_w + 1296])
def emit_row(chunks):
for _ops in chunks:
with b.repeat(1):
play(_ops)
⋯ 2 lines
# ---- 0/1. utility wires; encode: one affine repeat over the 192 links ------
Recode links 4,416 Toffolis · peak 1,788 qubits
with b.phase("encode"):
play([("alloc", _w) for _w in UTIL])
with b.repeat(192) as _i:
play(XF_FWD, wire=lambda w: LNK[9 * _i + w] if w < 9 else None)
# ---- 2. loops: the conveyor walk + 576 trace-adds ---------------------------
Open the walk 10,030 Toffolis · peak 1,859 qubits
with b.phase("loops/id-set"):
play([("alloc", _w) for _w in CHAINW] + [("x", accP[0]), ("x", Qm[0]), ("x", Rm[0])])
with b.phase("loops/row-open"):
for _ops in OPEN_CHUNKS:
play(_ops)
with b.repeat(3): # plane 0, rows y = 1..3
Shift the links 0 Toffolis · peak 1,800 qubits
with b.phase("loops/shift"):
y_shift()
Row program 122,472 Toffolis · peak 1,859 qubits
with b.phase("loops/row-y"):
emit_row(Y_CHUNKS)
with b.phase("loops/shift"):
y_shift()
with b.repeat(3): # planes z = 1..3
with b.phase("loops/shift"):
z_shift()
Plane program 30,618 Toffolis · peak 1,859 qubits
with b.phase("loops/row-z"):
emit_row(Z_CHUNKS)
with b.repeat(3):
with b.phase("loops/shift"):
y_shift()
with b.phase("loops/row-y"):
emit_row(Y_CHUNKS)
with b.phase("loops/shift"):
y_shift()
with b.phase("loops/shift"):
z_shift()
Close the walk 176 Toffolis · peak 1,809 qubits
with b.phase("loops/close"):
play(TAIL_OPS)
with b.phase("loops/id-clear"):
play([("x", accP[0]), ("x", Qm[0]), ("x", Rm[0])] + [("free", _w) for _w in CHAINW])
# ---- 3/4. decode: exact mirror of encode; release the utilities ------------
Restore links 4,416 Toffolis · peak 1,788 qubits
with b.phase("decode"):
with b.repeat(192) as _i:
play(XF_REV, wire=lambda w: LNK[1719 - 9 * _i + w] if w < 9 else None)
play([("free", _w) for _w in UTIL])