KernelMatrices.jl - A software package for working with hierarchical matrices

Package Announcements
21

If you end up comparing against multiple implementations of the ACA, would you mind benching the one in KernelMatrices against it and let me know how things look? At the time that I released the package it definitely wasn’t competitive, but the package has been quietly improving since then, so I am somewhat hopeful that now it might hold its own against more sophisticated implementations, particularly if you give the ACA a KernelMatrix object instead of a dense Matrix, although the function accepts both.

Example calls to KernelMatrices’ ACA:

using KernelMatrices

# Actually low rank:
truU  = randn(1024, 100)
truV  = randn(1024, 100)
A1    = truU*truV'
(U,V) = KernelMatrices.ACA(A1, 1.0e-8) # size(U,2) \approx 100. 

# Should be full rank, but for numerical reasons isn't (kernel is analytic everywhere/"entire").
# compared to LowRankApprox.jl:
using LowRankApprox
rand_pts = [randn(2) for _ in 1:1024]
A2  = [exp(-norm(x-y)^2) for x in rand_pts, y in rand_pts]
for rnk in (50, 150, 250)
  (fU,fV) = KernelMatrices.ACA(A2, 1.0e-16, rnk) 
  (pU, pS, pV) = psvd(A2, rank=rnk)
  println("Rank $rnk:")
  println("pSVD opnorm difference: $(round(opnorm(A2-pU*Diagonal(pS)*pV'), digits=5))")
  println("ACA  opnorm difference: $(round(opnorm(A2-fU*fV'), digits=5))\n\n")
end
2 Likes
22

Just stumbled upon this and was wondering if this package is differentiable? Is there a frule/rrule for H-matrices operations that we may need instead of differentiating through the source code directly?

23

Hi @KapilKhanal — thanks for the note! I don’t have any rules in here, no. But also, this package was my first real Julia project that I worked on maybe six months out of my undergrad degree. For a variety of reasons (primarily losing faith in H-matrices as the right tool for scaling up GP computations), I have kind of walked away from this package and haven’t kept up with it (and I don’t think the code quality is particularly good).

To you and any other people who stumble on this thread, I would recommend HMatrices.jl for your H-matrix needs. @maltezfaria really did excellent work and that code will probably take good care of you. I’m not sure if there are any diff rules in there, but if any package would be worth spending the time to contribute something like that to it would be that one.

4 Likes