qml.lie_closure often raises LinAlgError: Singular matrix

[Qottmann]

Running into this problem a lot when using sums of paulis, wonder how we can make this more reliable and stable (ideally with no significant performance regression).

In instances where the code breaks, I often find myself with large sums of Paulis with large coefficients, which makes me wonder whether picking up this old suggestion by @vincentmr could help?

An example:

n = 4
ZZs = [Z(i) @ Z(i+1) for i in range(n-1)]
Xs = [X(i) for i in range(n)]

H_TFIM = qml.dot(np.ones(2*n-1, dtype=float), ZZs + Xs) # np.random.rand(2*n-1)
H_TFIM = H_TFIM.pauli_rep

S = [Z(i) for i in range(n)]
S = [op.pauli_rep for op in S]

vspace = qml.pauli.PauliVSpace(S, dtype=float)

while True:
    len_before = len(vspace.basis)
    i = 1
    for m, Om in enumerate(vspace.basis.copy()):
        com = H_TFIM.commutator(Om)
        com.simplify()

        if len(com) == 0:  # skip because operators commute
            continue

        # result is always purely imaginary
        # remove common factor 2 with Pauli commutators
        for pw, val in com.items():
            com[pw] = val.imag / 2.

        vspace.add(com, tol=1e-8)
    
    print("----------------\n")
    print(f"After {i} epochs, the following candidate support:")
    print(vspace.basis)
    print("----------------\n")

    i += 1
    
    if len_before == len(vspace):
        break

LinAlgError: Singular matrix

[josh146]

Thanks @Qottmann! Just double checking, does #6075 aim to resolve this? It wasn't clear to me from the PR (as the description is blank)

[Qottmann]

Thanks @Qottmann! Just double checking, does #6075 aim to resolve this? It wasn't clear to me from the PR (as the description is blank)

Yes! Though it is currently not a full fix and we are still not sure how best to resolve it. Linear independence in a stable fashion is hard apparently 😢

[dwierichs]

@josh146 Sorry for the blank description, I thought this was going to be a very quick fix that is not even worth triaging.
However, it turned out to be more complex after discussions on the PR.

We need to re-triage this issue, in my opinion, as it was found that neither @Qottmann nor I can guarantee to work on this soon.

[josh146]

Thanks! I guess one question; does the open PR improve the behaviour even if it doesn't resolve it fully?

[dwierichs]

Yep I think so, because a previously failing case is fixed by it (see test addition).
However, as this is about numerical precision and the like, there might be a new edge case lingering somewhere that stops working due to the PR, I suppose 🤔 Although it would surprise me.

[josh146]

Got it, thanks @dwierichs!