Add schmidt_decomposition function to quantum_info

[diemilio]

Summary

Addresses #10087

Details and comments

Added function to quantum_info.states.utils that performs the Schmidt decomposition of a bipartite pure state. Following questions/items still need to be addressed:

• Currently if qargs = None, function raises an error. Would it be better if it returned a copy of the state? Only issue is that output format will be different to that of a properly factorized state. Same if user passes qargs equal to all subsystem positions.
• Sometimes, the svd function in numpy.linalg has rounding errors and returns very small singular values in the order of 1e-16. Are there any qiskit guidelines regarding the tolerance at which these errors should be rounded? In other words, to what precision are two numbers (let's say, prob amplitudes) considered to be equal?
• Is having the output be a list of tuples containing the Schmidt coefficients and their corresponding subsystem vectors OK? Is there some other preferred format, like for example three separate lists, one for coeffs, one for states of A, one for states of B?
• Regarding tests, where is the correct place to add them? Also, I have a few ideas of what to test, but please let me know if there is anything specific that is required.
• Add import to __init__.py under quantum_info and under quantum_info.states if needed.
• Will add release notes when items above are completed.

[qiskit-bot]

Thank you for opening a new pull request.

Before your PR can be merged it will first need to pass continuous integration tests and be reviewed. Sometimes the review process can be slow, so please be patient.

While you're waiting, please feel free to review other open PRs. While only a subset of people are authorized to approve pull requests for merging, everyone is encouraged to review open pull requests. Doing reviews helps reduce the burden on the core team and helps make the project's code better for everyone.

One or more of the the following people are requested to review this:

• @Qiskit/terra-core
• @ikkoham

[diemilio]

Hello @ikkoham, thanks so much for the interest in having this integrated into qiskit. Please let me know your thoughts regarding the items above, or any questions/concerns you might have about the code.

Thanks!

[coveralls]

Pull Request Test Coverage Report for Build 5315103594

• 24 of 27 (88.89%) changed or added relevant lines in 2 files are covered.
• 245 unchanged lines in 30 files lost coverage.
• Overall coverage increased (+0.04%) to 85.964%

---

Changes Missing Coverage	Covered Lines	Changed/Added Lines	%
qiskit/quantum_info/states/utils.py	23	26	88.46%

Files with Coverage Reduction	New Missed Lines	%
crates/accelerate/src/error_map.rs	1	81.63%
crates/accelerate/src/sabre_swap/neighbor_table.rs	1	87.18%
crates/accelerate/src/sabre_swap/sabre_dag.rs	1	97.14%
crates/qasm2/src/expr.rs	1	93.76%
crates/qasm2/src/lex.rs	1	91.39%
crates/qasm2/src/lib.rs	1	92.21%
qiskit/circuit/instruction.py	1	95.69%
qiskit/dagcircuit/dagcircuit.py	1	89.59%
qiskit/transpiler/passmanager_config.py	1	97.26%
crates/accelerate/src/edge_collections.rs	2	54.55%

Totals
Change from base Build 5094569308:	0.04%
Covered Lines:	71494
Relevant Lines:	83167

---

💛 - Coveralls

[ewinston]

Just a quick response,

• tests for quantum_info are under python/test/quantum_info
• toloerance defaults can be found from qiskit.quantum_info.operators.predicates import ATOL_DEFAULT, RTOL_DEFAULT

[diemilio]

Thanks @ewinston.

I have added some tests; hope these are enough. Do we test for cases when QiskitError is triggered?

[ikkoham]

Thank you very much. schmidt_decomposition is useful.

Currently if qargs = None, function raises an error. Would it be better if it returned a copy of the state? Only issue is that output format will be different to that of a properly factorized state. Same if user passes qargs equal to all subsystem positions.

Raising an error is fine. I can't think of a use case for this, but if some user has a request, we can extend it.

Is having the output be a list of tuples containing the Schmidt coefficients and their corresponding subsystem vectors OK? Is there some other preferred format, like for example three separate lists, one for coeffs, one for states of A, one for states of B?

Returning list of tuple is nice idea. We can guarantee that the lists are the same length.

[ikkoham]

Could you test to see if you are getting the correct decomposition, not just the return to the original.

[diemilio]

Hi @ikkoham. I had already tested that the correct decomposition was being returned, but working on the implementation of these tests made me realize something. As you probably know, the SVD is a unique decomposition but only up to the preservation of the sign parity of the singular vectors. So for example, the state |++⟩ has these two decompositions which are equally valid:
λ|u⟩|v⟩ = 1.0 (1.0 |+⟩) (1.0 |+⟩)
λ|u⟩|v⟩ = 1.0 (-1.0 |+⟩) (-1.0 |+⟩)

The global phase of the reconstructed total state remains the same, but the individual "global" phases of each singular vector can be different as long as they preserve the parity.

The sign selection for the singular vectors depends on the underlying algorithm used for the SVD which, in this case is numpy.linalg.svd. This really isn't a big problem, except that for simple cases of separable states (like the ones you can construct from the Statevector.from_label method), I am always getting singular vectors with the negative sign, so for the self.assertEqual to work, I am going to have to premultiply some terms by -1 in the test function. Again, there is nothing incorrect about this as long as it is done keeping in mind that the parity must be preserved, my worry is that if someone is going thru the code later on, it might be confusing.

What do you think?

[diemilio]

I went ahead and implemented the test as described above, but let me know if you would like to see it implemented in a different way. Thanks!

[ikkoham]

Thanks. AFAIK, there are no rules or conventions on that point. I'm fine with this implementation. (I'd like to hear other's opinion.)

[diemilio]

Hi @ikkoham. I found this very helpful report on the sign ambiguity of the SVD:
https://www.osti.gov/servlets/purl/920802

What they suggest to deal with the sign ambiguity, is:
"...in order to identify the sign of a singular vector, it is suggested that it be similar to the sign of the majority of vectors it is representing."

They propose a function to do this, but the issue is that they only consider real-valued matrices. I am not familiar with how SVD algorithms work, so I don't know if in the case of complex-valued data, the SVD could also return singular vectors with arbitrary phases as long long as the overall sign of the sum term is preserved. If this is the case, the implementation of this function can get complicated.

[ewinston]

Instead of checking for equality between individual statevectors perhaps it would be better to just verify that the result does indeed recompose and satisfy all properties as described in your documentation of a Schmidt decomposition.

The question of numerical stability of the decomposition due to the sign issue could be treated in a separate test since it addresses a somewhat distinct issue. I'm not yet familiar with what the best convention might be at the moment so I'd be ok with any convention, which might include unhandled, and clearly documenting the choice.

Also in your current test function it would be better to split it up into several tests along the lines of each commented block you already have.

[diemilio]

Hi @ewinston. Thanks for the feedback.

1. Currently there are tests for both: 1) check for equality between individual singular vectors and singular values, and 2) test to verify that the result recomposes. Is your suggestion to get rid of 1)?
2. Regarding sign ambiguity testing, do you mean we should just keep the test we have right now but move it under a different test function? Or are you referring to having a different set of tests altogether?
3. Regarding sign convention selection, I am not quite sure what algorithm is used in the numpy implementation of the svd, but in the paper I shared it is mentioned that if it is a Lanczos-based method, there is a random component to it, so the selected sign might not be guaranteed. The best solution would be to implement the sign correction process suggested in that paper, but again, here we're just trying to correct for some "global" phase factor which is less critical in quantum computing compared to other applications where the sign might play an important role.
4. I will split the test functions just as you suggested. I will probably wait to hear back from you regarding items 1 and 2 before doing so, so I can make the final changes all together.

[diemilio]

@ikkoham, thank you for your feedback. I will address these items asap.

[diemilio]

Hello @ikkoham. I completed all remaining tasks. Please let me know if you have any further comments or questions.

[ikkoham]

Looks good to me. It would be nice to have a second reviewer if possible.

[diemilio]

Thank you so much @ikkoham. I agree it will be nice to have a second review. Do you know anyone that can help?

[ShellyGarion]

@diemilio - thank you very much for your important contribution to Qiskit! The code is well-written and documented.

I have a few minor suggestions:

1. In your tests there are only examples where all the Schmidt coefficients are equal, can you perhaps add an example where they are not equal? Perhaps check some larger pseudo-random statevecors?

2. When running the tests, I see that the Schmidt coefficient is in fact 0.9999999999999999 and not 1.0. I wonder if we should add a parameter to round it. For example, in https://qiskit.org/documentation/stubs/qiskit.quantum_info.Statevector.probabilities.html there is a parameter decimals. @ikkoham - what do you think?

[diemilio]

Hi @ShellyGarion,

Thanks so much for your feedback!

Regarding 1. Yes!, that's a great idea I will add an example were the Schmidt Coefficients are not equal. When it comes to checking to pseudo-random statevectors, we can definitely add a test to check that the original state can be recovered, but I think that testing for the individual Schmidt coeffs and vectors will not be possible. Is this OK if we just check for the reconstructed state?

Regarding 2. Adding that feature should be fairly simple. Just wondering if it could cause confusion if not used correctly since once you approximate to certain decimal place, the original state cannot be recovered in its exact form. Also, should we also apply the approximation to the statevectors, or just the coefficients?

Let me know what you think and I will make the changes. I also have to wait hear back from @ewinston regarding these four items for the tests.

Thanks again for your help!

[diemilio]

Hi @ShellyGarion! I went ahead and added a test where different probability amplitudes are used (see test function test_schmidt_decomposition_entangled)

I also separated the tests into different functions as recommended by @ewinston.

Still waiting to hear back from @ikkoham regarding adding a parameter to round Schmidt coefficients. An easy solution is to actually just round to the value being used to discriminate if the coeffs should be kept or or not (ATOL_DEFAULT).

[ikkoham]

Thanks. I always wonder if we should add truncate/round/decimals. The reason is that it is not obvious what the use cases are and what the best rounding method is. I think it would be nice to add it later when it is needed, since it would be possible to add it without breaking the API.

[diemilio]

Hi @ikkoham. Do you have any other comments/changes regarding this PR? Happy to take care of anything that might be missing.

[ikkoham]

Thank you! LGTM