I am interested in differentiating things that are along the lines of what is described in section 2 of these notes where the linear systems are sparse. For several years, I have occasionally been checking out Julia’s AD ecosystem to see if there is a package that can do this sort of differentiation automatically. I haven’t yet had success. The code below shows a couple of my feeble attempts to get Zygote to do this for a very simple linear system. Am I doing something wrong, or is Zygote not ready for this yet? Is there another package that can automatically differentiate the code below? I am especially interested in the `g2`

example which is more representative of real problems that I am working with.

```
import LinearAlgebra
import SparseArrays
import Zygote
function g1(p)
n = 3
A = LinearAlgebra.diagm(-1=>fill(-p, n - 1), 0=>fill(2p, n), 1=>fill(-p, n - 1))
b = ones(n)
return (A \ b)[1]
end
@show g1(1.0)
@show g1'(1.0)#Mutating arrays is not supported
function g2(p::T) where {T}
n = 3
I = Int[]
J = Int[]
V = T[]
for i = 1:n
if i > 1
push!(I, i)
push!(J, i - 1)
push!(V, -p)
end
push!(I, i)
push!(J, i)
push!(V, 2p)
if i < n
push!(I, i)
push!(J, i + 1)
push!(V, -p)
end
end
b = ones(n)
return (SparseArrays.sparse(I, J, V) \ b)[1]
end
@show g2(1.0)
@show g2'(1.0)#segfault
```