The Julia Symmetric Matrix type stores the full matrix. I’m looking for a package that offers a way to store and to index into a symmetric matrix directly but only stores one triangle, but that is performant, operates as a matrix, and also allows me to change off-diagonal elements at will, e.g. setting X[2,3] = 4 is an acceptable operation.
What is my best bet here?
Background: storing Hessians for large-dimensional optimization problems.
Thanks!!
3 Likes
Google serves up StandardPacked.jl and RectangularFullPacked.jl. Note, as mentioned in the READMEs, that the memory savings come at the cost of some extra overhead and a slight loss of floating-point accuracy in many calculations.
3 Likes
Thanks a bunch @danielwe . I’ll have a look at these too. There’s a whole bunch of them, many dormant, now incompatible, or restrictive.
@danielwe Do you have a recommendation?
Nope, never used any of these.
Thanks. Guess I should have phrased my original question differently. I’m looking for which package to use out of the umpteen different ones that claim to do this.
I think the conventional wisdom here is that the factor-of-two savings of a packed format aren’t usually worth the performance hit that this incurs. Is a factor of two in memory really the difference between solvable and unsolvable for you?
Have you considered ways to avoid explicit storage of the whole n \times n Hessian? Note that you can compute Hessian–vector products much more efficiently without storing the matrix — in fact, you can save a factor of roughly n by forward-over-reverse-mode computation of Hessian–vector products (see section 8.4.1 of our matrix-calculus course notes)!
Given a fast Hessian–vector product, you can construct a low-rank approximation of the Hessian using e.g. randomized sketching methods in LowRankApprox.jl. (And L-BFGS optimization algorithms estimate low-rank Hessian approximations by a different method.)
6 Likes
In another direction, is it possible that your Hessian matrix is sparse?
1 Like
Thanks @stevengj . That’s helpful. No, it isn’t: waste not want not and all.
In fact, storing is what I’ll do since I suspect that they are not large by your standards.
I’m aware of the LBFGS algorithms with low-rank Hessians. I’m using Newton with Trust Regions which appears to be better behaved for the type of problem I’m working on.