# Odd behaviour regarding Adjoint Sensitivity Analysis
**URL:** https://discourse.julialang.org/t/odd-behaviour-regarding-adjoint-sensitivity-analysis/121141
**Category:** Numerics
**Tags:** question, sciml, differentialequation
**Created:** [October 10, 2024, 12:28pm UTC](https://discourse.julialang.org/t/odd-behaviour-regarding-adjoint-sensitivity-analysis/121141 "2024-10-10T12:28:58Z")
**Posts on this page:** 9
**Page:** 1
### Author:  [@nil](https://discourse.julialang.org/u/nil)
#### Post date: [October 10, 2024, 12:28pm UTC](https://discourse.julialang.org/t/odd-behaviour-regarding-adjoint-sensitivity-analysis/121141/1 "2024-10-10T12:28:58Z")
Hello dear community, I have a question regarding the following code I have written. I am trying to calculate the adjoint of the simple (linear) systems forward simulation with respect to the input u.
```julia
using LinearAlgebra
using DifferentialEquations
using SciMLSensitivity
using Interpolations
using Zygote
using ForwardDiff
using Plots
T = 4.0
N = 40
dt = T / N
np = N + 1
A = [0 1 0 0; 0 0 1 0; 0 0 0 1;-2 0 1 0]
B = [0,0,0,1]
uk = rand(np)
function f(x,u)
return A*x + B*u
end
function S_x(u)
function forward(x, p, t)
return f(x,p(t))
end
u0 = [1.0, 0.0, -1.0, 0.0]
tspan = (0.0, T)
prob = ODEProblem(forward, u0, tspan, linear_interpolation(0.0:dt:T,u))
sol = solve(prob,Midpoint(),saveat = T/N)
return sol.u
end
x_test = S_x(uk)
Zygote.jacobian(S_x, uk)
ForwardDiff.jacobian(S_x, uk)
```
Now I can evaluate the solution easily, however calculating the derivative w.r.t. u does not work with neither Zygote.jl nor Enzyme.jl. Zygote yields
> ERROR: MethodError: no method matching zero(::Type{Vector{Float64}})
ForwardDiff.jl works but yields weird results, a huge matrix filled with entries of type
Vector{ForwardDiff.Dual{ForwardDiff.Tag{typeof(S\_x), Float64}, Float64, 11}}
Maybe there is a better way for these kinds of problems, since my code is very naive right now.
Thanks in advance
---
### Author:  [@Vitor\_Patricio\_Canta](https://discourse.julialang.org/u/Vitor_Patricio_Canta)
#### Post date: [October 10, 2024, 8:22pm UTC](https://discourse.julialang.org/t/odd-behaviour-regarding-adjoint-sensitivity-analysis/121141/2 "2024-10-10T20:22:57Z")
Hi,
It is a bit hard to understand what you are trying to do there. The type to the result from ForwardDiff is correct. This weird type is what calculates the derivatives. Guess zygote is not compatible with your linear interpolation function.
I dont think you need to compute the adjoint by hand unless you want to deeply understand it. the solution should be differentiable (ScimlSensitivity documentation has more information on that). There are some lectures about adjoints in the SciML ecosystem around there, I can’t search it now, but you should take a look.
---
### Author:  [@nil](https://discourse.julialang.org/u/nil)
#### Post date: [October 10, 2024, 10:33pm UTC](https://discourse.julialang.org/t/odd-behaviour-regarding-adjoint-sensitivity-analysis/121141/3 "2024-10-10T22:33:06Z")
I will explain what I am trying to achieve. The function S\_x computes the solution of the differential equation
\dot x = Ax + Bu, ~~ x(0) = x\_0
which is a linear differential equation for now but could also be nonlinear.
S\_x takes as arguent u, which is the control function (a vector). I want to compute the derivative with respect to u. Now this can be easily writen down and calculated, even analytically utilizing the adjoint system, however I want to use automatic differentiation.
The result still seems odd to me, I need the matrix with the floating point values not the Dual type. This example is actually quite similar to the one in the SciMLSensitivitiy documentation, just that I use interpolation of the input paramters. I think Interpolation is not the problem, I can easily differentiate thorugh interpolations and it gives matrices with numeric floating point values.
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [October 10, 2024, 10:54pm UTC](https://discourse.julialang.org/t/odd-behaviour-regarding-adjoint-sensitivity-analysis/121141/4 "2024-10-10T22:54:35Z")
> [@nil](#):
>
> `prob = ODEProblem(forward, u0, tspan, linear_interpolation(0.0:dt:T,u))`
I’m not sure if this is going to work. `p` is expected to be a SciMLStructure, generally an Array. A `linear_interpolation` from Interpolations.jl is not an Array.
Instead you’d need to do `prob = ODEProblem(forward, u0, tspan, u)` and
```julia
function forward(x, p, t)
return f(x,linear_interpolation(0.0:dt:T,p)(t))
end
```
should work. Check that first. Caching that computation can come second.
---
### Author:  [@nil](https://discourse.julialang.org/u/nil)
#### Post date: [October 11, 2024, 11:10am UTC](https://discourse.julialang.org/t/odd-behaviour-regarding-adjoint-sensitivity-analysis/121141/5 "2024-10-11T11:10:58Z")
Hello Chris, thank you for your answer.
The code changes did not result in a solution, the error messages still persist for Zygote and ForwardDiff calculates weird results. Here is the part of the changed code
```julia
uk = rand(np)#zeros(np)
function f(x,u)
return A*x + B*u
end
function S_x(u)
function forward(x, p, t)
return f(x,linear_interpolation(0.0:dt:T,p)(t))
end
u0 = [1.0, 0.0, -1.0, 0.0]
tspan = (0.0, T)
prob = ODEProblem(forward, u0, tspan, u)
sol = solve(prob,Midpoint(),saveat = T/N)
return sol.u
end
S_x(uk) #works fine
ForwardDiff.jacobian(S_x,uk) #weird result
Zygote.jacobian(S_x,uk) #error
```
The output for ForwardDiff is still the arrangement of type
> Matrix{Vector{ForwardDiff.Dual{ForwardDiff.Tag{typeof(S\_x), Float64}, Float64, 11}}}:
while Zygote yields
> ERROR: MethodError: no method matching zero(::Type{Vector{Float64}})
Thanks for your help
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [October 11, 2024, 11:49am UTC](https://discourse.julialang.org/t/odd-behaviour-regarding-adjoint-sensitivity-analysis/121141/6 "2024-10-11T11:49:49Z")
What is the full stacktrace?
---
### Author:  [@nil](https://discourse.julialang.org/u/nil)
#### Post date: [October 11, 2024, 11:59am UTC](https://discourse.julialang.org/t/odd-behaviour-regarding-adjoint-sensitivity-analysis/121141/7 "2024-10-11T11:59:25Z")
> ERROR: MethodError: no method matching zero(::Type{Vector{Float64}})
>
> Closest candidates are:
> zero(::Type{Union{}}, Any…)
> @ Base number.jl:310
> zero(::Type{LibGit2.GitHash})
> @ LibGit2 ~/.julia/juliaup/julia-1.10.5+0.x64.linux.gnu/share/julia/stdlib/v1.10/LibGit2/src/oid.jl:221
> zero(::Type{Dates.Date})
> @ Dates ~/.julia/juliaup/julia-1.10.5+0.x64.linux.gnu/share/julia/stdlib/v1.10/Dates/src/types.jl:439
> …
>
> Stacktrace:
> [1] zeros(::Type{Vector{Float64}}, dims::Tuple{Int64, Int64})
> @ Base ./array.jl:637
> [2] Matrix{Vector{Float64}}(s::UniformScaling{Bool}, dims::Tuple{Int64, Int64})
> @ LinearAlgebra ~/.julia/juliaup/julia-1.10.5+0.x64.linux.gnu/share/julia/stdlib/v1.10/LinearAlgebra/src/uniformscaling.jl:407
> [3] Matrix{Vector{Float64}}(s::UniformScaling{Bool}, m::Int64, n::Int64)
> @ LinearAlgebra ~/.julia/juliaup/julia-1.10.5+0.x64.linux.gnu/share/julia/stdlib/v1.10/LinearAlgebra/src/uniformscaling.jl:414
> [4] \_eyelike(y::Vector{Vector{Float64}})
> @ Zygote ~/.julia/packages/Zygote/Tt5Gx/src/lib/grad.jl:164
> [5] withjacobian(f::Function, args::Vector{Float64})
> @ Zygote ~/.julia/packages/Zygote/Tt5Gx/src/lib/grad.jl:148
> [6] jacobian(f::Function, args::Vector{Float64})
> @ Zygote ~/.julia/packages/Zygote/Tt5Gx/src/lib/grad.jl:128
> [7] top-level scope
> @ REPL[88]:1
Julia was installed with juliaup and the versioninfo() yields
> Julia Version 1.10.5
> Commit 6f3fdf7b362 (2024-08-27 14:19 UTC)
> Build Info:
> Official [https://julialang.org/](https://julialang.org/) release
> Platform Info:
> OS: Linux (x86\_64-linux-gnu)
> CPU: 8 × Intel(R) Core™ i5-8250U CPU @ 1.60GHz
> WORD\_SIZE: 64
> LIBM: libopenlibm
> LLVM: libLLVM-15.0.7 (ORCJIT, skylake)
> Threads: 1 default, 0 interactive, 1 GC (on 8 virtual cores)
> Environment:
> JULIA\_LOAD\_PATH = /usr/share/gmsh/api/julia:
Edit: Updating to 1.11.0 right now using juliaup
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [October 11, 2024, 12:21pm UTC](https://discourse.julialang.org/t/odd-behaviour-regarding-adjoint-sensitivity-analysis/121141/8 "2024-10-11T12:21:44Z")
You were just returning a vector of vectors. AD libraries do not know how to represent Jacobians of such outputs. Instead, just return a normal array and it’s fine:
```julia
using LinearAlgebra
using OrdinaryDiffEq
using SciMLSensitivity
using Interpolations
using Zygote
using ForwardDiff
using Plots
T = 4.0
N = 40
dt = T / N
np = N + 1
A = [0 1 0 0; 0 0 1 0; 0 0 0 1;-2 0 1 0]
B = [0,0,0,1]
uk = rand(np)#zeros(np)
function f(x,u)
return A*x + B*u
end
function S_x(u)
function forward(x, p, t)
return f(x,linear_interpolation(0.0:dt:T,p)(t))
end
u0 = [1.0, 0.0, -1.0, 0.0]
tspan = (0.0, T)
prob = ODEProblem(forward, u0, tspan, u)
sol = solve(prob,Midpoint(),saveat = T/N)
return Array(sol)
end
S_x(uk)
ForwardDiff.jacobian(S_x,uk)
Zygote.jacobian(S_x,uk)
```
---
### Author:  [@nil](https://discourse.julialang.org/u/nil)
#### Post date: [October 11, 2024, 9:49pm UTC](https://discourse.julialang.org/t/odd-behaviour-regarding-adjoint-sensitivity-analysis/121141/9 "2024-10-11T21:49:13Z")
Sorry for the late reply. This indeed does the trick! Thank you for your solution!