1. Compilation/
  2. Diagonal unitary decomposition

Diagonal unitary decomposition

Here, we show how the decomposition of a diagonal unitary operator comes about. We will first look at the single-qubit case and then generalize. We will also outline the recursion that leads to the complete decomposition.

Single-qubit diagonal unitary

A single-qubit unitary operator is described by a 2×22\times 2 matrix UU with the property UU=I2UU^\dagger=\mathbb{I}_2. As a result, a diagonal unitary operator Δ\Delta on one qubit has a diagonal with two entries, aa and bb, that satisfy

ΔΔ=(a00b)(a00b)=(1001)    a=eiφ0,b=eiφ1.\Delta\Delta^\dagger = \begin{pmatrix} a & 0 \\ 0 & b \end{pmatrix} \begin{pmatrix} a^\ast & 0 \\ 0 & b^\ast \end{pmatrix} = \begin{pmatrix} 1 & 0 \\ 0 & 1 \end{pmatrix} \ \ \Rightarrow\ \ a = e^{-i\varphi_0}, b=e^{-i\varphi_1}.

We will decompose this diagonal operator into a product of global phase and an RZR_Z rotation. As both are also diagonal unitary operators, they take the same form as Δ\Delta, but with φ0=φ1\varphi_0=\varphi_1 and φ0=φ1\varphi_0=-\varphi_1 for qml.GlobalPhase and qml.RZ, respectively. Denoting the global phase angle as ϕ\phi and the RZR_Z rotation angle as θ\theta, we find for their product:

(ei2θ00ei2θ)(eiϕ00eiϕ)=(ei(ϕ+12θ)00ei(ϕ12θ))\begin{pmatrix} e^{-\frac{i}{2}\theta} & 0 \\ 0 & e^{\frac{i}{2}\theta} \end{pmatrix} \begin{pmatrix} e^{-i\phi} & 0 \\ 0 & e^{-i\phi} \end{pmatrix} = \begin{pmatrix} e^{-i(\phi + \tfrac{1}{2}\theta)} & 0 \\ 0 & e^{-i(\phi - \tfrac{1}{2}\theta)} \end{pmatrix}

If we want to find those angles ϕ\phi and θ\theta that reproduce the diagonal unitary operator Δ\Delta, we just need to equate the matrices and solve for the angles:

θ=φ0φ1,\theta = \varphi_0 - \varphi_1,
ϕ=12(φ0+φ1).\phi = \frac{1}{2}(\varphi_0 + \varphi_1).

To conclude this section, we draw the single-qubit decomposition:

Decomposition of a single-qubit diagonal unitary operator

Generalizing to multiple qubits

To generalize the decomposition from one to multiple qubits, we take two steps: first, we establish the general form this decomposition should assume, and second, we calculate the corresponding angles.

The first step results from the so-called Multiplexer Extension Property (MEP) [1], which states that attaching multiplexer nodes to a circuit identity yields a valid circuit identity again. That is, given a circuit identity

0: ─A─┤ = ─B─┤,

for parametrized circuits AA and BB, we may conclude that

0: ─╭◻─┤ = ─╭○───╭●──┤ = ─╭○───╭●──┤ = ─╭◻─┤
1: ─╰A─┤   ─╰A0──╰A1─┤   ─╰B0──╰B1─┤   ─╰B─┤,

where A0A_0 and A1A_1 are two instances of the parametrized circuit AA, and B0B_0 and B1B_1 are the equivalent instances of BB. Applying multiple multiplexing nodes follows the same logic and results in valid circuit identities again.

Note that a multi-qubit diagonal operator can be interpreted as a multiplexed single-qubit diagonal operator, which generally takes the form of a block diagonal matrix. This is because a diagonal matrix is a special case of a block-diagonal matrix. Therefore, we can directly apply the MEP to the single-qubit decomposition from above and add n1n-1 multiplexing nodes to obtain the following decomposition of a multi-qubit diagonal unitary:

Decomposition of a multi-qubit diagonal unitary operator

Where we used the fact that a multiplexed global phase is the same as a diagonal operator on the controlling qubits.

For the second step, i.e., the calculation of the parameters in the decomposition, we may again generalize the single-qubit scenario. For each basis state on the control qubits, the multiplexed phase, which then becomes a smaller diagonal operator, applies the same phase to two basis states of the target qubit. Similarly, the multiplexed RZR_Z rotation applies phases with opposite sign to the target, with individual phases chosen for each control state. The matrix equation for the decomposition above then reads

diag(eiφ0eiφ1eiφ2eiφ3eiφ2n1)=diag(ei2θ0ei2θ0ei2θ1ei2θ1ei2θ2n11)diag(eiϕ0eiϕ0eiϕ1eiϕ1eiϕ2n11).\operatorname{diag}\begin{pmatrix} e^{-i\varphi_0} \\ e^{-i\varphi_1} \\ e^{-i\varphi_2} \\ e^{-i\varphi_3} \\ \vdots \\ e^{-i\varphi_{2^{n}-1}} \end{pmatrix} = \operatorname{diag}\begin{pmatrix} e^{-\tfrac{i}{2}\theta_0} \\ e^{\tfrac{i}{2}\theta_0} \\ e^{-\tfrac{i}{2}\theta_1} \\ e^{\tfrac{i}{2}\theta_1} \\ \vdots \\ e^{\tfrac{i}{2}\theta_{2^{n-1}-1}} \end{pmatrix} \operatorname{diag}\begin{pmatrix} e^{-i\phi_0} \\ e^{-i\phi_0} \\ e^{-i\phi_1} \\ e^{-i\phi_1} \\ \vdots \\ e^{-i\phi_{2^{n-1}-1}} \end{pmatrix}.

Now, we can apply the same strategy as in the single-qubit case to each pair of matrix entries. That is, we can compute the parameters θp\theta_p and ϕp\phi_p from the ppth pair of phases {exp(iφ2p),exp(iφ2p+1)}\{\exp(-i\varphi_{2p}), \exp(-i\varphi_{2p+1})\} from the original diagonal that we aim to implement:

ei(ϕp+12θp)=eiφ2p+1e^{-i(\phi_p+\tfrac{1}{2}\theta_p)} = e^{-i\varphi_{2p\phantom{+1}}}
ei(ϕp12θp)=eiφ2p+1.e^{-i(\phi_p-\tfrac{1}{2}\theta_p)} = e^{-i\varphi_{2p+1}}.

For this, we extract the phase angles φ\vec{\varphi} (e.g., using np.angle), and compute

ϕp=12(φ2p+φ2p+1),θp=φ2pφ2p+1,\phi_p = \tfrac{1}{2}(\varphi_{2p} + \varphi_{2p+1}),\quad \theta_p = \varphi_{2p} - \varphi_{2p+1},

to find the required parameters. This means that we just need to chain together np.angle, some indexing, np.mean and subtraction in order to compute the parameters. This does not require much computational effort, which is fortunate because a diagonal operator on nn qubits has 2n2^n entries and quickly becomes expensive to handle. An implementation of this can be found in PennyLane's qml.DiagonalQubitUnitary.

Recursive decomposition

To conclude, we can apply the decomposition from above recursively to break down a multi-qubit diagonal unitary operator into multiplexed RZR_Z rotations (and a global phase):

Recursive decomposition of a multi-qubit diagonal into multiplexed rotations

The multiplexed RZR_Z rotations, qml.SelectPauliRot in PennyLane, can then be decomposed individually. For more, see the Select-U(2) compilation page.