Hi! I’m new to using `NonlinearSolve.jl`

and trying to understand how to provide an analytical jacobian to the problem. I realize that this is unnecessary if we have autodiff in the mix since that will be enough - but in cases where I have say registered functions (and corresponding registered jacobians) I want to understand how to explicitly supply the jacobian to the problem.

Here’s what I have as an example:

```
# Starting with a simple example
function f(u, p)
@. u*u - p
end
u0 = [1.0, 1.0, 1.0]
p = [2.0, 1.0, 4.0]
prob = NonlinearProblem(f, u0, p)
```

```
julia> sol = solve(prob, NewtonRaphson())
retcode: Success
u: 3-element Vector{Float64}:
1.4142135623730951
1.0
2.000000000000002
```

Now if I try to give it the jacobian and setup a `NonlinearFunction`

(was not very sure how else to setup a `NonlinearProblem`

with a jacobian than doing this):

```
df(u,p) = 2.0 .* u
h = NonlinearFunction(f, jac = df)
prob = NonlinearProblem(h, u0, p)
```

If I try solving this problem now I get this:

```
julia> sol = solve(prob, NewtonRaphson())
retcode: Unstable
u: 3-element Vector{Float64}:
1.0
1.0
1.0
```

And looking at this I tried specifying to the `NewtonRaphson`

algorithm `concrete_jac = true`

but that still does not seem to help.

```
julia> sol = solve(prob, NewtonRaphson( ;concrete_jac = true))
retcode: Unstable
u: 3-element Vector{Float64}:
1.0
1.0
1.0
```

All this is using the latest `NonlinearSolve`

version:

```
[8913a72c] NonlinearSolve v3.9.1
```

Any thoughts or suggestions would be very appreciated!