Add one_qubit_decomposition function

[trbromley]

Feature details

A single-qubit unitary matrix $U$ can be decomposed into X, Y, and Z rotations. In PennyLane, it is possible to perform the zyz_decomposition or the xyx_decomposition. We would like to add a single one_qubit_decomposition function, similar to two_qubit_decomposition that allows a user to decompose any single-qubit unitary matrix $U$ into chosen rotations. In addition, we would like to add a zxz_decomposition into only Z and X rotations.

Implementation

The new one_qubit_decomposition function should have the signature:

one_qubit_decomposition(U: tensor/array, rotations: string = "ZYZ", wire) -> list[Operation]

In addition, when rotations="ZXZ", a new zxz_decomposition function will be called:

zxz_decomposition(U: tensor/array, wire, return_global_phase=False) -> list[Operation]

These additions should live within the transforms module, and should include documentation, examples, and tests.

How important would you say this feature is?

1: Not important. Would be nice to have.

Additional information

No response

[jiggychauhi]

I would like to be assigned on this

[trbromley]

Hi @jiggychauhi! That's awesome, you can check out our development guide or just post here if you have any questions. Good luck!

[abhishekabhishek]

hi @trbromley just wondering how many rotation strings should we have support for i.e. should we consider all the standard 6: {XYX, XZX, YXY, YZY, ZXZ, ZYZ} ?

[trbromley]

Hi @abhishekabhishek. The minimum we're looking for is ZYZ and XYX (which currently exist) and the addition of ZXZ. Though, I'm open to more combinations being added. Thanks!

[tnemoz]

I'd have some questions on this issue.

• From what i've seen, the zyz_decomposition and xyx_decomposition work quite differently. In particular zyz_decomposition returns a single Rot gate, while xyx_decomposition returns a succession of RX and RY gates. My guess is that one_qubit_decomposition should return a list of gates, even if rotations is set to "ZYZ", is that correct?

• From my understanding, the xyx_decomposition function returns the gates in the wrong order. The documentation states

  Returns a list of of gates, an RX, an RY and another RX gate, which when applied in the order of appearance in the list is equivalent to the unitary U up to global phase.

  However, take the following code:

>>> import pennylane as qml
>>> from pennylane import numpy as np
>>> def rx(x):
...     return np.array([[np.cos(x/2), -np.sin(x/2)*1j], [-np.sin(x/2)*1j, np.cos(x/2)]])
... 
>>> def ry(x):
...     return np.array([[np.cos(x/2), -np.sin(x/2)], [np.sin(x/2), np.cos(x/2)]])
... 
>>> U = np.array([
...     [-0.28829348-0.78829734j,  0.30364367+0.45085995j],
...     [ 0.53396245-0.10177564j,  0.76279558-0.35024096j]
... ])
>>> decomp = qml.transforms.xyx_decomposition(U, 0, return_global_phase=True)
>>> decomp = [gate.parameters[0] for gate in decomp]
>>> rx(decomp[0]) @ ry(decomp[1]) @ rx(decomp[2]) * decomp[3]
tensor([[-0.28829348-0.78829734j,  0.30364367+0.45085995j],
        [ 0.53396244-0.10177564j,  0.76279558-0.35024096j]], requires_grad=True)
>>> rx(decomp[2]) @ ry(decomp[1]) @ rx(decomp[0]) * decomp[3]
tensor([[ 0.76279558-0.35024096j,  0.30364367+0.45085995j],
        [ 0.53396244-0.10177564j, -0.28829348-0.78829734j]], requires_grad=True)

We can see that the correct output is obtained for rx(decomp[0]) @ ry(decomp[1]) @ rx(decomp[2]) * decomp[3]. However, this doesn't correspond to what the documentation states, since this amounts to applying rx(decomp[2]) first, not rx(decomp[0]). Thus, if this issue is valid, should it be the documentation that should be corrected, or the code?

[Newtech66]

I think would be better to have xyx_decomposition return the gates in the correct order because it enables directly calling the function to add the decomposition to a circuit. Also, it would be great to make zyz_decomposition work the same way as xyx_decomposition, then we can return the global phase in one_qubit_decomposition as well.

Upon further inspection, the following depend on the current implementation of zyz_decomposition:

• ctrl_decomp_zyz (but it is easily fixed)
• unitary_to_rot

[Newtech66]

Should the decomposition functions support batched unitaries? The documentation does not make it clear but there is explicit support for batched unitaries in zyz_decomposition (not so for xyx_decomposition).

[Newtech66]

Am I making a mistake or does xyx_decomposition not give the same answer as qml.Rot for the following case?

>>> import pennylane as qml
>>> import pennylane.transforms as tr
>>> from pennylane import numpy as np
>>> from functools import reduce
>>> U=qml.Rot(2.3,2.3,2.3,wires=0).matrix()
>>> obtained_gates = tr.xyx_decomposition(U,0,return_global_phase=True)
>>> print(obtained_gates)
[RX(array(-2.81949056), wires=[0]), RY(array(3.04714574), wires=[0]), RX(array(-3.46369475), wires=[0]), (1+6.3335721407800714e-18j)*(Identity(wires=[0]))]
>>> print(U)
[[-0.27216539-0.30461121j -0.91276394-0.j        ]
 [ 0.91276394+0.j         -0.27216539+0.30461121j]]
>>> obtained_mat = obtained_gates[0].matrix()@obtained_gates[1].matrix()@obtained_gates[2].matrix()@obtained_gates[3].matrix()
>>> print(obtained_mat)
[[-0.04720591-3.16208337e-01j -0.94751459-1.66434329e-18j]
 [ 0.94751459+1.03379607e-17j -0.04720591+3.16208337e-01j]]

[Newtech66]

It seems there is indeed a bug in the code for xyx_decomposition, it calculates $\theta$ incorrectly. The equivalent Qiskit function returns the correct value of $\theta$ for this unitary.

>>> import qiskit as qt
>>> Q = qt.quantum_info.OneQubitEulerDecomposer(basis='XYX')
>>> Q.angles_and_phase(U)
[2.590307346392943, -2.8194905595576767, 2.8194905595576767, -3.141592653589793]
>>> def rx(x):
...     return np.array([[np.cos(x/2), -np.sin(x/2)*1j], [-np.sin(x/2)*1j, np.cos(x/2)]])
...

>>> def ry(x):
...     return np.array([[np.cos(x/2), -np.sin(x/2)], [np.sin(x/2), np.cos(x/2)]])
...
>>> rx(-2.8194905595576767)@ry(2.590307346392943)@rx(2.8194905595576767)*np.exp(-3.141592653589793j)
tensor([[-0.27216539-3.04611214e-01j, -0.91276394-1.18720238e-16j],
        [ 0.91276394+1.04842450e-16j, -0.27216539+3.04611214e-01j]], requires_grad=True)

[Newtech66]

The bug arises because of a faulty conditional: if math.allclose(lams_plus_phis, 0):

It neglects that lams_plus_phis could also be some multiple of $\pi$. Accounting for this makes xyx_decomposition return a correct answer. I'll fix it in my PR.

[Newtech66]

Please don't forget to assign me to this issue for UnitaryHack.

[trbromley]

Hi @Newtech66! Congrats on getting this merged and thank you for being a contributor to PennyLane! 🚀

As part of our 0.31 release, we're planning to promote UnitaryHack contributions to PennyLane on social media. This is completely optional, but would you like to be mentioned in our social media posts? Specifically, we'd need either your Twitter profile, your LinkedIn profile, or your name.

[Newtech66]

@trbromley Sure! My name and LinkedIn can be found on my GitHub profile.