Module idea for helping calculate large sparse Jacobians

Specific Domains Optimization (Mathematical)
1

I’ve put together a minimal, tiny draft of a project idea I had. It’s about helping to work with large structured models, the kind of thing you might use g2o for (and I imagine similar things might exist for FEM). Basically we have formulas that need to be evaluated for many indices, and the idea is to help doing that and also calculating the Jacobian.

I would love to hear some comments! Check it out at: https://github.com/nlw0/SparseBlockJacobians.jl

3 Likes
2

Just FYI, in Julia master we now have AbstractSparseMatrixCSC interface (as an example usage, see also https://github.com/tkf/SparseXX.jl/pull/2). If you want to bundle block information in a custom sparse matrix type, it could be useful. (Although simple function-based API is the best approach sometimes.)

Note that this interface is still not “official.” The details would be finalized in:

https://github.com/JuliaLang/julia/pull/33054

2 Likes
3

Thanks, I definitely have to learn more about what is going on regarding sparse arrays. I’m not so sure how this connects to what I’m looking for, though. Right now I’m just filling up a SparseMatrixCSC and at most would like a CSR version of it. The essence of the project is that we have models defined by a few formulas that involve few variables and they are applied over long vectors to define the full model. The main questions in my head right now are: what’s a good, declarative way to represent the model and its functions, how to ensure an efficient Jacobian update, and how to interface nicely with AD tools. Perhaps my liberal use of the term “block” is misleading about the emphasis?

4

I hope there’s a bunch of iteration and coloring ideas in that interface. See https://github.com/JuliaDiffEq/ArrayInterface.jl

5

I saw you are passing block sizes to update_jacobian_values! and update_jacobian_indices!. I just thought it may make sense to embed this information in the matrix so that you only need to pass functions. Of course, as I said, a function-based API that manipulates SparseMatrixCSC directly totally makes sense as well.

6

It’s just a dead simple CSC matrix. I think the only thing that is relevant to ArrayInterface is that you can write something like has_sparsestruct(x::Type{<:AbstractSparseMatrixCSC}) instead of has_sparsestruct(x::Type{<:SparseMatrixCSC}).