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

[WIP] Bosonic backend support for threshold measurement #604

Closed
wants to merge 15 commits into from

Conversation

JeyKelly
Copy link

@JeyKelly JeyKelly commented Jul 1, 2021

Context: Currently, the Bosonic backend does not support the MeasureFock and MeasureThreshold measurement operations. This pull request aims to add support for it. Ref : issue #571

Description of the Change: Adding support for measure_threshold and measure_fock

Benefits:

Possible Drawbacks: N/A

Related GitHub Issues:

@josh146
Copy link
Member

josh146 commented Jul 2, 2021

Thanks @JeyKelly for this PR! Feel free to ask any questions here, and you can comment directly on the code itself 🙂

@JeyKelly
Copy link
Author

JeyKelly commented Jul 2, 2021

@josh146 I've checked the result of the test and it seems like the test was designed to raise a non-implemented error when calling the MeasureThreshold operation (~strawberryfields\tests\backend\test_bosonic_backend.py line 664). Assuming the Fock state case works and post-measurement is updated properly, should I update the pull request to include removing the line of code and adding two extra modes to test_measurement verifying the proper sample shapes/states after applying MeasureFock and MeasureThreshold?

Comment on lines 845 to 859
##measure photon number in every gaussian mode
for i in range(len(cov)):
cov_i = cov[i]
reduced_cov_i = cov_i[np.ix_(modes_idxs, modes_idxs)]
tor_i = torontonian_sample_state(reduced_cov_i, shots)
measurements[i] = tor_i[0][0]

index = np.array([np.random.choice(list(range(len(measurements))),1,p=np.abs(weights))])
samples = measurements[index]
if measurements[index] == 1:
print("wip : updating the state in the case there is a click")
else:
print("wip : updating the state in the case there is no click")
return samples

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

dear @JeyKelly : my first suggestion would be to use the method fidelity_vacuum. This method will tell you the probability p0 that the detector does not click. The probability that it does click is given by 1-p0. You can then flip a coin with this bias (using np.random.choice) and then based on that one will have to update the vectors of coefficients, means, and covariance matrices. Maybe I'd suggest you start by writing the logic deciding whether to click or not. I will try to give you more detailed info on how to update the state based on the two different outcomes.

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@nquesada This is definitely a more efficient solution, thank you. It works properly except in the case of vacuum where I used is_vacuum to circumvent the rounding error (p0>1). Is there a similarly simple solution for the Fock measurement or is sampling each gaussian mode individually like I put in the pull request a decent approach? It still has the issue of using non-normalized weights, so in the case of a pure Fock state, a measurement can't be done using the current definition.

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hi @JeyKelly --- As promised here are some notes explaining the algorithm:

bosonic_threshold.pdf

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Also, this update rule was implemented in Fortran in this repo here https://github.com/XanaduAI/torontonian-sampling/blob/master/src/torontonian_samples.f90

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hi @JeyKelly --- As promised here are some notes explaining the algorithm:

bosonic_threshold.pdf

I took a look at this, thanks for explaining it. I was wondering though: wouldn't the rewritten cov matrix \sigma^{(l)} have (2M - 2) for matrices dimensions excluding the last mode in \sigma_{(B)} (and similarly for \vec{r}_{(A)})? I'm asking because I have a hard time understanding why the cov matrices have a (M - 2) x (M - 2) sizes post no-click measurements.

The way I interpret this is we're cutting out the mode completely post-measurement, and so we end up with modified matrices for the remaining states - however, I don't quite understand why in the case where we have a click, the remaining cov/weights/means matrices have sizes which reflect this extra mode that we have measured? Initially, I thought we were replacing the mode with vacuum and adjusting the remaining ones to reflect this depending on the measurement outcome, but as far as I understand it's not what is done here. Do you have an idea where I'm not understanding properly?

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Also, this update rule was implemented in Fortran in this repo here https://github.com/XanaduAI/torontonian-sampling/blob/master/src/torontonian_samples.f90

@nquesada I am having trouble interfacing this in Python; it's telling me I'm missing the dependencies (use kinds, vars, structures) in the torontonian_samples module. Would this ultimately be translated in python code to be implemented in Strawberryfields? If I can manage to make it run on my end I would gladly rewrite it in a python file that can be used directly in the measure_threshold method.

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hi @JeyKelly --- I am not suggesting you use that Fortran code verbatim. I suggested it because it can give you and idea of how to implement the update rules. Also, I made typo, in the notes, the updated matrices are of size 2M-2 x 2M-2. Sorry about this!

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

also to answer your question @JeyKelly about why they are 2M-2 x 2M-2, it is because they only describe the modes which were not measured. Since optical measurements are typically demolition measurements, the measured mode is typically reset to vacuum, as you said. This can be done afterwards, just be adding 2 extra rows/cols to the covariance matrix in the positions corresponding to the measured modes (just with hbar/2 on the diagonal and zero everywhere else). And adding 2 extra zeros to the updated 2M-2 length displacement vector.

@JeyKelly
Copy link
Author

@nquesada I have revisited the measure_threshold method in the way you suggested. I have ran into a few issues when testing, and this is why I didn't update the test file yet. If I ran the circuit with a Fock() state parameter of 7, the circuit returns nan values for the weights (hence can't even run without a measurement), and if I put more, the fidelity of vacuum is bigger than one. This doesn't seem to be related to these modifications, but I thought I would mention it before creating an automatic test for the measurement.

That being said, the measure_fock method still remains and I'm not sure on how to approach updating the state (and even if the measurement is done properly). Do you have a similar piece theory that would help achieve this one?

@nquesada
Copy link
Collaborator

Hi @JeyKelly --- The PR is looking really nice. For tests I'd suggest a few:

  • Prepare the vacuum and check that when you use threshold measurements you never get a click.
  • Prepare the odd cat states (with p=1) and also verify that you never get click.
  • This somewhat trivial, but check that the shapes of the covariances and vectors of means are correct when you get a click/no click.

@JeyKelly
Copy link
Author

Hi @JeyKelly --- The PR is looking really nice. For tests I'd suggest a few:

* Prepare the vacuum and check that when you use threshold measurements you never get a click.

* Prepare the odd cat states (with `p=1`) and also verify that you never get click.

* This somewhat trivial, but check that the shapes of the covariances and vectors of means are correct when you get a click/no click.

@nquesada Hi, thanks for the response. Just to be sure, do you mean for automatic tests?

@nquesada
Copy link
Collaborator

Hi @JeyKelly --- Yes, I mean automatic tests. You can start by using the backend independent tests here:

https://github.com/XanaduAI/strawberryfields/blob/master/tests/backend/test_threshold_measurement.py

To test the bosonic backend it should be as easy as adding "bosonic:" in this line:

@pytest.mark.backends("gaussian")

@nquesada
Copy link
Collaborator

Then you could add tests specific to the bosonic backend in here https://github.com/XanaduAI/strawberryfields/blob/master/tests/bosonic_files/test_backend_bosoniccircuit.py

This could be the ones I described before involving cat states and the ones checking the shapes of the updated coefficients, vectors of means and covariance matrices are correct.

@JeyKelly
Copy link
Author

JeyKelly commented Jul 13, 2021

To test the bosonic backend it should be as easy as adding "bosonic:" in this line:

@pytest.mark.backends("gaussian")

Hi Again @nquesada, by using the test test_two_mode_squeezed_measurements in this file, it made me realize that running the circuit with two identical squeezed state (on two uninteracting modes) and measuring them gives me a different result. It seems that the self.fidelity_vacuum() is the culprit. Multiple calls to the function in one program increments its value when the state gets updated post-measurement (namely the self.weights) on uninteracting modes, when there is "no click" ("clicks" measurements decrease the weights so it doesn't seem to run into problems). I'm not sure if the problem lies in the use of it in the code I submitted here or if the function itself has a flaw, but this prevents running the test because the vacuum fidelity goes above one in this case. Like I've mentioned previously, the fidelity is also problematic (above one) when calculated on a Fock state of more than 7 photons, even on a clean version of strawberryfields.

@nquesada
Copy link
Collaborator

Hi @JeyKelly --- That is because two-mode squeezing is different from applying single-mode squeezing in two modes. As for your second question, you are correct that Fock states have a number of numerical issues. For the moment I suggest you set them aside (both in terms of preparing them or measuring them) and focus in threshold detection. Let's worry about Fock measurements in a future PR.

@JeyKelly
Copy link
Author

JeyKelly commented Jul 13, 2021

Hi @JeyKelly --- That is because two-mode squeezing is different from applying single-mode squeezing in two modes.

@nquesada Yes, but I ran into fidelity problems by running the single-mode squeezing in two modes while trying to see where the problem was for two-mode squeezing (sorry, it might have been superfluous to mention it). Vacuum fidelity of the remaining mode increases (above 1) when there is a "no click" event in the first single-mode squeezed state.

@nquesada
Copy link
Collaborator

Can you share a minimal non-working example of the issue you are finding with the two single-mode squeezed states @JeyKelly ?

@JeyKelly
Copy link
Author

JeyKelly commented Jul 13, 2021

Can you share a minimal non-working example of the issue you are finding with the two single-mode squeezed states @JeyKelly ?

import strawberryfields as sf

hbar = 2
nmodes = 2

prog_cat_fock = sf.Program(nmodes)
with prog_cat_fock.context as q:
    sf.ops.Catstate(alpha=0.1)   | q[0]    
    sf.ops.Catstate(alpha=0.1)   | q[1]         
    sf.ops.MeasureThreshold()    | q[0]   
    sf.ops.MeasureThreshold()    | q[1]   
 
eng = sf.Engine("bosonic", backend_options={"hbar": hbar})
result = eng.run(prog_cat_fock)

@nquesada This will do the trick, or any combination of state close to unity vacuum fidelity (for example, by replacing both states with sf.ops.Squeezed(r=0.1)). When the weights are adjusted here post "no-click" with the c' = c/p0 line, vacuum_fidelity exceeds 1 when calculated for the q[1] mode.

@nquesada nquesada changed the title [WIP] Bosonic backend support for threshold and Fock measurement [WIP] Bosonic backend support for threshold measurement Aug 3, 2021
@codecov
Copy link

codecov bot commented Aug 6, 2021

Codecov Report

Merging #604 (da7cbd7) into master (45309dc) will not change coverage.
The diff coverage is n/a.

@@           Coverage Diff           @@
##           master     #604   +/-   ##
=======================================
  Coverage   98.56%   98.56%           
=======================================
  Files          77       77           
  Lines        8904     8904           
=======================================
  Hits         8776     8776           
  Misses        128      128           

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 45309dc...da7cbd7. Read the comment docs.

@nquesada
Copy link
Collaborator

nquesada commented Aug 10, 2021

@elib20 : I take the point that it should be equivalent to heterodyne post selected on the \alpha = 0 outcome. But this implementation is still correct. I guess it is hard to see how one is equivalent to the other since heterodyne calls generaldyne. Is your point that more tests should be added? Or that I should use reuse the functionality from heterodyne for the no click case here?

@nquesada nquesada closed this Aug 23, 2021
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
enhancement New feature or request WIP
Projects
None yet
Development

Successfully merging this pull request may close these issues.

Add photon number and threshold detection to the Bosonic backend
4 participants