LSMR docs not added to solvers.md

[jpbrodrick89]

Should probably be added to the least squares solvers section with a reference in the iterative solvers but could be the other way round.

[patrick-kidger]

Yeah, I'm slightly wary of adding this in as at least for now it keeps failing in tests. :/ For example.

The failures seems to happen when doing a jvp with a singular matrix, which isn't hugely encouraging for a method that is supposed to handle singular matrices! :D Most test suite runs have one or two failing tests each time, but sometimes we get lucky and have zero, which is probably how this got in in the first place.

I haven't entirely made up my mind which way to go on this yet – remove it as something I'm not able to maintain / put in the time to try and debug this / etc. So at least for now it's intentionally left out of the documentation as only semi-public API.

(CC @PTNobel )

[f0uriest]

Is there any way to figure out what rng values are causing it to fail? I'd be interested in trying to debug.

One thought that occurs if it only happens on a jvp of a singular matrix, is if the derivative effectively changes the range of the matrix? Under the hood we're basically solving the normal equations:

$$ A^T A x = A^T b $$

When A isn't full rank this system is still compatible, bc $A^T b$ projects out the part of $b$ that lies in the nullspace of $A$, so there is still some $x$ that satisfies that equation. However the jvp looks like

$$A^T A dx = dA^T b + A^T db - dA^T A x - A^T dA x$$

Because of the dA terms I'm not sure the RHS of this still lies in the range of $A^T A$ if $dA$ has a different nullspace than $A$

[patrick-kidger]

I'd be interested in trying to debug.

Awesome! Thank you 🙏

Is there any way to figure out what rng values are causing it to fail?

Yup! At the start of each failing test the local variable getkey = GetKey(seed=411397021, ... is reported. You can manually re-run the test locally with the environment variable EQX_GETKEY_SEED=411397021 in order to set this value.

That said I don't know to what extent this ends up depending on hardware etc, as finickity floating-point things sometimes do.

In terms of changing the range, I don't think this should matter. In principle LSMR should return the pseudoinverse solution (I am told).