[BUG] Precision error in two-qubit decomposition of qml.QubitUnitary about pennylane HOT 10 OPEN

Hi @abhishekabhishek, thank you for opening this bug. What precision were you looking for? Is it more a matter of being aware of the precision beforehand? Or do you need a precision of 1e-9 or similar?

from pennylane.

Hi @CatalinaAlbornoz, I was expecting a precision of at least 1e-8 since that's what np.allclose uses for matrix equivalence checking by default. This is what led to some tests failing for my use case above. However, since I don't have a strict requirement on the precision, 1e-7 suffices but I think it is worth mentioning in the documentation in case other users do have a strict requirement.

from pennylane.

Ah good catch @abhishekabhishek ! Yes I agree that mentioning it in the docs is important. Thanks for the feedback!

from pennylane.

The problem seems to lie in the optimization we are using to reduce the number of CNOTs. At least the code above works out if one deactivates this optimization and always decomposes to the most general circuits containing 3 CNOTs...

from pennylane.

I took a closer look at this (as the original author, I feel responsible for the issues here! 😅). I tried with/without the small perturbation in precision that was added in #5448, but this was not the problem, and there is precision loss in both cases.

I think what it boils down to is, the matrix used here is a special case (controlled-phase gate) whose decomposition requires only 2 CNOTs. You can check that the following circuit identity holds, up to a global phase:

Applying some simple circuit identities on the single-qubit rotations, we can see that the two_qubit_decomposition function produces a very convoluted version of this with many extra rotations, some of which are 0. The large amount of extra angles, each with some small loss of precision, contributes to the overall loss of precision.

However, I'm a bit stumped about why it works with the 3 CNOT case which has just as many angles. It could be a matter of how the rotations are derived from the matrix products computed internally by functions like _extract_su2su2_prefactors; in the 2-CNOT case for this controlled phase, many of the matrix entries are very close to 0.

from pennylane.

Indeed, looking at all these extra rotations and setting the angles that are very close to $\pi$ or $3\pi$ to np.pi fixes the precision loss (to standard numpy precision). In particular, the middle angle of the single-qubit rotations into ZYZ seems to cause precision loss 🤔

edit: I was able to solve the issue using arccos of the diagonal entry rather than arcsin of the off-diagonal (in _zyz_get_rotation_angles). However, I'd expect that both will suffer from precision loss at extremal points, so this is not really a general fix.
Maybe we can come up with a generalized routine that smartly chooses how to extract this type of angle information from multiple (theoretically) redundant matrix entries, picking the numerically favourable one.

from pennylane.

Perhaps a simple check like

if qml.math.isclose(np.sin(theta), 0):
    theta = np.pi

in those angle extraction functions would be useful for "resetting" things that are close-to-multiples-of-pi to pi? It could be a "good first issue". (Though TBH I could have sworn I have implemented something like that in the past, but maybe that was in somewhere in The Ionizer)

from pennylane.

The above will no longer be differentiable at these special points (an issue I've been fighting with elsewhere (#5715) as well 😅 )

from pennylane.

Right, of course... 🤦‍♀️

from pennylane.

these kind of numerical-value-based checks/"fixes" can also sometimes lead to problems with jitting

from pennylane.