GULPS (Global Unitary Linear Programming Synthesis) is the first open tool that robustly compiles arbitrary two-qubit unitaries optimally into non-standard instruction sets.

Most existing compilers only target CNOT gates. Analytical rules exist for a few special cases like fractional CNOT (XX family), Berkeley (B), and $\sqrt{\text{iSWAP}}$, but nothing more general. Numerical methods can in principle handle arbitrary gates, but they are slow, unreliable, and do not scale as instruction sets grow. GULPS fills this gap by combining linear programming with lightweight numerics to achieve:

• Support for fractional, continuous, or heterogeneous gate sets.
• Scalability to larger ISAs, unlike black-box numerical methods.
• A fast, practical tool integrated with Qiskit if you study gate compilation from two-body Hamiltonians or parameterized unitary families.

pip install gulps

Optional extras:

If your backend's ISA is already defined in a Qiskit Target, GULPS works as a drop-in translation stage plugin:

from qiskit import transpile

output_qc = transpile(input_qc, target=my_target, translation_method="gulps")

For full control, define your ISA manually. Gate costs must be additive (e.g., normalized durations where fractional gates cost proportionally to their basis gate).

from qiskit.circuit.library import iSwapGate
from gulps import GulpsDecomposer
from gulps.core.isa import DiscreteISA

isa = DiscreteISA(
    gate_set=[iSwapGate().power(1 / 2), iSwapGate().power(1 / 3)],
    costs=[1 / 2, 1 / 3],
    names=["sqrt2iswap", "sqrt3iswap"],
)
decomposer = GulpsDecomposer(isa=isa)

Once initialized, call the decomposer with a Qiskit Gate or a 4x4 np.ndarray:

from qiskit.quantum_info import random_unitary

u = random_unitary(4, seed=0)
v = decomposer(u)
v.draw()

To compile a full QuantumCircuit, use the TransformationPass. Because GULPS leaves single-qubit gates unsimplified, append Optimize1qGatesDecomposition to rewrite them:

from gulps import GulpsDecompositionPass
from qiskit.transpiler import PassManager
from qiskit.transpiler.passes import Optimize1qGatesDecomposition
from qiskit.circuit.random import random_circuit

input_qc = random_circuit(4, 4, max_operands=2)
pm = PassManager(
    [
        GulpsDecompositionPass(decomposer),
        Optimize1qGatesDecomposition(basis="u3"),
    ]
)
output_qc = pm.run(input_qc)
output_qc.draw("mpl")

The decomposition begins by identifying the cheapest feasible basis gate sentence (a sequence of native gates sufficient to construct the target unitary). We use monodromy polytopes to describe the reachable space of canonical invariants for each sentence in the ISA.

For example, this ISA:

from gulps.core.isa import DiscreteISA

isa = DiscreteISA(
    gate_set=[iSwapGate().power(1 / 2), iSwapGate().power(1 / 3)],
    costs=[1 / 2, 1 / 3],
    names=["sqrt2iswap", "sqrt3iswap"],
    precompute_polytopes=True,
)

has the following coverage set:

from gulps.core.coverage import coverage_report

coverage_report(isa.coverage_set)

Once a sentence is chosen, a linear program is used to determine a trajectory of intermediate invariants. These represent the cumulative two-qubit nonlocal action after each gate in the sentence, starting from the identity and ending at the target.

from gulps.core.invariants import GateInvariants
from gulps.viz.invariant_viz import plot_decomposition

example_input = random_unitary(4, seed=31)
target_inv = GateInvariants.from_unitary(example_input)
constraint_sol = decomposer._best_decomposition(target_inv=target_inv)
plot_decomposition(
    constraint_sol.intermediates, constraint_sol.sentence, decomposer.isa
);

In this example, the optimal sentence is composed of 2 $\sqrt[3]{\texttt{iSWAP}}$ gates and 1 $\sqrt[2]{\texttt{iSWAP}}$. That is, the resulting circuit falls into a parameterized ansatz like this:

Unlike other decomposition techniques, the linear program contains additional information about the intermediate points used to reduce the problem into simpler subproblems, each corresponding to a depth-2 circuit segment. In this case, the circuit has three segments, although the first red segment (beginning at Identity is trivial). That leaves two segments requiring synthesis:

Red(2)	Blue

We solve for the local one-qubit gates in each segment using a Gauss-Newton solver on the Makhlin invariants, followed by a Weyl-coordinate polish. The solver (implemented in Rust) is tuned to work well across a broad range of ISAs, but there is no one-size-fits-all for every possible gate set, so edge-case performance may vary.

solutions = decomposer._local_synthesis._solve_segments(
    constraint_sol.sentence,
    constraint_sol.intermediates,
    n_inner=len(constraint_sol.sentence) - 1,
)
solutions[0]
# SegmentSolution(u0=..., u1=..., weyl_residual=1.2e-16, max_residual=4.8e-09, success=True)

After solving the individual segments, we apply a final stitching step to handle orientation between segments and to promote local equivalence into global unitary equivalence:

decomposer._local_synthesis.synthesize_segments(
    gate_list=constraint_sol.sentence,
    invariant_list=constraint_sol.intermediates,
    target=target_inv,
).draw("mpl")

See more:

• https://quantum-journal.org/papers/q-2020-03-26-247/
• https://quantum-journal.org/papers/q-2022-04-27-696/
• https://threeplusone.com/pubs/on_gates.pdf
• https://chromotopy.org/latex/papers/xx-synthesis.pdf
• qiskit-advocate/qamp-spring-23#33
• Qiskit/qiskit#9375
• https://weylchamber.readthedocs.io/en/latest/readme.html