Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Adds find_classical_subsystem #193

Merged
merged 6 commits into from
Sep 3, 2020
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
27 changes: 27 additions & 0 deletions thewalrus/quantum.py
Original file line number Diff line number Diff line change
Expand Up @@ -114,6 +114,7 @@
is_valid_cov
is_pure_cov
is_classical_cov
find_classical_subsystem
total_photon_num_dist_pure_state
gen_single_mode_dist
gen_multi_mode_dist
Expand Down Expand Up @@ -947,6 +948,7 @@ def is_classical_cov(cov, hbar=2, atol=1e-08):
Args:
cov (array): a covariance matrix
hbar (float): value of hbar in the uncertainty relation
atol (float): the absolute tolerance parameter used in `np.allclose`

Returns:
(bool): whether the given covariance matrix corresponds to a classical state
Expand All @@ -962,6 +964,31 @@ def is_classical_cov(cov, hbar=2, atol=1e-08):
return False


def find_classical_subsystem(cov, hbar=2, atol=1e-08):
"""Find the largest integer ``k`` so that subsystem in modes ``[0,1,...,k-1]`` is a classical state.


Args:
cov (array): a covariance matrix
hbar (float): value of hbar in the uncertainty relation
atol (float): the absolute tolerance parameter used when determining if the state is classical

Returns:
int: the largest k so that modes ``[0,1,...,k-1]`` are in a classical state.
"""
n, _ = cov.shape
nmodes = n // 2
if is_classical_cov(cov, hbar=hbar, atol=atol):
return nmodes
k = 0
mu = np.zeros(n)
is_classical = True
while is_classical:
_, Vk = reduced_gaussian(mu, cov, list(range(k + 1)))
is_classical = is_classical_cov(Vk, hbar=hbar, atol=atol)
k += 1
return k - 1

def gen_single_mode_dist(s, cutoff=50, N=1):
"""Generate the photon number distribution of :math:`N` identical single mode squeezed states.

Expand Down
38 changes: 37 additions & 1 deletion thewalrus/tests/test_quantum.py
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,7 @@
import numpy as np
from scipy.stats import poisson

from thewalrus.symplectic import rotation, squeezing, interferometer, two_mode_squeezing, beam_splitter, loss
from thewalrus.symplectic import rotation, squeezing, interferometer, two_mode_squeezing, beam_splitter, loss, expand

from thewalrus.random import random_covariance, random_interferometer

Expand All @@ -45,6 +45,7 @@
pure_state_amplitude,
state_vector,
is_classical_cov,
find_classical_subsystem,
total_photon_num_dist_pure_state,
gen_single_mode_dist,
fock_tensor,
Expand Down Expand Up @@ -1299,3 +1300,38 @@ def test_photon_number_expectation_two_mode_squeezed(r, phi, hbar):
for modes, expected in zip(mode_list, expected_vals):
val = photon_number_squared_expectation(means, cov, modes=modes, hbar=hbar)
assert np.allclose(val, expected)

@pytest.mark.parametrize("hbar", [0.5, 1, 2, 1.6])
def test_find_classical_susbsytem_tmsq(hbar):
"""Test that for a system of 2*n squeezed vacua the first n are in a classical state"""
n = 10
nmodes = 2 * n
S = np.identity(2 * nmodes)
for i in range(n):
S = expand(two_mode_squeezing(1.0, 0), [i, i + n], nmodes) @ S
cov = S @ S.T * hbar / 2
k = find_classical_subsystem(cov, hbar=hbar)
assert k == n


@pytest.mark.parametrize("hbar", [0.5, 1, 2, 1.6])
def test_find_classical_subsystem_product_sq(hbar):
"""Tests that for a product state of squeezed vacua the classical system size is 0"""
nmodes = 10
r = 1
vals = np.ones([nmodes]) * hbar / 2
cov = np.diag(np.concatenate([np.exp(2 * r) * vals, vals * np.exp(-2 * r)]))
k = find_classical_subsystem(cov, hbar=hbar)
assert k == 0


@pytest.mark.parametrize("hbar", [0.5, 1, 2, 1.6])
def test_find_classical_subsystem_thermal(hbar):
"""Tests that for a multimode thermal state the whole system is classical"""
nmodes = 20
diags = (2 * np.random.rand(nmodes) + 1) * hbar / 2
cov = np.diag(np.concatenate([diags, diags]))
O = interferometer(random_interferometer(nmodes))
cov = O @ cov @ O.T
k = find_classical_subsystem(cov, hbar=hbar)
assert k == nmodes