Fixed-point Gaussian

The Gaussian e−x² in fixed point. The baseline squares its input first; the discovered circuit never forms the square, and skips most of its table, because beyond |x| ≈ 4 the Gaussian is smaller than the allowed error.

In brief

The baseline (Häner et al., 2018) computes x2 with a full multiplication, then evaluates cubics on 32 intervals of x2.

The discovered circuit evaluates a quadratic in a local coordinate on each of 1,024 equal cells of |x|. Folding negative inputs costs no Toffolis, because a shortcut that is off by one input step stays inside the error budget. Cells where the Gaussian is below the error budget return zero, and the table lookup skips every subtree whose words are all equal.

TQ = 5·10⁵TQ = 10⁶TQ = 2·10⁶TQ = 5·10⁶2003004002,0005,00010⁴2·10⁴peak qubits Q →Toffoli count T →Discovered circuitHäner et al. (2018)Live qubits as Toffolis execute24k1.6k362139
Removing the square removes one large product and its register.
1

Registers, wiring, and the fold

init, fold

The circuit keeps two 39-qubit registers: the input x, a two's-complement number with 32 fractional bits, and the accumulator k, which receives round(232e−x²) modulo 239. Fixed slices of x do the rest of the indexing. Ten bits addr_msb select one of 1,024 equal segments of width 2−7, the low 25 bits Dq are the coordinate within the segment, and Uq holds the top eleven bits of that coordinate.

The Gaussian is even, so the circuit folds negative inputs onto positive ones. Thirty-eight CNOTs from the sign bit complement the other bits. An exact absolute value would also add one; the circuit omits that carry chain. The folded value then sits one input unit below |x|, which moves the output by at most 0.858 output units, inside the error budget.

After the fold, bits 35 to 37 of x are zero for every input in [−8, 8). The circuit lends them to a wire pool, and the first three slope bits of the table word bind to them. This part of the circuit has no Toffolis.

b = Builder()
_r = declare_registers(b, [("x", N_X), ("k", KW)])
XQ, KQ = _r["x"], _r["k"]
b.claim("x", XQ)
b.claim("k", KQ)
# ------------------------------------------------------------ wiring
Dq = XQ[0:25]
Uq = XQ[25 - M : 25]
addr_msb = [XQ[34 - i] for i in range(10)]
low_msb = [XQ[28], XQ[27], XQ[26], XQ[25]]
hi_msb = [XQ[34], XQ[33], XQ[32], XQ[31], XQ[30], XQ[29]]
# ---- fold: sign-fold x, pool the certified-zero spares, allocate words
Fold 0 Toffolis · peak 78 qubits
with b.phase("fold/sign"):
for i in range(38):
b.cx(XQ[38], XQ[i])
wp = WirePool(b, SPARES)
b.release("x", SPARES) # the spares leave x while lent to the pool
with b.phase("fold/words"):
# The pool bank empties into the slope immediately: qFw[0..2] bind
# the spares here (pooled allocs emit no op; the wires are the
# persistent certified-zero x bits). Every other packed-word wire
# materializes at the first leaf that sets its bit, and a bit no
# address ever sets fails the post-walk coverage assert instead of
# riding as a dead wire.
TW = FW + EWID + A0BITS
T = [None] * TW
for _j in range(3):
T[_j] = wp.qubit()
assert set(T[:3]) == set(SPARES) # the pool bank empties into the slope
b.claim("slope", T[:3])
2

A lookup that stops at uniform subtrees

lookup

The segment address selects a 55-bit word: a 25-bit constant a0, an 11-bit signed slope correction c, and 19 stored slope bits s. The two lowest slope bits are always zero, so the circuit wires them as implied zeros. A binary tree over the ten address bits forms one AND per descent and keeps at most nine live. At a leaf, leaf_main writes the word with conditioned CNOTs.

The walk prunes by content. A subtree whose words are all zero is skipped (nz_main), and a subtree whose words are all equal gets one conditioned write at its root (unif_main). Of the 1,024 segments, 455 hold quadratics, 59 hold constants, and 510 return zero, because the Gaussian there is below the error budget. Pruning lowers the walk from 1,022 Toffolis to 495.

A word wire comes into existence at the first leaf that sets its bit, so no wire is live before it holds data. The lookup reaches the peak of the whole circuit: 78 input and accumulator wires, 52 fresh word wires, and nine branch ANDs.

01234567801e^(−x²)|x| →skippedskippedskipped455 quadratic59 constant510 zero
Each divider bounds a 128-cell subtree; the lookup skips the last three subtrees, which are all zero.
def _matT(j):
if T[j] is None:
T[j] = wp.qubit()
b.claim("slope" if j < FW else "correction" if j < FW + EWID else "constant", [T[j]])
return T[j]
# ---- lookup: content-pruned QROM walk writes the packed word --------
CUM = [0] * (NADDR + 1)
for a in range(NADDR):
CUM[a + 1] = CUM[a] + (1 if WORDS[a] else 0)
def nz_main(lo, span):
return CUM[lo + span] - CUM[lo] > 0
def unif_main(lo, span):
w = WORDS[lo]
if w == 0:
return False
for a in range(lo + 1, lo + span):
if WORDS[a] != w:
return False
return True
def leaf_main(a, wire):
w = WORDS[a]
j = 0
while w:
if w & 1:
if wire is None:
b.x(_matT(j))
else:
b.cx(wire, _matT(j))
w >>= 1
j += 1
def lmark(lo, span):
if span == 128:
b.set_phase(f"lookup/b{lo // 128}")
# descend = the walk levels above the first 128-address subtree; each
# b<n> mark then owns that subtree's leaves plus the pops back to the
# next one
Lookup 495 Toffolis · peak 139 qubits
with b.phase("lookup/descend"):
pruned_walk(
b, addr_msb, 0, NADDR, leaf_main, nz_main, uniform=unif_main, mark=lmark, pool=wp
)
3

Add the constant, then measure it away

a0add

The constant enters the accumulator directly: k ← k + a0·27. The argument gl sets how many low carry positions are erased by measurement (one Toffoli and one extra live qubit each); the remaining positions use in-place carries (two Toffolis, no extra qubit).

The circuit no longer needs the constant word, so it measures each bit in the X basis and returns the wires to the pool. Each outcome leaves a known sign that depends on the address. The circuit keeps the outcomes in sbA0 and repairs those signs at the end, together with the signs of the other two fields.

Add constant 53 Toffolis · peak 139 qubits
with b.phase("a0add/add"):
add(b, KQ[A0SH:], qA0w, gl=GL0, pool=wp)
with b.phase("a0add/meas"):
sbA0 = []
for q in qA0w:
sbA0.append(b.hmr(q))
b.release("constant", qA0w) # demolished: ownership ends before the frees
for q in reversed(qA0w):
wp.recycle(q)
4

Correct the slope in place

mid_mul

The table stores the quadratic of each segment as a constant and a slope that changes linearly across the segment. Six radix-4 Booth rows add ⌊ut c 22t−15⌋ to the slope register, where ut is digit t of Uq. Because Uq approximates the coordinate divided by 214, the sum supplies the correction cD/229.

The update changes the slope register in place; no product register exists.

for t in range(NT):
Slope update 198 Toffolis · peak 137 qubits
with b.phase(f"mid_mul/t{t}"):
mid_row(t, False)
5

One Booth row

acc_row, mid_row, booth_row

Radix-4 Booth recoding reads overlapping bit pairs of the multiplier as signed digits in {−2, …, 2}. In mid_row, CNOTs form m1, which marks magnitude one, and z. One AND in acc_row turns them into m2, which marks magnitude two. At most one of the two controls is active, so the one-times and two-times partial products of a column share one ancilla, built by booth_columns. A negative digit complements the target around the same addition instead of building a negated operand.

After the in-place addition, booth_columns_release erases the column ANDs by measurement, newest first. In the slope rows, one shared AND, the sign wire, extends the signed correction through all high columns. With invert set, the row complements the slope around the addition and so subtracts; that is how the circuit later restores the slope.

booth_row uses the same skeleton for the main product. It keeps only the columns at or above drop, and it picks gl, the number of measurement-erased carries, as large as the qubit ceiling of the lookup allows.

def acc_row(m1, partner, flip, cols, wrap, T, colneg=None):
for u, v in flip:
b.cx(u, v)
b.x(m1)
m2 = wp.and_compute(m1, partner)
b.x(m1)
for u, v in flip:
b.cx(u, v)
ts, reps = booth_columns(b, m1, m2, cols, neg=colneg, pool=wp)
A, akw, unwrap = wrap(m2, ts)
add(b, T, A, pool=wp, **akw)
unwrap()
booth_columns_release(b, ts, reps, neg=colneg, pool=wp)
for u, v in flip:
b.cx(u, v)
b.x(m1)
wp.and_release(m2, m1, partner)
b.x(m1)
for u, v in flip:
b.cx(u, v)
⋯ 6 lines
def mid_row(t, invert):
j0, j1, j2 = 2 * t - 1, 2 * t, 2 * t + 1
SH1 = M + 4 - 2 * t
qlo = max(FLOW, EQ - SH1 + 1)
alloc_m1 = j0 >= 0
if alloc_m1:
m1 = wp.qubit()
b.cx(Uq[j1], m1)
b.cx(Uq[j0], m1)
else:
m1 = Uq[j1]
z = wp.qubit()
b.claim("recode", ([m1] if alloc_m1 else []) + [z])
if j2 < M:
b.cx(Uq[j2], z)
b.cx(Uq[j1], z)
neg = Uq[j2] if j2 < M else None
q2top = EQ + EWID - SH1
cols = [(qEw[q + SH1 - EQ], qEw[q + SH1 - EQ - 1]) for q in range(qlo, q2top)]
def wrap(m2, ts):
# sign wire (covers columns q2top .. FW-1)
b.cx(m1, m2)
sw = wp.and_compute(m2, qEw[EWID - 1])
b.cx(m1, m2)
ext = FW - (q2top - FLOW)
ext = max(ext, 0)
for qq in qFw:
if invert:
b.x(qq)
if neg is not None:
b.cx(neg, qq)
def unwrap():
for qq in qFw:
if neg is not None:
b.cx(neg, qq)
if invert:
b.x(qq)
cbit = b.hmr(sw)
b.cx(m1, m2)
b.cz(m2, qEw[EWID - 1], cond=cbit)
b.cx(m1, m2)
wp.recycle(sw)
return ts + [sw] * ext, {"shared": sw}, unwrap
acc_row(m1, z, [], cols, wrap, qFw[qlo - FLOW :])
if j2 < M:
b.cx(Uq[j2], z)
b.cx(Uq[j1], z)
b.release("recode", [z])
wp.recycle(z)
if alloc_m1:
b.cx(Uq[j1], m1)
b.cx(Uq[j0], m1)
b.release("recode", [m1])
wp.recycle(m1)
⋯ 8 lines
def booth_row(i, drop):
sh = FBITS - 2 * i
f1w = qFw[2 * i - FLOW]
have2 = 2 * i + 1 < FBITS
f2w = qFw[2 * i + 1 - FLOW] if have2 else None
Aw = qFw[2 * i - 1 - FLOW]
b.cx(f1w, Aw)
neg = f2w if have2 else None
cols = [
(Dq[sh + q] if sh + q <= 24 else None, Dq[sh + q - 1]) for q in range(drop, 2 * i + 5)
]
m = len(cols) + 1
# [GID]-prefix length: the largest that keeps this row's runway
# moment (row wires + carries + remaining headroom term) under
# CAPRAW
gl = max(0, min(m - 1, CAPRAW - ROWBASE - m - (33 - 2 * i)))
def wrap(g, ts):
b.cx(Aw, g)
if neg is not None:
b.cx(neg, g)
def unwrap():
if neg is not None:
b.cx(neg, g)
b.cx(Aw, g)
return [g] + ts, {"neg": neg, "gl": gl}, unwrap
acc_row(
Aw,
f2w if have2 else f1w,
[(f1w, f2w)] if have2 else [],
cols,
wrap,
KQ[drop - 1 :],
colneg=neg,
)
b.cx(f1w, Aw)
6

Multiply into the accumulator

final_rows

The corrected slope times the coordinate gives both the linear and the quadratic term of the segment, and the circuit subtracts that product from k directly: k ← k − ⌊s′D/221⌋. No register holds D2 or the full product. Complementing k before and after the rows turns the subtraction into additions.

At full precision the product would need more rows; nine rows can reach the retained bits of k. Each row keeps at most 18 columns and leaves the five lowest bits of k untouched (drop). An offline search over the circuit's cost chose these cutoffs, and an exhaustive scan over every coordinate checks the error bound.

Product into k 589 Toffolis · peak 139 qubits
with b.phase("final_rows/neg"):
for q in KQ:
b.x(q)
for i in range(I0, 11):
with b.phase(f"final_rows/i{i}"):
drop = max(RB, (2 * i + 5) - CAPW)
booth_row(i, drop)
with b.phase("final_rows/unneg"):
for q in KQ:
b.x(q)
7

Restore the slope

mid_unmul

The six slope rows run again in reverse order with invert set, which subtracts each correction. The slope register then holds the looked-up slope exactly. The measurement-based erasure that follows needs this: each measured bit must equal a stored table bit, so that the sign it leaves is a known function of the address.

for t in reversed(range(NT)):
Slope restored 198 Toffolis · peak 137 qubits
with b.phase(f"mid_unmul/t{t}"):
mid_row(t, True)
8

Measure the words and repair all signs in one pass

unlook, unfold

Reversing the lookup would cost another 495 Toffolis. The circuit instead measures the slope and correction fields in the X basis, appends those outcomes to the ones from the constant, and repairs every recorded sign in one walk.

The three borrowed input wires return to x once the pool is empty again, and the same 38 CNOTs undo the fold, leaving x unchanged.

fixup(qFw + qEw, sbA0, WORDS, low_msb, hi_msb, 16)
assert len(wp._free) == 3 and set(wp._free) == set(SPARES), "spares in flight"
b.claim("x", SPARES) # the certified-zero spares re-join x
# ---- unfold ---------------------------------------------------------
Unfold 0 Toffolis · peak 78 qubits
with b.phase("unfold"):
for i in range(38):
b.cx(XQ[38], XQ[i])
9

The repair walk

fixup

A measured table bit with outcome m multiplies the amplitude of address j by (−1)m b(j), where b(j) is the stored bit. The outcome is classical, so the repair is a known Boolean function of the address. All three fields are measured before the repair starts.

The repair groups addresses into blocks of sixteen. build_onehot expands the four low address bits into a rail with exactly one active wire, in layers of 2, 4, 8, and 16 wires. A pruned walk over the six high bits then reaches each block that holds a nonzero word, and its leaf applies, for each set bit of a stored word, one CZ conditioned on that bit's measurement outcome. The rail costs 14 ANDs and the walk 36, so the whole repair takes 50 Toffolis. un_onehot then clears the rail.

def build_onehot(bits_msb):
hb = bits_msb[0]
w1 = wp.qubit()
b.cx(hb, w1)
w0 = wp.qubit()
b.x(w0)
b.cx(hb, w0)
b.claim("onehot", [w0, w1])
cur = [w0, w1]
levs = []
for bit in bits_msb[1:]:
lev = []
new = []
for w in cur:
hi = wp.and_compute(w, bit)
lo = wp.qubit()
b.cx(w, lo)
b.cx(hi, lo)
new.extend([lo, hi])
lev.append((w, bit, lo, hi))
b.claim("onehot", new)
levs.append(lev)
cur = new
return cur, ((hb, w0, w1), levs)
⋯ 1 line
def un_onehot(meta):
(hb, w0, w1), levs = meta
for lev in reversed(levs):
b.release("onehot", [q for (_w, _bit, lo, hi) in lev for q in (lo, hi)])
for w, bit, lo, hi in reversed(lev):
b.cx(hi, lo)
b.cx(w, lo)
wp.recycle(lo)
wp.and_release(hi, w, bit)
b.release("onehot", [w0, w1])
b.cx(hb, w0)
b.x(w0)
wp.recycle(w0)
b.cx(hb, w1)
wp.recycle(w1)
⋯ 15 lines
def fixup(qubits, pre_sb, wbits, low_msb, hi_msb, lam):
# The parameters are coupled to the closure: the demolition order
# must match the slope + correction releases, and lam must be the
# rail size the leaf indexes.
assert qubits == qFw + qEw and lam == 1 << len(low_msb)
Lookup erased 50 Toffolis · peak 110 qubits
with b.phase("unlook/meas"):
sb = []
for q in qubits:
sb.append(b.hmr(q))
b.release("slope", qFw) # demolished: ownership ends
b.release("correction", qEw) # before the frees
for q in reversed(qubits):
wp.recycle(q)
sb = sb + list(pre_sb)
with b.phase("unlook/onehot"):
oh, meta = build_onehot(low_msb)
assert len(oh) == lam
def nzf(lo, span):
a0 = lo * lam
a1 = min((lo + span) * lam, NADDR)
for a in range(a0, a1):
if wbits[a]:
return True
return False
def leaf(h, wire):
for l in range(lam):
a = h * lam + l
if a >= NADDR:
break
w = wbits[a]
j = 0
while w:
if w & 1:
b.cz(wire, oh[l], cond=sb[j])
w >>= 1
j += 1
def wmark(lo, span):
if span == 8:
b.set_phase(f"unlook/w{lo}")
pruned_walk(b, hi_msb, 0, 1 << len(hi_msb), leaf, nzf, mark=wmark, pool=wp)
with b.phase("unlook/unhot"):
un_onehot(meta)
Start (init): 0 Toffolis, peak 78 qubitsFold (fold): 0 Toffolis, peak 78 qubitsLookup (lookup): 495 Toffolis, peak 139 qubitsAdd constant (a0add): 53 Toffolis, peak 139 qubitsSlope update (mid_mul): 198 Toffolis, peak 137 qubitsProduct into k (final_rows): 589 Toffolis, peak 139 qubitsSlope restored (mid_unmul): 198 Toffolis, peak 137 qubitsLookup erased (unlook): 50 Toffolis, peak 110 qubitsUnfold (unfold): 0 Toffolis, peak 78 qubits13978peak live qubits per phase, above the 78 I/O qubits1,583 Toffolis →
b = Builder()
_r = declare_registers(b, [("x", N_X), ("k", KW)])
XQ, KQ = _r["x"], _r["k"]
b.claim("x", XQ)
b.claim("k", KQ)
# ------------------------------------------------------------ wiring
Dq = XQ[0:25]
Uq = XQ[25 - M : 25]
addr_msb = [XQ[34 - i] for i in range(10)]
low_msb = [XQ[28], XQ[27], XQ[26], XQ[25]]
hi_msb = [XQ[34], XQ[33], XQ[32], XQ[31], XQ[30], XQ[29]]
# ---- fold: sign-fold x, pool the certified-zero spares, allocate words
Fold 0 Toffolis · peak 78 qubits
with b.phase("fold/sign"):
for i in range(38):
b.cx(XQ[38], XQ[i])
wp = WirePool(b, SPARES)
b.release("x", SPARES) # the spares leave x while lent to the pool
with b.phase("fold/words"):
# The pool bank empties into the slope immediately: qFw[0..2] bind
# the spares here (pooled allocs emit no op; the wires are the
# persistent certified-zero x bits). Every other packed-word wire
# materializes at the first leaf that sets its bit, and a bit no
# address ever sets fails the post-walk coverage assert instead of
# riding as a dead wire.
TW = FW + EWID + A0BITS
T = [None] * TW
for _j in range(3):
T[_j] = wp.qubit()
assert set(T[:3]) == set(SPARES) # the pool bank empties into the slope
b.claim("slope", T[:3])
def _matT(j):
if T[j] is None:
T[j] = wp.qubit()
b.claim("slope" if j < FW else "correction" if j < FW + EWID else "constant", [T[j]])
return T[j]
# ---- lookup: content-pruned QROM walk writes the packed word --------
CUM = [0] * (NADDR + 1)
for a in range(NADDR):
CUM[a + 1] = CUM[a] + (1 if WORDS[a] else 0)
def nz_main(lo, span):
return CUM[lo + span] - CUM[lo] > 0
def unif_main(lo, span):
w = WORDS[lo]
if w == 0:
return False
for a in range(lo + 1, lo + span):
if WORDS[a] != w:
return False
return True
def leaf_main(a, wire):
w = WORDS[a]
j = 0
while w:
if w & 1:
if wire is None:
b.x(_matT(j))
else:
b.cx(wire, _matT(j))
w >>= 1
j += 1
def lmark(lo, span):
if span == 128:
b.set_phase(f"lookup/b{lo // 128}")
# descend = the walk levels above the first 128-address subtree; each
# b<n> mark then owns that subtree's leaves plus the pops back to the
# next one
Lookup 495 Toffolis · peak 139 qubits
with b.phase("lookup/descend"):
pruned_walk(
b, addr_msb, 0, NADDR, leaf_main, nz_main, uniform=unif_main, mark=lmark, pool=wp
)
⋯ 8 lines
Add constant 53 Toffolis · peak 139 qubits
with b.phase("a0add/add"):
add(b, KQ[A0SH:], qA0w, gl=GL0, pool=wp)
with b.phase("a0add/meas"):
sbA0 = []
for q in qA0w:
sbA0.append(b.hmr(q))
b.release("constant", qA0w) # demolished: ownership ends before the frees
for q in reversed(qA0w):
wp.recycle(q)
⋯ 96 lines
for t in range(NT):
Slope update 198 Toffolis · peak 137 qubits
with b.phase(f"mid_mul/t{t}"):
mid_row(t, False)
⋯ 44 lines
Product into k 589 Toffolis · peak 139 qubits
with b.phase("final_rows/neg"):
for q in KQ:
b.x(q)
for i in range(I0, 11):
with b.phase(f"final_rows/i{i}"):
drop = max(RB, (2 * i + 5) - CAPW)
booth_row(i, drop)
with b.phase("final_rows/unneg"):
for q in KQ:
b.x(q)
# ---- mid_unmul: exact reverse of mid_mul returns F to B1 ------------
for t in reversed(range(NT)):
Slope restored 198 Toffolis · peak 137 qubits
with b.phase(f"mid_unmul/t{t}"):
mid_row(t, True)
⋯ 110 lines
fixup(qFw + qEw, sbA0, WORDS, low_msb, hi_msb, 16)
assert len(wp._free) == 3 and set(wp._free) == set(SPARES), "spares in flight"
b.claim("x", SPARES) # the certified-zero spares re-join x
# ---- unfold ---------------------------------------------------------
Unfold 0 Toffolis · peak 78 qubits
with b.phase("unfold"):
for i in range(38):
b.cx(XQ[38], XQ[i])