density matrix measurements, comp basis measurements, circuit run

[sarsid]

Checklist
Thank you for contributing to QuTiP! Please make sure you have finished the following tasks before opening the PR.

• Contributions to qutip should follow the pep8 style.
  You can use pycodestyle to check your code automatically
• Please add tests to cover your changes if applicable.
• If the behavior of the code has changed or new feature has been added, please also update the documentation and the notebook. Feel free to ask if you are not sure.

Delete this checklist after you have completed all the tasks. If you have not finished them all, you can also open a Draft Pull Request to let the others know this on-going work and keep this checklist in the PR description.

Description
This PR is supposed to mark on-going additions to a new class qutip.qip.circuit.Measurement as well as an elementary qutip.qip.circuit.Measurement.QubitCircuit.run() function to exactly simulate quantum circuits. The major functions added right now are:

• Measurement.density_measurement(self, measurement_ops, state) : This adds generalized state/density matrix measurements given a list of observables. Still needs addition of checks to see if the list of observables is valid
• Measurement.measurement_comp_basis(self, state) : This adds measurement statistics for specific qubit measurements in the computational basis and also returns the collapsed states (without changing the dimensions).
• QubitCircuit.run(self, state, cbits): Given a ket input to the circuit, applies gates and measurements from the circuit and returns the resultant ket.

Right now, the measurement module is very bare-bones and splintered between the previous code
from @hodgestar and the new code added by me. The eigenstate based measurements are slightly different than the measurements usually used on circuits (or I might be missing something). Ofc, in the case of measurements on all qubits, they amount to the same result. I have not added tests yet but I have added a basic quantum teleportation notebook. It will be nice to get feedback on testing methods as well as what kind of notebooks the community is interested in. Also, any other feature requests on this module would also be appreciated.

PS: I have also added measurement code from Simon's PR here so maybe we can either add him as a contributor to this PR or figure out something else.

Update:

1. Addressing some concerns about the PR containing multiple diverging pieces of code, I think classical registers, measurement on the circuit and the classically-controlled gates are best grouped together as they all rely on each other in some ways. The valid concern is that of the Measurement class encapsulating all projective style measurements. I also want to point out that it was raised by @BoxiLi that we may want to move the Gate class and Measurement class completely to qutip.qip.operations to which I agree but I believe that should be addressed in a clean-up PR with naming changes suggested by @ajgpitch.

Related issues or PRs
Adds to #1090

Changelog
Adds partial and complete measurements to state vectors/density matrices and a basic run function

[BoxiLi]

This might be helpful in adding @hodgestar as co-authors in a commit. https://help.github.com/en/github/committing-changes-to-your-project/creating-a-commit-with-multiple-authors

[hodgestar]

Hi @sarsid! I'd like to recommend that we land #1090 separately. It's a stand alone piece of work and coupling PRs together into big PRs where it can be avoided generally makes things harder to review and land.

I can show you how to merge that branch into this one if you need to it make progress before #1090 lands.

I would also recommend making the formatting changes you have made (I assume using black or your editor?) in a separate PR to this one. Lumping them into one PR also makes reviewing harder since one doesn't know which changes are intentional or meaningful.

I'm also happy to help with this if you need!

I've started work on the docs for #1090, so hopefully it can land this week.

[sarsid]

I think I agree that the two PRs should be kept separate. Initially, I wanted to make them the same since I was trying to use the code in #1090 but right now the two don't interact much.

Wrt the formatting changes, they were put inadvertently and I'll make sure to revert them but thanks for pointing them out. Let me know what you think of the code in this PR and if you think there's some redundancy between the two (I think there should be but I can't quite put a finger on it).

[ajgpitch]

I have now merged #1090.
@sarsid @hodgestar good to see you two working together on this.
We would like as much as possible that is generic to measurement of quantum systems to be in the qutip core, with then only what is specific to QIP circuits being in the qip subpackage. Perhaps you two can continue to develop on this?
We can have qubit specific functions in the qutip core if that is necessary to preserve execution speed.

[sarsid]

I think I am pretty happy now with all my functions and variable names (for now). I would say we should leave the variable name changing and moving gates and measurement to separate files for a different PR (it can just be a restructuring one !). I was wondering how to check if it passes current tests (I am not sure what changes were there on Travis).

[BoxiLi]

You can either run it on your local machine by the command pytest or turn it from a draft PR into a ready-for-reveiw PR, which will trigger the Travis CI tests.

[ajgpitch]

I will reiterate @BoxiLi suggestion. If you find that tests pass locally, then switch from draft to review ready. You can also set the 'ready for review' label with your group membership permissions. I will take a look over the code now, but will not fully review it until you mark it as ready for such. Would be best if tests were in place before this.

Take a look at the codeclimate issues. They don't have to fixed, they are just suggestions.

[ajgpitch]

I am happy with the general style of the code.
It seems like this PR implements classical control as well as measurement, or am I missing something? As far as is possible PR should be self contained to one additional feature on an existing sub-package.
I am also not sure about having the Measurement class in the circuit.py file. We should try to keep to smaller code module files where possible.

[sarsid]

Some Updates:

1. I have added some tests for both the Measurement (in qutip/tests/test_circuit_measurement.py) portion and the circuit run (qutip/tests/test_qubitcircuit.py) portions of the code. There is a slight hiccup with the Python 3.8 run of one of the tests (test_run_statistics). I am not sure why it is only failing in the Python 3.8 case.
2. I have converted all the tests in qutip/tests/test_qubitcircuit.py to pytest style tests. Maybe, @jakelishman can comment if further changes are required.
3. Finally, since measurement_density and measurement_ket allow for POVM style measurements, it might be better to move those functions to qutip/measurement.py but I am not sure. Maybe @hodgestar can comment on that.

[BoxiLi]

I am not sure why it is only failing in the Python 3.8 case.

Well, it's a statistical average, so it can deviate. I don't think it has anything to do with Python 3.8. But probably it can be avoided? I guess its the same question as #1268

I'm thinking, what about adding a targets parameters to QubitCircuit.run_statistics. Usually, people don't interested in the full output state of the circuit because many of them are ancillary qubits. They can use targets to specify what is the qubits they want to look at. In the function, we can use ptrace(state, targets) to trace out the ancilla.

[sarsid]

I am not sure why it is only failing in the Python 3.8 case.

Well, it's a statistical average, so it can deviate. I don't think it has anything to do with Python 3.8. But probably it can be avoided? I guess its the same question as #1268

I'm thinking, what about adding a targets parameters to QubitCircuit.run_statistics. Usually, people don't interested in the full output state of the circuit because many of them are ancillary qubits. They can use targets to specify what is the qubits they want to look at. In the function, we can use ptrace(state, targets) to trace out the ancilla.

Seems like a good idea, maybe we can add it as optional parameter to both QubitCircuit.run_statistics and QubitCircuit.run.

Also, regarding Python 3.8, I was saying because both times it only failed in that test run, also never failed in any of my runs. Should I make tolerance higher ?

[BoxiLi]

This seems to be quite vulnerable to floating-point error and global phase. The ideal way would be to check if the fidelity of two states is smaller than some threshold.

[sarsid]

Good idea, will do ! I guess the state output can just be the first instance (I'm assuming it doesn't matter much what the small deviations are)

[BoxiLi]

Also, regarding Python 3.8, I was saying because both times it only failed in that test run, also never failed in any of my runs. Should I make tolerance higher ?

Ok.... Let me try restart the CI test once

[jakelishman]

Consider what you're doing from a statistical sense here. If everything works correctly, then you're effectively trying a two-tail hypothesis test that your test binomial distribution has a probability of 0.25, given you made n observations of it. The analytic success rate of your test is then

import scipy.stats
def success(n, p, tol):
    dist = scipy.stats.binom(n, p)
    return dist.cdf(n * (p+tol)) - dist.cdf(n * (p-tol))

where success(4096, 0.25, 0.02) is 99.67%. Since there are 5 independent runs on Travis, that's equivalent to a failure rate of 1.6% on every single CI run, which is several orders of magnitude too high. It's probably just luck that it was Python 3.8 on both occasions - it's a one-in-five chance that in two failing CI runs, it was the same setup both times.

You have to consider what is most appropriate to fix the test, and what you can reasonably test in ~1 second of runtime. You don't want to increase the tolerance too much, because then you can get a lot of false positives even if something is broken.

If you can bump the number of runs up to 100,000 and set the tolerance at 0.01, you'll have a per-test failure rate of ~3e-13, which is more like what we'd want. If you keep it at 4096, the tolerance should be more like 0.05 (which is pretty big tbh).

[sarsid]

Consider what you're doing from a statistical sense here. If everything works correctly, then you're effectively trying a two-tail hypothesis test that your test binomial distribution has a probability of 0.25, given you made n observations of it. The analytic success rate of your test is then

import scipy.stats
def success(n, p, tol):
    dist = scipy.stats.binom(n, p)
    return dist.cdf(n * (p+tol)) - dist.cdf(n * (p-tol))

where success(4096, 0.25, 0.02) is 99.67%. Since there are 5 independent runs on Travis, that's equivalent to a failure rate of 1.6% on every single CI run, which is several orders of magnitude too high. It's probably just luck that it was Python 3.8 on both occasions - it's a one-in-five chance that in two failing CI runs, it was the same setup both times.

You have to consider what is most appropriate to fix the test, and what you can reasonably test in ~1 second of runtime. You don't want to increase the tolerance too much, because then you can get a lot of false positives even if something is broken.

If you can bump the number of runs up to 100,000 and set the tolerance at 0.01, you'll have a per-test failure rate of ~3e-13, which is more like what we'd want. If you keep it at 4096, the tolerance should be more like 0.05 (which is pretty big tbh).

Thanks for the excellent analysis. I don't think it's feasible to do 100,000 given the current. efficiency. It did give me incentive to make it somewhat faster. What is a good guideline for the maximum time a test can take. It seems like I can maybe do 150 runs in ~ 1s. In any case, it doesnt seem feasible to do even 4096 claims without taking quite a bit of time. Any ideas on how to structure the test differently ?

[BoxiLi]

Hmmm... 4000/150~25 second. That's not very short actually. The total time for a whole round of qutip CI test costs 15min~20min on Travis.

Just ideas. Tests here seem to be two-folded:

• Test classical controlled gates in a circuit
• Test the measurement functions

So maybe we can split it. The teleportation circuit test can be done without measurement. Just check the final state tracing out the ancillary qubits. The test for measurement can probably be done for single or two Qobj along without circuit. But run_circuit_statistics... Probably we then only need to check if the number of elements in the result is correct after some 50 runs?

On a different matter, I'm wondering if running the statistics takes so long, for such a small teleportation circuit, is it still advantageous in any case? Since we are doing simulation and have the full quantum state, one can actually calculate all 4 possible final states, classically mix them into a density matrix with the corresponding measurement probability and calculate the exact statistics distribution. I doubt that will be slower than 25s, although coding will be harder I suppose.

[sarsid]

Hmmm... 4000/15025 second. That's not very short actually. The total time for a whole round of qutip CI test costs 15min20min on Travis.

Just ideas. Tests here seem to be two-folded:

• Test classical controlled gates in a circuit
• Test the measurement functions

So maybe we can split it. The teleportation circuit test can be done without measurement. Just check the final state tracing out the ancillary qubits. The test for measurement can probably be done for single or two Qobj along without circuit. But run_circuit_statistics... Probably we then only need to check if the number of elements in the result is correct after some 50 runs?

On a different matter, I'm wondering if running the statistics takes so long, for such a small teleportation circuit, is it still advantageous in any case? Since we are doing simulation and have the full quantum state, one can actually calculate all 4 possible final states, classically mix them into a density matrix with the corresponding measurement probability and calculate the exact statistics distribution. I doubt that will be slower than 25s, although coding will be harder I suppose.

I added the teleportation circuit because it seemed like a simple enough example to test both classically controlled gates and measurements. I have some separate (non-circuit based) examples in the other file. Maybe the run_statistics test can be not on the teleportation circuit and something even simpler ?

Re: the idea for run_statistics, do you mean tracking the various probability elements during each measurement (along with the state) ? Seems like a decent idea !

[sarsid]

Updates:

• I have switched the order of state and op/ops.
• There is now an optional targets arguments on each of the measurement functions to enable applying ops on particular
  qubits although this is inherently delicate, especially when dealing with non-qubit type dimensions. Right now, it's restricted to qubit style expansions (in terms of dimensions). This indicated in the docstring. I will add some more information on how to use this in the user guide but I think it works best for more simpler use cases !

[sarsid]

Again, I am not sure why some unrelated tests are failing !

[BoxiLi]

These failing tests are different from those random failing we've been seeing before. I also have this in my scheduler PR, but it is completely unrelated to the PR. It starts to appear yesterday, but nothing was merged in the last three weeks. Does anyone have a clue?

The error seems to come from the core data part @jakelishman @Ericgig. Scipy made a release 4 days ago and we are using one of the private attributes.

[jakelishman]

scipy 1.5 has changed some of its private attributes for matrix multiplication (to a much better name!), and so matrix-matrix multiplication is completely broken for us at the moment: see scipy/scipy@53fac7a.

There's also been some changes to how Hermitian eigenvalues and vectors are calculated (which actually is good news for us in general), which may has a bit of a knock-on for some of the zheevr tests, which is going to be a bit of a nuisance to fix. To avoid polluting this PR with (any more) off-topic discussion, I've opened #1299.

The new data layer types will fix the matrix multiplication issue permanently, because we'll not be duplicating/reusing large tracts of scipy private code - we have our own Cython versions that operate faster on more optimised types. In the meantime, #1298 is a patch to catch the renames.

[jakelishman]

@sarsid, @BoxiLi: the scipy issues are all fixed and merged to master.

[BoxiLi]

@sarsid, @BoxiLi: the scipy issues are all fixed and merged to master.

Many thanks!

[sarsid]

@quantshah @BoxiLi @hodgestar I have wrapped up the two measurement functions in the measure and measurement_statistics function which calls the required functions based on whether ops is a list of Qobjs or a Qobj. Do we want the individual functions (especially POVM) to be more granular? Moreover, the docstring for the wrapper function is really clunky (given the case-wise output types). Any ideas on how to make that simpler ?(one idea is to not specify return types and only have parameters given it's bound to be clunky!)

ps. I'll fix the tests all at once after we decide on the api

[BoxiLi]

The only concern I have for this API is that now measure(zero_state, plus_state) and measure(zero_state, [plus_state]) will get completely different behavior. The first is observable measurement while the later the usual one. This might be confusing and was also my concern before. Did you come to some conclusion last Thursday? @quantshah @hodgestar @sarsid

[BoxiLi]

Please highlight the difference also here or in the function documentation:

• if ops is list, it is measurement povm or projective
• if ops is qobj, it is measurement observerable

Probably also give an extra emphasis in the rst doc.

[BoxiLi]

Same here

[BoxiLi]

Is this a good way to do it? @nathanshammah @ajgpitch @quantshah @hodgestar

[nathanshammah]

• please add a line descriptor under each variable, e.g., case (which is repeated twice on line 450 and 456).
• using :class:qutip.Qobj would make it clickable and take one to the API. Or even just :class:.Qobj.

[BoxiLi]

I suggest we merge it after the tests and the codeclimate issues are (partially) fixed. We can make an update PR for doc if necessary.

[sarsid]

Addressed all the above concerns, made measurement documentation to be as Alex indicated in email!

[BoxiLi]

I find it's good to go. If there is no other big suggestions, I'll merge. Nice work!

[nathanshammah]

Looks fine to me, great work! I added some suggestions that you @sarsid should be able to commit directly. You need to add ticks like "" and "" just after :class: and :func: to be picked up by Sphinx. This is a general feature of restructured text.

[sarsid]

Looks fine to me, great work! I added some suggestions that you @sarsid should be able to commit directly. You need to add ticks like "" and "" just after :class: and :func: to be picked up by Sphinx. This is a general feature of restructured text.

Thanks for the fixes!

[nathanshammah]

Great work @sarsid. (turns out that everytime I was writing `, it did not render here.)

[jakelishman]

If useful, you can escape backticks in markdown with a \, so you can write \` to get one. If you need a single backtick in an inline code block, you can use two backticks to start the code block, so to get ` you type `` ` ``.