Hi Community,

I’m working on a program that requires solving A\B repeatedly, where A is a sparse matrix that looks like this (it’s from finite-difference of an 2d PDE).

I was wondering is there any method that can utilize its “banded-block-banded” structure as well as its sparsity to achieve higher speed? I tried to solve it as a `bandedblockbandedmatricies`

but it’s ten times slower than lu from umfpack. I guess it’s because it’s too sparse.

I’m not super satisfied with the current method because I’m converting this program from Fortran, which also uses umfpack but runs 5 times faster than my Julia version so I was wondering any chance there is some black magic that can help.