# Hessian is not jacobian-compose-gradient in Symbolics.jl

**URL:** <https://discourse.julialang.org/t/hessian-is-not-jacobian-compose-gradient-in-symbolics-jl/128921>\
**Category:** General Usage\
**Tags:** question, symbolics, autodiff\
**Created:** [May 12, 2025, 3:02pm UTC](https://discourse.julialang.org/t/hessian-is-not-jacobian-compose-gradient-in-symbolics-jl/128921 "2025-05-12T15:02:28Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![WalterMadelim](https://avatars.discourse-cdn.com/v4/letter/w/3e96dc/32.png) [@WalterMadelim](https://discourse.julialang.org/u/WalterMadelim)\
**Post date:** [May 12, 2025, 3:02pm UTC](https://discourse.julialang.org/t/hessian-is-not-jacobian-compose-gradient-in-symbolics-jl/128921/1 "2025-05-12T15:02:28Z")

</div>

I wonder why is (the title).  
My results is like this. Did I use its API correctly?

```julia
julia> using Symbolics

julia> @variables p[1:10]
1-element Vector{Symbolics.Arr{Num, 1}}:
 p[1:10]

julia> f = bobj(e, p) # some scalar function
-1.5log(7 + p[4]) - 1.5log(5 + p[3]) - 1.5log(3 - p[9]) - 1.5log(3 - p[7]) - 1.5log(7 + p[1]) - 1.5log(8 - p[2]) - 1.5log(9 - p[4]) - 1.5log(9 - p[1]) - 1.5log(3 + p[8]) - 1.5log(6 + p[6]) - 1.5log(6 - p[3]) - 1.5log(6 - p[6]) - 1.5log(5 - p[8]) - 1.5log(5 + p[9]) - 1.5log(3 + p[7]) - 1.5log(5 + p[5]) - 1.5log(5 + p[2]) - 1.5log(5 - p[5]) - 1.5log(6 + p[10]) - 1.5log(6 - p[10]) - 0.08p[1]*p[6] + 0.09p[1]*p[7] + 0.07p[1]*p[8] - 0.02p[10]*p[3] - 0.07p[10]*p[4] - 0.02p[10]*p[5] - 0.03p[2]*p[7] - 0.07p[2]*p[8] - 0.04p[2]*p[9] + 0.02p[3]*p[8] + 0.06p[3]*p[9] + 0.08p[4]*p[9]

julia> Symbolics.jacobian(Symbolics.gradient(f, p), p)
10×10 OffsetArray(::Matrix{Num}, 1:10, 1:10) with eltype Num with indices 1:10×1:10:
 1.5 / ((7 + p[1])^2) + 1.5 / ((9 - p[1])^2) 0 … 0 0
                                           0 1.5 / ((8 - p[2])^2) + 1.5 / ((5 + p[2])^2) -0.04 0
                                           0 0 0.06 -0.02
                                           0 0 0.08 -0.07
                                           0 0 0 -0.02
                                       -0.08 0 … 0 0
                                        0.09 -0.03 0 0
                                        0.07 -0.07 0 0
                                           0 -0.04 1.5 / ((3 - p[9])^2) + 1.5 / ((5 + p[9])^2) 0
                                           0 0 0 1.5 / ((6 + p[10])^2) + 1.5 / ((6 - p[10])^2)

julia> Symbolics.hessian(f, p)
10×10 Matrix{Num}:
 0 0 0 0 0 0 0 0 0 0
 0 0 0 0 0 0 0 0 0 0
 0 0 0 0 0 0 0 0 0 0
 0 0 0 0 0 0 0 0 0 0
 0 0 0 0 0 0 0 0 0 0
 0 0 0 0 0 0 0 0 0 0
 0 0 0 0 0 0 0 0 0 0
 0 0 0 0 0 0 0 0 0 0
 0 0 0 0 0 0 0 0 0 0
 0 0 0 0 0 0 0 0 0 0

```

~~(If you need the source code, I’ll append then.)~~ It was already in the post (though a bit lengthy).

---

<div class="post-metadata">

**Author:** ![gdalle](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/gdalle/32/27854_2.png) [@gdalle](https://discourse.julialang.org/u/gdalle)\
**Post date:** [May 12, 2025, 4:44pm UTC](https://discourse.julialang.org/t/hessian-is-not-jacobian-compose-gradient-in-symbolics-jl/128921/2 "2025-05-12T16:44:53Z")

</div>

The source code would be helpful indeed

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [May 12, 2025, 8:28pm UTC](https://discourse.julialang.org/t/hessian-is-not-jacobian-compose-gradient-in-symbolics-jl/128921/3 "2025-05-12T20:28:21Z")

</div>

I think this is a bug with symbolic arrays. Check `Symbolics.hessian(f, collect(p))`, since `p` is not an array of symbolic values but a symbolic array, it might be hitting the wrong dispatch.

---

<div class="post-metadata">

**Author:** ![WalterMadelim](https://avatars.discourse-cdn.com/v4/letter/w/3e96dc/32.png) [@WalterMadelim](https://discourse.julialang.org/u/WalterMadelim)\
**Post date:** [May 12, 2025, 11:22pm UTC](https://discourse.julialang.org/t/hessian-is-not-jacobian-compose-gradient-in-symbolics-jl/128921/4 "2025-05-12T23:22:15Z")

</div>

It works. I’ve opened an [issue](https://github.com/JuliaSymbolics/Symbolics.jl/issues/1535).

> [@ChrisRackauckas](#):
>
> it might be hitting the wrong dispatch.

This is unreasonable. Since it invents the specialized vector, I guess there are some motivations (for speed etc.) Then the specialized method should have been more reliable. Otherwise we don’t need to invent new types of containers.

As a comparison

```julia
julia> using JuMP

julia> m = JuMP.Model(); @variable(m, x[1:3]) # the default container
3-element Vector{VariableRef}:
 x[1]
 x[2]
 x[3]

julia> exit()
PS K:\uc24> julia

julia> using Symbolics

julia> @variables x[1:3]
1-element Vector{Symbolics.Arr{Num, 1}}:
 x[1:3]

julia> typeof(x) # the special container
Symbolics.Arr{Num, 1}

julia> collect(x)
3-element Vector{Num}:
 x[1]
 x[2]
 x[3]

```

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [May 13, 2025, 12:03am UTC](https://discourse.julialang.org/t/hessian-is-not-jacobian-compose-gradient-in-symbolics-jl/128921/5 "2025-05-13T00:03:12Z")

</div>

> [@WalterMadelim](#):
>
> This is unreasonable. Since it invents the specialized vector, I guess there are some motivations (for speed etc.) Then the specialized method should have been more reliable. Otherwise we don’t need to invent new types of containers.

Yes, it has an O(1) representation. But there isn’t necessarily an O(1) representation of the Hessian, so the safe thing to do is just unroll here. If we end up with block vector things and someone wants to specialize it to improve performance then we can get to that, but for now [Fix and better test derivatives w.r.t. array variables by ChrisRackauckas · Pull Request #1536 · JuliaSymbolics/Symbolics.jl · GitHub](https://github.com/JuliaSymbolics/Symbolics.jl/pull/1536) should fix this case. For the record, I put some tests on Gradient and Jacobian as well and those seemed to work, so it’s just Hessian that got missed for some reason, so this now tests it more comprehensively.