How to NLsolve with sparsity but no jacobian?

Specific Domains Optimization (Mathematical)
, , ,
1

From the NLsolve-docs i know that one can pass a jacobian to nlsolve and make the solver treating the jacobian as sparse.
But what can I do when I don’t want to write down a jacobian of a PDE discretisation, and only know its sparsity pattern?
I tried to pass this pattern to DiffEqDiffTools.forwarddiff_color_jacobian!. The resulting jacobian than to nlsolve. And surprise… specifing a jacobian that way ended up with nlsolve taking for ever to get a wrong / NaN solution.

Here is what I have tired.

using SparseDiffTools
using DiffEqDiffTools
using NLsolve

f! = (dx, x) ->  rhs!(dx, x, some_parameters, 0)  # Discretized PDE

jac = some_sparsity_pattern_ of_J_from_f!
colors = matrix_colors(jac)
 j = (jac, x) -> forwarddiff_color_jacobian!(jac, f!, x, colorvec = colors)

x0 = some_initial_guess 
#then I call
nlsolve(f, j, x0)  # slow
nlsolve(f , x0)  # fast and works (but want to have it faster)

Is there any way exploit the sparsity of the jacobian to speed up nlsolve ?

2

I’ve just figured out that the problem caused by the colorvec argument. Without it works but no gain in terms of speed is observed.

1 Like
3

Using just the colorvec will build the compressed Jacobian. You need to make sure you decompress it, and for that you need to pass the sparsity pattern to sparsity. So in most use cases you’ll want to pass both colorvec and sparsity.

2 Likes
4

Ahhh … The inplace variant does not provide such an sparsity argument ( but the out of place does). I looks like that nlsolve simply offers dense target variables for the Jacobians, so that the sparsity information is lost, this I totally have missed…thanks so far. OK tomorrow I will elaboate on this…

5

Yes it does. It’s the argument called sparsity.

6

Indeed, thanks for that hint again. Now it works quite nicely

f! = (dx, x) ->  rhs!(dx, x, some_parameters, 0)  # Discretized PDE

jac_sp = some_sparsity_pattern_ of_J_from_f!
colors = matrix_colors(jac_sp)
j!(jac,x)=   forwarddiff_color_jacobian!(jac,f, x;
                                        colorvec = colors,
                                        sparsity = jac_sp)

x0 = some_initial_guess 

nlsolve(f!, x0)  # 13 s
nlsolve(f!, j!, x0)  # 11 s

# ask nlsolve to manipulate sparse matrices instead of full ones
df =  OnceDifferentiable(f!, j!, u0,copy(u0), jac_sp)
nlsolve(df, u0) # 2 s

that’s quite nice (650% faster in my case)… @ChrisRackauckas : might it also be worth to wrap something like that i to SSRootfind?

2 Likes
7

Indeed it would be good to do that in SSRootfind. But also, I think NLsolve should just be enhanced so that the colorvec = colors argument is exposed and then sparsity should always be set to the sparse Jacobian if you’ve passed that to NLsolve. Then this would be fairly automatic.

Then, Optim.jl should also take the colorvec and sparsity arguments, and everything with sparse handling would be cohesive.

@pkofod see the example above.

1 Like
8

For future reference on what was going on here with the compression/decompression, consult the following notes:

https://mitmath.github.io/18337/lecture9/stiff_odes

1 Like