Transpilation benchmarks

Note

Qiskit has no synthesis method that does what GULPS does, so neither comparison is like for like. The default synthesis uses one kind of entangling gate per block, which is not optimal when the Target offers several. XXDecomposer finds least-cost sequences for these gates, but it runs in Python and has not been ported to Rust.

These benchmarks transpile QFT, EfficientSU2, quantum-volume, and random circuits of 4 to 64 qubits for an all-to-all Target with U and \(R_{ZZ}(\pi/2)\), \(R_{ZZ}(\pi/4)\), and \(R_{ZZ}(\pi/6)\). Duration is the sum of the gate durations, without scheduling.

Against Qiskit’s default synthesis

Qiskit’s default synthesis produces circuits 1.5 to 2 times longer than GULPS on every circuit except EfficientSU2. It synthesizes each block with TwoQubitBasisDecomposer over \(R_{ZZ}(\pi/2)\) alone, which is locally equivalent to CX, while GULPS also uses the two shorter gates.

import json
from pathlib import Path

rows = json.loads(Path("docs/_scripts/fullcircuit.json").read_text())
qubits = sorted({r["qubits"] for r in rows})
families = list(dict.fromkeys(r["circuit"] for r in rows))
ratio = {
    (r["circuit"], r["qubits"]): r["qiskit_duration_us"] / r["gulps_duration_us"]
    for r in rows
}
print("Default / GULPS summed duration")
print(f"{'Circuit':<14}" + "".join(f"{f'{n} qubits':>11}" for n in qubits))
for family in families:
    print(f"{family:<14}" + "".join(f"{ratio[family, n]:>11.2f}" for n in qubits))

Default / GULPS summed duration
Circuit          4 qubits   8 qubits  16 qubits  32 qubits  64 qubits
QFT                  1.70       1.83       1.87       1.89       1.90
EfficientSU2         1.00       1.00       1.00       1.00       1.00
QV                   1.60       1.59       1.60       1.58       1.58
Random               2.01       1.79       1.65       1.65       1.64

EfficientSU2 entangles with CX, so both outputs use \(R_{ZZ}(\pi/2)\) only.

Against XXDecomposer

GULPS transpiles tens to hundreds of times faster than XXDecomposer, with the same summed two-qubit duration in 19 of the 20 cases. To make the two minimize the same cost, each \(R_{ZZ}\) gate gets the fidelity (1 - ERR_BASE) ** p, where p is its duration relative to \(R_{ZZ}(\pi/2)\), so XXDecomposer’s largest fidelity product is the least summed duration.

XXDecomposer’s transpile time divided by GULPS’s. Each point is the geometric mean of two rounds, and the bars span both rounds.

The one case with different duration

In the random circuit at 64 qubits, XXDecomposer’s output is shorter by one sixth of the \(R_{ZZ}(\pi/2)\) duration. The difference comes from one block with Weyl coordinates \((11/24,11/24,0)\). XXDecomposer synthesizes this block with one \(R_{ZZ}(\pi/6)\) and three \(R_{ZZ}(\pi/4)\), whose interaction strengths sum to \(11/12\), the sum of the block’s coordinates. The block therefore lies on the boundary of the region of this cheaper sentence. XXDecomposer’s reachability tolerance accepts the block, and its output has process infidelity up to \(9.5\times10^{-13}\) against it. GULPS uses two \(R_{ZZ}(\pi/2)\) gates.

Measurement setup

All pipelines ran with Qiskit 2.5.2 and GULPS built from this repository, pinned to one core of an AMD Ryzen 5 5600X under WSL2. BLAS, OpenMP, and Rayon each used one thread. Each measurement has two rounds with fresh decomposers and pass managers, and the pipeline order changes between rounds.

GATES = [
    ("zz", np.pi / 2, 1.0),
    ("sq2_zz", np.pi / 4, 1 / 2),
    ("sq3_zz", np.pi / 6, 1 / 3),
]
DUR_BASE, DUR_1Q, ERR_BASE = 500e-9, 52.3e-9, 1e-3
ROUNDS = 2


def build_target(n: int) -> Target:
    """An all-to-all Target on ``n`` qubits with the three RZZ gates and U."""
    target = Target()
    pairs = list(permutations(range(n), 2))
    for name, angle, scale in GATES:
        props = InstructionProperties(
            duration=DUR_BASE * scale, error=1 - (1 - ERR_BASE) ** scale
        )
        target.add_instruction(RZZGate(angle), {p: props for p in pairs}, name=name)
    u_err = 1 - (1 - ERR_BASE) ** (DUR_1Q / DUR_BASE)
    u_props = {
        (q,): InstructionProperties(duration=DUR_1Q, error=u_err) for q in range(n)
    }
    target.add_instruction(
        UGate(Parameter("t"), Parameter("p"), Parameter("l")), u_props
    )
    return target

class XXPass(TransformationPass):
    """Qiskit's mixed-strength XXDecomposer on every consolidated two-qubit block."""

    def __init__(self) -> None:
        """Build one decomposer with the three RZZ strengths and their fidelities."""
        super().__init__()
        embodiments, fidelities = {}, {}
        for _, angle, scale in GATES:
            qc = QuantumCircuit(2)
            qc.h([0, 1])
            qc.rzz(angle, 0, 1)
            qc.h([0, 1])
            # XXDecomposer maximizes the product of fidelities. With fidelity
            # (1 - ERR_BASE)**scale it minimizes the summed scale, as GULPS does.
            embodiments[angle], fidelities[angle] = qc, (1 - ERR_BASE) ** scale
        self.xx = XXDecomposer(fidelities, euler_basis="U", embodiments=embodiments)

    def run(self, dag: DAGCircuit) -> DAGCircuit:
        """Replace each two-qubit unitary node with its decomposition."""
        for node in dag.op_nodes():
            if node.op.name == "unitary" and node.op.num_qubits == 2:
                dag.substitute_node_with_dag(
                    node, self.xx(node.op.to_matrix(), approximate=False, use_dag=True)
                )
        return dag

def pass_managers(target: Target) -> dict[str, PassManager]:
    """The GULPS, XXDecomposer, and default UnitarySynthesis pipelines."""
    unroll, merge1q = (
        Unroll3qOrMore(target=target),
        Optimize1qGatesDecomposition(target=target),
    )
    consolidate = ConsolidateBlocks(force_consolidate=True)
    # The decomposer is built from the same gates and durations as the Target.
    gulps = GulpsDecomposer(
        [RZZGate(angle) for _, angle, _ in GATES], [scale for *_, scale in GATES]
    )
    return {
        "gulps": PassManager([unroll, GulpsDecompositionPass(gulps), merge1q]),
        "xx": PassManager([unroll, consolidate, XXPass(), merge1q]),
        "qiskit": PassManager(
            [unroll, consolidate, UnitarySynthesis(target=target), merge1q]
        ),
    }

Each pipeline runs once on a small circuit before it is timed.

def benchmark() -> list[dict]:
    """Transpile every circuit at every size with the three pipelines."""
    rows = []
    for n in QUBITS:
        target = build_target(n)
        warm = QuantumCircuit(n)
        warm.unitary(np.eye(4), [0, 1])
        warm.unitary(Operator(RZZGate(0.3)), [0, 1])
        for i, (family, qc) in enumerate(circuits(n).items()):
            row = dict(circuit=family, qubits=n, gulps_s=[], xx_s=[], qiskit_s=[])
            for rnd in range(ROUNDS):
                # Fresh pipelines per round: GULPS keeps each block class's sentence
                # across calls, so a repeated circuit would be timed warm.
                managers = pass_managers(target)
                for pm in managers.values():
                    pm.run(warm)
                keys = list(managers)
                shift = (rnd + i) % len(keys)
                for key in keys[shift:] + keys[:shift]:
                    start = time.perf_counter()
                    out = managers[key].run(qc)
                    row[key + "_s"].append(time.perf_counter() - start)
                    if rnd == 0:
                        if n <= 8:
                            assert Operator(out).equiv(Operator(qc)), (key, family, n)
                        row[key + "_duration_us"] = duration_us(out)
                        row[key + "_2q_duration_us"] = duration_us(out, one_qubit=False)
            row["same_cost"] = math.isclose(
                row["xx_2q_duration_us"],
                row["gulps_2q_duration_us"],
                rel_tol=SAME_COST_RTOL,
            )
            row.update(audit_blocks(qc, target))
            rows.append(row)
            print(
                f"{family:12s} {n:2d}q  GULPS {min(row['gulps_s']):.4f}-{max(row['gulps_s']):.4f}s"
                f"  XX {min(row['xx_s']):.4f}-{max(row['xx_s']):.4f}s"
                f"  default {min(row['qiskit_s']):.4f}-{max(row['qiskit_s']):.4f}s"
                f"  duration GULPS {row['gulps_duration_us']:.4f} XX {row['xx_duration_us']:.4f}"
                f" default {row['qiskit_duration_us']:.4f} us"
                f"  two-qubit GULPS {row['gulps_2q_duration_us']:.4f}"
                f" XX {row['xx_2q_duration_us']:.4f} us"
                f"  XX inexact blocks {row['blocks_xx_inexact']}",
                flush=True,
            )
    return rows

Plotting code

def plot(
    rows: list[dict], out: Path = HERE.parent / "_static" / "fullcircuit.svg"
) -> None:
    """Draw the XXDecomposer/GULPS ratio of transpile time."""
    with plt.rc_context({"svg.hashsalt": "fullcircuit"}):
        fig, ax = plt.subplots(figsize=(6.5, 3.4), constrained_layout=True)
        extent = []
        for family, marker in zip(FAMILIES, MARKERS):
            fam = sorted(
                (r for r in rows if r["circuit"] == family and r["xx_s"]),
                key=lambda r: r["qubits"],
            )
            qubits = [r["qubits"] for r in fam]
            rounds = np.array([np.divide(r["xx_s"], r["gulps_s"]) for r in fam])
            ratio = np.exp(np.log(rounds).mean(axis=1))
            extent += [rounds.min(), rounds.max()]
            ax.errorbar(
                qubits,
                ratio,
                yerr=(ratio - rounds.min(axis=1), rounds.max(axis=1) - ratio),
                marker=marker,
                capsize=1.5,
                elinewidth=0.8,
                label={"QV": "Quantum volume"}.get(family, family),
            )
        # The y range ends at the powers of two just outside the data.
        lo = int(np.floor(np.log2(min(extent))))
        hi = int(np.ceil(np.log2(max(extent))))
        if lo <= 0:
            ax.axhline(1, color="0.5", lw=0.8, ls="--", zorder=0)
        ax.set_xscale("log", base=2)
        ax.set_yscale("log", base=2)
        ax.set_ylim(2.0**lo, 2.0**hi)
        ax.set_xticks(QUBITS, [str(q) for q in QUBITS])
        ticks = [2**k for k in range(lo, hi + 1)]
        ax.set_yticks(ticks, [f"{t:g}" for t in ticks])
        ax.minorticks_off()
        ax.set_xlabel("Qubits")
        ax.set_ylabel("Transpile time, XXDecomposer / GULPS")
        fig.legend(loc="outside upper center", ncol=4)
        fig.savefig(out, metadata={"Date": None})

Reproduce the measurements

docs/_scripts/fullcircuit.py writes docs/_scripts/fullcircuit.json and the figure; --plot redraws the figure from the JSON file. Pinned to one core:

taskset -c 2 python docs/_scripts/fullcircuit.py

The run takes about three minutes, most of it in XXDecomposer at 32 and 64 qubits.