Output truncation error for truncated decompositions

[Jutho]

No description provided.

[lkdvos]

Overall looks good!

I am slightly on the fence about whether we should compute the truncation error separately, or do that from inside the find_truncated functions. There's definitely some strategies where that information is either already available or more easily computable within these functions. In any case it probably makes sense to have a compute_truncerr function available for some of these strategies to use.

[Jutho]

Overall looks good!

I am slightly on the fence about whether we should compute the truncation error separately, or do that from inside the find_truncated functions. There's definitely some strategies where that information is either already available or more easily computable within these functions. In any case it probably makes sense to have a compute_truncerr function available for some of these strategies to use.

The reason I did not quite liked this in the Copilot PR is that find_truncated is also used for left_null. But there it doesn't make sense at all to do this. But then the number of output arguments of find_truncated needs to depend on the specific function it receives as first argument, which is not great.

[Jutho]

This is also ready from my side, if at least the current strategy is acceptable.

[lkdvos]

I rebased onto the current main branch and made some small updates to the docs, which I think I overlooked in the previous PR.

Additionally, I tried adding a compute_truncerr that does not destroy the input, and depending on the case at hand actually turns out to be more or less as performant. This is currently not being used, I just wanted to test the waters and see how you feel about this. I think that the accuracy loss in the fallback case is probably acceptable, am I correct in thinking this would mostly be an issue for heavy truncation, where norm(A)^2 - norm(Strunc)^2 is a difference of a large and small number? Given that we already have to include this for the decompositions that don't have access to the full spectrum, I don't think this is an unfair concession.

Let me know what you think though, I can try and add this in and properly test it, but if you think this is not worth the hassle I'm happy with just removing it and merging this with the input-destroying function.

[github-actions]

After the build completes, the updated documentation will be available here

[Jutho]

I think that the accuracy loss in the fallback case is probably acceptable, am I correct in thinking this would mostly be an issue for heavy truncation, where norm(A)^2 - norm(Strunc)^2 is a difference of a large and small number? Given that we already have to include this for the decompositions that don't have access to the full spectrum, I don't think this is an unfair concession.

No inaccuracies would appear in the case of very little truncation, where norm(A) and norm(Strunc) are both big, but almost equal.

Let me know what you think though, I can try and add this in and properly test it, but if you think this is not worth the hassle I'm happy with just removing it and merging this with the input-destroying function.

I am certainly ok with the non-destructive compute_truncerr implementations, so for me that is ok. If you keep this, can you also remove the copy in the pullback rules which then becomes unnecessary?

[lkdvos]

I updated the JET compat here as well, and used the non-destructive implementation everywhere.
As a continuation on the naming discussions, I renamed the function to truncation_error to avoid random abbreviations, and to avoid using compute_x, which I seem to remember you don't usually like that much. Happy to change it to anything else though.
I would expect that this is a useful function to overload externally as well, so I marked it as public.

[lkdvos]

It turns out that my initial approach of sorting the vector of indices in-place breaks the pullbacks: since we are first computing ind, then slicing with ind, and only then computing the truncation error I cannot sort!(ind) there.

I currently patched it by allocating a new ind when the inds aren't sorted, but that leads to another allocation. I can come up with some solutions, although I'm not sure which is best: either using the inaccurate fallback or change the findtruncated implementations to return sorted inds in the first place. Finally, I could also compute the truncation error within the truncate functions, therefore doing that before the slicing actually happens which means I can modify ind from within the truncation_error function.

Ideas and opinions would be very welcome.

[Jutho]

If a copy of ind is necessary, then I would just compute the equivalent of setdiff(eachindex(values), ind), possibly also specialized for ranges or BitVector. Or just take a copy of values, set copy_values[ind] = 0 and compute the norm of the result, i.e. my original approach (but possibly with an extra copy if you don't like the mutating part).

[lkdvos]

I'll just change it back to your implementation, it was worth a shot but I feel like it isn't working out too much, so your approach was probably best :)

[lkdvos]

I'm a bit confused about the failing test, I would have thought this now all is computing the same thing in more or less the same way

[codecov]

Codecov Report

❌ Patch coverage is 87.09677% with 4 lines in your changes missing coverage. Please review.

Files with missing lines	Patch %	Lines
ext/MatrixAlgebraKitChainRulesCoreExt.jl	75.00%	4 Missing ⚠️

Files with missing lines	Coverage Δ
src/implementations/eig.jl	99.02% <100.00%> (+<0.01%)	⬆️
src/implementations/eigh.jl	94.48% <100.00%> (+0.04%)	⬆️
src/implementations/svd.jl	92.36% <100.00%> (+2.22%)	⬆️
src/implementations/truncation.jl	93.42% <100.00%> (-0.95%)	⬇️
src/interface/eig.jl	85.71% <ø> (ø)
src/interface/eigh.jl	80.00% <ø> (ø)
src/interface/svd.jl	80.00% <ø> (ø)
src/interface/truncation.jl	68.18% <ø> (ø)
ext/MatrixAlgebraKitChainRulesCoreExt.jl	81.96% <75.00%> (-0.49%)	⬇️

... and 2 files with indirect coverage changes

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[Jutho]

I've modified the tests a bit. The rtol = sqrt(eps) specification for the kept values shouldn't do anything, as that is the default already. For the truncation error, it makes sense that we need an atol specification, as this is supposed to be about small values the small values that are thrown away with respect to the large values. If the tests pass, this is good to be merged for me.