TwoQubitWeylDecomposition.circuit() has the wrong global phase for about half of all unitaries #17401

Environment

• Qiskit version: 2.6.1
• Python version: 3.12
• Operating system: macOS

What is happening?

TwoQubitWeylDecomposition(U).circuit() reproduces U only up to a global phase, and the phase is wrong by exactly pi (the circuit equals -U) for roughly half of Haar-random 4x4 unitaries. The decomposition's own attributes (K1l, K1r, K2l, K2r, a, b, c, global_phase) reassemble U exactly, so only circuit() is affected.

How can we reproduce the issue?

import numpy as np
from qiskit.quantum_info import Operator, random_unitary
from qiskit.synthesis import TwoQubitWeylDecomposition

bad = 0
for seed in range(2000):
    u = random_unitary(4, seed=seed).data
    v = Operator(TwoQubitWeylDecomposition(u).circuit()).data
    bad += np.abs(u - v).max() > 1e-9
print(bad, "of 2000")          # 1017 of 2000

# one failing case
u = random_unitary(4, seed=5).data
v = Operator(TwoQubitWeylDecomposition(u).circuit()).data
print(np.angle(np.trace(v.conj().T @ u)))   # -3.14159...  (u = -v)
print(np.abs(u + v).max())                  # ~1e-15

What should happen?

Operator(TwoQubitWeylDecomposition(U).circuit()) equals U including its global phase, as the attributes do.

Any suggestions?

In crates/synthesis/src/two_qubit_decompose/weyl_decomposition.rs, TwoQubitWeylDecomposition::circuit. The phase of each single-qubit Euler sequence is added to global_phase, but the code adds c2r.global_phase a second time where c1r.global_phase is meant (line ~1143), and c1l.global_phase is never added:

        for gate in c1r.gates { ... }
        global_phase += c2r.global_phase;     // should be c1r.global_phase
        ...
        for gate in c1l.gates { ... }
        gate_sequence.set_global_phase_f64(global_phase);   // c1l.global_phase is missing