How to actually take hessian in Zygote

Specific Domains Numerics
1

for a simple function like this, why can’t

hessian(((a, b),) -> a[1]*b[2] + a[2]*b[1], [[2,1], [3,4]])

It throws BoundsError: attempt to access ForwardDiff.Dual{Nothing, Any, 2} at index [2]
biut for hessian(((a, b,c,d),) -> a*b + c*d, [2,1, 3,4]) it works. So for large array inputs, they should be unpacked and repacked ? that seems like a lot of work. Is there any other better solutions?

2

From the documentation of hessian, it seems it only acceps a single argument which must be an array or a number. Any other behavior, even if it appears to work like your second example, is undocumented and therefore not guaranteed.

One possible solution would be to use ComponentArrays.jl to encapsulate several inputs into one. However, this code throws another bug that I don’t understand:

using ComponentArrays, Zygote
x = ComponentVector(a = [2, 1], b = [3, 4])
f(x) = x.a[1] * x.b[2] + x.a[2] * x.b[1]

The Jacobian works but the Hessian fails:

julia> jacobian(f, x)
([4 3 1 2],)

julia> hessian(f, x)
ERROR: type Array has no field a
Stacktrace:
  [1] adjoint
    @ ~/.julia/packages/Zygote/jxHJc/src/lib/lib.jl:229 [inlined]
  [2] _pullback
    @ ~/.julia/packages/ZygoteRules/M4xmc/src/adjoint.jl:67 [inlined]
  [3] f
    @ ~/Downloads/quadgk.jl:4 [inlined]
  [4] _pullback(ctx::Zygote.Context{false}, f::typeof(f), args::Vector{ForwardDiff.Dual{Nothing, Int64, 4}})
    @ Zygote ~/.julia/packages/Zygote/jxHJc/src/compiler/interface2.jl:0
  [5] pullback(f::Function, cx::Zygote.Context{false}, args::Vector{ForwardDiff.Dual{Nothing, Int64, 4}})

Probably because ForwardDiff.jl discards the ComponentVector