Bloch-Messiah returns incorrect results #728
Description
Before posting a bug report
- I have searched exisisting GitHub issues to make sure the issue does not already exist.
Expected behavior
Consider
import numpy as np
from thewalrus.random import random_symplectic
from thewalrus.quantum import is_symplectic
from strawberryfields.decompositions import bloch_messiah
r = np.array([1,1,1,1,2,2,0,0,0,0,0,0])
Sdiag= np.diag(np.concatenate([np.exp(r), np.exp(-r)]))
num_modes = len(r)
O1 = random_symplectic(num_modes, passive=True)
O2 = random_symplectic(num_modes, passive=True)
S = O1 @ Sdiag @ O2
assert is_symplectic(S)
# Do the Bloch-Messiah from SF
from strawberryfields.decompositions import bloch_messiah
O1, S, O2 = bloch_messiah(S)
assert is_symplectic(Sdiag)
assert is_symplectic(O1)
assert is_symplectic(O2)Actual behavior
The returned matrices O1 and O2 in the example above are not symplectic.
Reproduces how often
Always
System information
Strawberry Fields: a Python library for continuous-variable quantum circuits.
Copyright 2018-2020 Xanadu Quantum Technologies Inc.
Python version: 3.9.7
Platform info:
Installation path:
Strawberry Fields version: 0.22.0-dev
Numpy version: 1.21.5
Scipy version: 1.8.0
SymPy version: 1.10
NetworkX version: 2.6.3
The Walrus version: 0.20.0-dev
Blackbird version: 0.4.0
XCC version: 0.1.2
TensorFlow version: 2.6.2Source code
No response
Tracebacks
No response
Additional information
A fixed Bloch-Messiah can be found in https://github.com/polyquantique/thewalrus/pull/2