Fix compatibility issues with SciPy 1.14

[jakelishman]

Summary

The main change here is that SciPy 1.14 on macOS now uses the system Accelerate rather than a bundled OpenBLAS by default. These have different characteristics for several LAPACK drivers, which caused numerical instability in our test suite. Fundamentally, these problems existed before; it was always possible to switch out the BLAS/LAPACK implementation that SciPy used, but in practice, the vast majority of users (and our CI) use the system defaults.

The modification to Operator.power to shift the branch cut was suggested by Lev. As a side-effect of how it's implemented, it fixes an issue with Operator.power on non-unitary matrices, which Sasha had been looking at.

The modification to the Choi-to-Kraus conversion reverts back from a Schur-based decomposition to an eigh-based one. This was noted in a code comment that it was causing instability, and tracing the git history through to fdd5603 (gh-3884) suggests that this was around the time of Scipy 1.1 to 1.3, with OpenBLASes around 0.3.6. The bundled version of OpenBLAS with SciPy 1.13 was 0.3.27, so several releases on. If eigh problems re-appear, we can look at changing the decomposition back to something else, but the current eigh seems to be more stable than the Schur form for a wider range of inputs now.

The test changes to the Kraus and Stinespring modules are to cope with the two operators only being defined up to a global phase, which the test previously did not account for. The conversion to Kraus-operator form happens to work fine with OpenBLAS, but caused global-phase differences on macOS Accelerate. A previous version of this commit attempted to revert the Choi-to-Kraus conversion back to using eigh instead of the Schur decomposition, but the eigh instabilities noted in fdd5603 (gh-3884) (the time of Scipy 1.1 to 1.3, with OpenBLASes around 0.3.6) remain with Scipy 1.13/1.14 and OpenBLAS 0.3.27.

Details and comments

Let's see what happens in the full test matrix, especially Windows (which iirc is usually the source of eigh problems). It worked with both SciPy 1.13 and 1.14 on my Intel Mac.

edit: it remained unstable, so it's reverted.

Fix #12655
Fix #13305
Fix #13307

This should hopefully unblock #13309.

This obsoletes the bugfix of #13319 as a side-effect of solving #13305, but the additional method / keyword argument to power there is still likely very useful (though will need adapting to use the same technique as here).

[qiskit-bot]

One or more of the following people are relevant to this code:

• @Qiskit/terra-core

[coveralls]

Pull Request Test Coverage Report for Build 11502063338

Details

• 10 of 10 (100.0%) changed or added relevant lines in 2 files are covered.
• 37 unchanged lines in 6 files lost coverage.
• Overall coverage decreased (-0.04%) to 88.65%

---

Files with Coverage Reduction	New Missed Lines	%
crates/accelerate/src/two_qubit_decompose.rs	1	92.09%
crates/qasm2/src/expr.rs	1	94.02%
crates/qasm2/src/lex.rs	2	93.73%
qiskit/transpiler/passes/synthesis/unitary_synthesis.py	2	61.59%
qiskit/synthesis/two_qubit/xx_decompose/decomposer.py	7	90.77%
crates/qasm2/src/parse.rs	24	96.23%

Totals
Change from base Build 11481266566:	-0.04%
Covered Lines:	74873
Relevant Lines:	84459

---

💛 - Coveralls

[levbishop]

It might be nice to explicitly test eigenvalues barely either side of the intended branch cut, but perhaps that would be a pain to write?

[jakelishman]

I could adapt the exact same test, I think? We'd just need to take care not to couple the test to the particular choice of the rotation epsilon.

[levbishop]

Yeah, it might get confusing though with the different epsilons involved. There's eps_bc that defines the intended position of the branch cut, eps_fp thats potential numerical imprecision during the test, and then eps_test that's how close the test comes either side of the branch cut. I guess you want eps_test > eps_fp and then check eigenvalues at eps_bc±eps_test away from the axis.

[jakelishman]

I was hoping to make eps_bc an implementation detail and specifically avoid testing anything about its value, but actually, maybe we make it a keyword arg on Operator.power to allow the user to move it around? Then it's easier to test as well.

[levbishop]

That seems in line with how we deal with other tolerance-ish things. I don't have a strong opinion.

[jakelishman]

Seems like eigh remains unstable on macOS. I'll switch that Choi-to-Kraus conversion over to something else.

[jakelishman]

Ok, I updated this to revert the superoperator conversion back to the Schur form (with a couple of comments to explain the validity of it), and to really fix their tests underlying problems, which were that they were failing to account for the operators only being defined up to a global phase. The Kraus representation would also potentially need to account for degeneracies in an eigensubspace of the Choi matrix, but that's not relevant for the problematic test here.

I exposed a branch_cut_rotation parameter to Operator.power, so the user can choose where the branch cut happens to stabilise their own roots for the types of operators they typically care about.

[levbishop]

Looks good