Hessian sparsity detection

[shashi]

Hessian sparsity detection works by propagating linearity information. See linearity.jl for various definitions derived from DiffRules.jl (found by differentiating twice and checking for zero-ness on a number of primitive functions). See combine_terms for how every linearity case is handled.

Demo:

julia> hsparsity(x->x[1]+x[2], [1,2]) |> spy
     Sparsity Pattern
     ┌─────┐
   1 │⠀⠀⠀⠀⠀│ > 0
   2 │⠀⠀⠀⠀⠀│ < 0
     └─────┘
     1     2
     nz = 0

julia> hsparsity(x->x[1]*x[2], [1,2]) |> spy
     Sparsity Pattern
     ┌─────┐
   1 │⠀⠀⠀⠄⠀│ > 0
   2 │⠀⠐⠀⠀⠀│ < 0
     └─────┘
     1     2
     nz = 2

julia> hsparsity(x->x[1]^2*x[2], [1,2]) |> spy
     Sparsity Pattern
     ┌─────┐
   1 │⠀⠠⠀⠄⠀│ > 0
   2 │⠀⠐⠀⠀⠀│ < 0
     └─────┘
     1     2
     nz = 3

julia> hsparsity(x->x[1]/x[2], [1,2]) |> spy
     Sparsity Pattern
     ┌─────┐
   1 │⠀⠀⠀⠄⠀│ > 0
   2 │⠀⠐⠀⠂⠀│ < 0
     └─────┘
     1     2
     nz = 3

julia> hsparsity(x->x[1]\x[2], [1,2]) |> spy
     Sparsity Pattern
     ┌─────┐
   1 │⠀⠠⠀⠄⠀│ > 0
   2 │⠀⠐⠀⠀⠀│ < 0
     └─────┘
     1     2
     nz = 3


julia> hsparsity(x->sin(x[1]+x[2]), [1,2]) |> spy
     Sparsity Pattern
     ┌─────┐
   1 │⠀⠠⠀⠄⠀│ > 0
   2 │⠀⠐⠀⠂⠀│ < 0
     └─────┘
     1     2
     nz = 4

julia> hsparsity(x->sin(x[1]+x[2])*x[3], [1,2,3]) |> spy
     Sparsity Pattern
     ┌─────┐
   1 │⠀⠄⠠⠠⠀│ > 0
   3 │⠀⠅⠨⠈⠀│ < 0
     └─────┘
     1     3
     nz = 8

julia> hsparsity(x->sin(x[1]+x[2])*x[3], [1,2,3]) |> spy
     Sparsity Pattern
     ┌─────┐
   1 │⠀⠄⠠⠠⠀│ > 0
   3 │⠀⠅⠨⠈⠀│ < 0
     └─────┘
     1     3
     nz = 8
julia> hsparsity(x->sum(x[1:5])/sum(x[6:10]), rand(10)) |> spy
      Sparsity Pattern
      ┌─────┐
    1 │⠀⠀⢰⣶⣶│ > 0
      │⣤⣤⣼⣿⣿│ < 0
   10 │⠿⠿⠿⠿⠿│
      └─────┘
      1    10
      nz = 75

TODO:

• Add tests
• Throw when a tainted leaf-node function call is of unknown linearity
• Use take-all-branches cassette pass

[ChrisRackauckas]

comment that