Error fixes for #81

[cgeoga]

Hey @maltezfaria---

I got around to taking a look at what's going on here. The following code was hitting two edge cases:

using StaticArrays, HMatrices

pts = rand(SVector{2,Float64}, 1000)
km  = KernelMatrix((x,y)->Float64(x==y), pts, pts)
hk  = assemble_hmatrix(km; rtol=1e-10)
hkf = cholesky(Hermitian(hk))

1. The check in the ACA in the case where no suitable pivot used a getindex in the check to decide whether to throw a warning. But since that getindex throws an error, that wasn't working. For the moment I've just wrapped it in a try-catch, but I assume you'll have a more thoughtful suggestion about how to do that check.
2. The _mul_dense! routine was missing routines for adjoints, so I just added that to the branches on Adata. Maybe this one is actually a decent solution? But like above, very happy to do something else if you think it is more suitable.

Thanks again for the great package!

[maltezfaria]

I would rather not have a try/catch block in this performance critical part of the code. Why is it needed? Is there a way to handle the exceptional case without running into the exception (I don't remember what the exception was).

[cgeoga]

Hey @maltezfaria---this is bullet point 1 above. Both of those all statements will throw an error when you are doing a low-rank approximation of a MulLinearOp. This seems to come up when you are factorizing a matrix and you are re-compressing a matrix product/sum where one of the blocks is basically all zeros.

Maybe an alternative fix would be to have a type-level check? Like, something that does the all(...) checks when K isa KernelMatrix and doesn't attempt to do so when K isa MulLinearOp or something?

[cgeoga]

Just pushed up a commit to change to that if-else check.

Another thing to consider here, though, is what to do about those warnings. If you run the example I have above with the identity matrix, you get a warning for every single admissible block. I'm not really familiar with ways that the ACA can fail in the sense of failing to select a pivot when the approximation isn't machine-precision accurate, so maybe that is possible. But maybe some less noisy way of communicating that would be good.

[maltezfaria]

Would adding maxlog=1 to @warn resolve the "communication issue"? At some point I should look closer into the problem to understand the proper way to handle it, but since I won't have time right now perhaps we can settle for this hacky alternative?

[cgeoga]

Good idea! Done.

[cgeoga]

Ah, and one other comment here---would you mind merging #81 and tagging a release? At the moment, as I tinker with this in various codebases I am having issues with the cholesky methods not being found sometimes just because the versioning is a little complicated.

[maltezfaria]

LGTM, thanks!

Minor suggestions in case you want to fix it before the PR is merged (feel free to ignore them as well):

• minor version bump
• add maxlog to the @warn

[maltezfaria]

Do you want to append maxlog=1 here to avoid the message appearing multiple times? Maybe modify the message to say that further instances of this issue will be silenced?

[cgeoga]

Done, thanks!

[codecov]

Codecov Report

❌ Patch coverage is 50.00000% with 3 lines in your changes missing coverage. Please review.
✅ Project coverage is 76.43%. Comparing base (88a154f) to head (d654f9d).
⚠️ Report is 5 commits behind head on main.

Files with missing lines	Patch %	Lines
src/compressor.jl	50.00%	2 Missing ⚠️
src/multiplication.jl	50.00%	1 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main      #82      +/-   ##
==========================================
- Coverage   76.54%   76.43%   -0.11%     
==========================================
  Files          14       14              
  Lines        1475     1477       +2     
==========================================
  Hits         1129     1129              
- Misses        346      348       +2     

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

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

[cgeoga]

@maltezfaria, thanks again! One other comment here for potential future stuff. I wonder if an elegant way around all of these warnings would be to have a struct like

struct RankZeroMatrix
  n::Int
  m::Int
end

and then define a few special methods for MulLinearOp when R or P is a RankZeroMatrix or something? Then the ACA compressor could say that if it can't find a pivot on the first search it just returns a RankZeroMatrix of the appropriate size, and re-compressing MulLinearOps with RankZeroMatrix objects could also just immediately return a RankZeroMatrix. That would also save some potentially significant passes over totally degenerate blocks anyway.

Maybe I'm not entirely understanding the problem, but that's the first idea that comes to mind for me!

[maltezfaria]

@maltezfaria, thanks again! One other comment here for potential future stuff. I wonder if an elegant way around all of these warnings would be to have a struct like

struct RankZeroMatrix
  n::Int
  m::Int
end

and then define a few special methods for MulLinearOp when R or P is a RankZeroMatrix or something? Then the ACA compressor could say that if it can't find a pivot on the first search it just returns a RankZeroMatrix of the appropriate size, and re-compressing MulLinearOps with RankZeroMatrix objects could also just immediately return a RankZeroMatrix. That would also save some potentially significant passes over totally degenerate blocks anyway.

Maybe I'm not entirely understanding the problem, but that's the first idea that comes to mind for me!

I need to understand the issue better, but why not simply have an RkMatrix with A and B of sizes m x 0 and n x 0? These rank 0 matrices are simply the zero matrix, and some logic can be added to handle the special (dynamic) case where rank(R)=0.