From 96ecdda87b5804816a2cca7f6176c613bb4a9d5b Mon Sep 17 00:00:00 2001 From: Jakob Mannstadt Date: Tue, 28 Jan 2025 20:01:38 +0100 Subject: [PATCH 1/3] changed the decomposition.triangular function to return (None, diag(localV), t_list. and added the If BS1 not None: to the interferometer_decompose method to match the decomposition with the interferometer function. I also changed the return in tests/frontend/test_decompositions.py for the triaungular decomposition to the new convention. Furthermore, I added a rudimental test for the triangular interferometer in test_interferometer_triangular. --- strawberryfields/decompositions.py | 4 +- strawberryfields/ops.py | 35 ++++++------ tests/frontend/test_decompositions.py | 4 +- .../test_interferometer_triangular.py | 54 +++++++++++++++++++ 4 files changed, 77 insertions(+), 20 deletions(-) create mode 100644 tests/frontend/test_interferometer_triangular.py diff --git a/strawberryfields/decompositions.py b/strawberryfields/decompositions.py index f4f3d340d..f2a6c79aa 100644 --- a/strawberryfields/decompositions.py +++ b/strawberryfields/decompositions.py @@ -612,7 +612,7 @@ def triangular(V, tol=1e-11): :math:`|VV^\dagger-I| \leq` tol Returns: - tuple[array]: returns a tuple of the form ``(tlist,np.diag(localV), None)`` + tuple[array]: returns a tuple of the form ``(None,np.diag(localV), tlist)`` where: * ``tlist``: list containing ``[n,m,theta,phi,n_size]`` of the T unitaries needed @@ -630,7 +630,7 @@ def triangular(V, tol=1e-11): tlist.append(nullT(nsize - j - 1, nsize - i - 2, localV)) localV = T(*tlist[-1]) @ localV - return list(reversed(tlist)), np.diag(localV), None + return None, np.diag(localV), tlist def M(n, sigma, delta, m): diff --git a/strawberryfields/ops.py b/strawberryfields/ops.py index 68ff22505..9276d0bdf 100644 --- a/strawberryfields/ops.py +++ b/strawberryfields/ops.py @@ -2664,26 +2664,29 @@ def _decompose(self, reg, **kwargs): decomp_fn = getattr(dec, mesh) BS1, R, BS2 = decomp_fn(self.p[0], tol=tol) - for n, m, theta, phi, _ in BS1: - theta = theta if np.abs(theta) >= _decomposition_tol else 0 - phi = phi if np.abs(phi) >= _decomposition_tol else 0 - if "symmetric" in mesh: - # Mach-Zehnder interferometers - cmds.append( - Command( - MZgate(np.mod(theta, 2 * np.pi), np.mod(phi, 2 * np.pi)), - (reg[n], reg[m]), + if BS1 is not None: + + for n, m, theta, phi, _ in BS1: + theta = theta if np.abs(theta) >= _decomposition_tol else 0 + phi = phi if np.abs(phi) >= _decomposition_tol else 0 + + if "symmetric" in mesh: + # Mach-Zehnder interferometers + cmds.append( + Command( + MZgate(np.mod(theta, 2 * np.pi), np.mod(phi, 2 * np.pi)), + (reg[n], reg[m]), + ) ) - ) - else: - # Clements style beamsplitters - if not (drop_identity and phi == 0): - cmds.append(Command(Rgate(phi), reg[n])) + else: + # Clements style beamsplitters + if not (drop_identity and phi == 0): + cmds.append(Command(Rgate(phi), reg[n])) - if not (drop_identity and theta == 0): - cmds.append(Command(BSgate(theta, 0), (reg[n], reg[m]))) + if not (drop_identity and theta == 0): + cmds.append(Command(BSgate(theta, 0), (reg[n], reg[m]))) for n, expphi in enumerate(R): # local phase shifts diff --git a/tests/frontend/test_decompositions.py b/tests/frontend/test_decompositions.py index 9df20787d..01632b2d0 100644 --- a/tests/frontend/test_decompositions.py +++ b/tests/frontend/test_decompositions.py @@ -351,7 +351,7 @@ def test_identity(self, tol): n = 20 U = np.identity(n) - tlist, diags, _ = dec.triangular(U) + _, diags, tlist = dec.triangular(U) qrec = np.diag(diags) @@ -373,7 +373,7 @@ def test_random_unitary(self, tol): n = 20 U = haar_measure(n) - tlist, diags, _ = dec.triangular(U) + _, diags, tlist = dec.triangular(U) qrec = np.diag(diags) diff --git a/tests/frontend/test_interferometer_triangular.py b/tests/frontend/test_interferometer_triangular.py new file mode 100644 index 000000000..8bb03b3ab --- /dev/null +++ b/tests/frontend/test_interferometer_triangular.py @@ -0,0 +1,54 @@ +import strawberryfields as sf +from strawberryfields.ops import * +import numpy as np + +M = 4 +# Create a random complex matrix +random_matrix = np.random.rand(M, M) + 1j * np.random.rand(M, M) + +# Perform QR decomposition to get a unitary matrix +q, r = np.linalg.qr(random_matrix) + +# Normalize to ensure unitary property (q should already be unitary, but we'll double-check) +unitary_matrix = q / np.linalg.norm(q, axis=0) + + +boson_sampling_triangular = sf.Program(4) +#initiate a starting state |1101> +with boson_sampling_triangular.context as q: + # prepare the input fock states + Fock(1) | q[0] + Fock(1) | q[1] + Vac | q[2] + Fock(1) | q[3] + + # apply the triangular interferometer + Interferometer(unitary_matrix, mesh = 'triangular', tol = 1e-10) | q + +#run the simulation +boson_sampling_triangular.compile(compiler="fock") +eng = sf.Engine(backend="fock", backend_options={"cutoff_dim": 5}) +results = eng.run(boson_sampling_triangular) +# extract the joint Fock probabilities +probs_triangular = results.state.all_fock_probs() + +#now try with mesh = 'rectangular' +boson_sampling_rectangular = sf.Program(4) + +with boson_sampling_rectangular.context as q: + # prepare the input fock states + Fock(1) | q[0] + Fock(1) | q[1] + Vac | q[2] + Fock(1) | q[3] + + # apply the rectangular interferometer + Interferometer(unitary_matrix, mesh = 'rectangular', tol = 1e-10) | q + +boson_sampling_rectangular.compile(compiler="fock") +eng = sf.Engine(backend="fock", backend_options={"cutoff_dim": 5}) +results = eng.run(boson_sampling_rectangular) +# extract the joint Fock probabilities +probs_rect = results.state.all_fock_probs() +#all probabilities should be the same here! Sp the output should be True +print(f"{np.allclose(probs_triangular, probs_rect)}") \ No newline at end of file From d50ab2c14faef345819f408444558616f29272f3 Mon Sep 17 00:00:00 2001 From: Jakob Mannstadt Date: Tue, 28 Jan 2025 20:46:55 +0100 Subject: [PATCH 2/3] changed TestTriangularDecompsition to match new convention --- tests/frontend/test_decompositions.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/frontend/test_decompositions.py b/tests/frontend/test_decompositions.py index 01632b2d0..3f81f53dc 100644 --- a/tests/frontend/test_decompositions.py +++ b/tests/frontend/test_decompositions.py @@ -355,7 +355,7 @@ def test_identity(self, tol): qrec = np.diag(diags) - for i in tlist: + for i in reversed(tlist): qrec = dec.Ti(*i) @ qrec assert np.allclose(U, qrec, atol=tol, rtol=0) From a7ae05f3ea807699abeba182725511d3c53e343e Mon Sep 17 00:00:00 2001 From: Jakob Mannstadt Date: Wed, 12 Feb 2025 10:23:09 +0100 Subject: [PATCH 3/3] changed the interfoermeter check to hopefully pass the test now --- .../test_interferometer_triangular.py | 44 +++++++++---------- 1 file changed, 21 insertions(+), 23 deletions(-) diff --git a/tests/frontend/test_interferometer_triangular.py b/tests/frontend/test_interferometer_triangular.py index 8bb03b3ab..58352fea1 100644 --- a/tests/frontend/test_interferometer_triangular.py +++ b/tests/frontend/test_interferometer_triangular.py @@ -2,30 +2,31 @@ from strawberryfields.ops import * import numpy as np -M = 4 +""" +Test to check that the triangular and rectangular meshes give the same results in the +interferometer opreation. +The tests creates a random M by M unitary matrix and applies it to a M-mode interferometer +with the rectangular and triangular meshes. The test checks that the output probabilities. +""" + +M = 5 # Create a random complex matrix random_matrix = np.random.rand(M, M) + 1j * np.random.rand(M, M) - # Perform QR decomposition to get a unitary matrix q, r = np.linalg.qr(random_matrix) - # Normalize to ensure unitary property (q should already be unitary, but we'll double-check) unitary_matrix = q / np.linalg.norm(q, axis=0) - -boson_sampling_triangular = sf.Program(4) -#initiate a starting state |1101> +boson_sampling_triangular = sf.Program(M) with boson_sampling_triangular.context as q: # prepare the input fock states - Fock(1) | q[0] - Fock(1) | q[1] - Vac | q[2] - Fock(1) | q[3] - - # apply the triangular interferometer + for i in range(M): + if i % 2 == 0: + Fock(1) | q[i] + else: + Vac | q[i] Interferometer(unitary_matrix, mesh = 'triangular', tol = 1e-10) | q -#run the simulation boson_sampling_triangular.compile(compiler="fock") eng = sf.Engine(backend="fock", backend_options={"cutoff_dim": 5}) results = eng.run(boson_sampling_triangular) @@ -33,16 +34,14 @@ probs_triangular = results.state.all_fock_probs() #now try with mesh = 'rectangular' -boson_sampling_rectangular = sf.Program(4) +boson_sampling_rectangular = sf.Program(M) with boson_sampling_rectangular.context as q: - # prepare the input fock states - Fock(1) | q[0] - Fock(1) | q[1] - Vac | q[2] - Fock(1) | q[3] - - # apply the rectangular interferometer + for i in range(M): + if i % 2 == 0: + Fock(1) | q[i] + else: + Vac | q[i] Interferometer(unitary_matrix, mesh = 'rectangular', tol = 1e-10) | q boson_sampling_rectangular.compile(compiler="fock") @@ -50,5 +49,4 @@ results = eng.run(boson_sampling_rectangular) # extract the joint Fock probabilities probs_rect = results.state.all_fock_probs() -#all probabilities should be the same here! Sp the output should be True -print(f"{np.allclose(probs_triangular, probs_rect)}") \ No newline at end of file +assert(np.allclose(probs_triangular, probs_rect)) \ No newline at end of file