# Forward differentiation and differential equations

**URL:** <https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002>\
**Category:** Modelling & Simulations\
**Tags:** diffeq, differentiation, optim, optimization\
**Created:** [February 20, 2019, 7:08am UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002 "2019-02-20T07:08:51Z")\
**Posts on this page:** 14\
**Page:** 1

<div class="post-metadata">

**Author:** ![moesphere](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/moesphere/32/3348_2.png) [@moesphere](https://discourse.julialang.org/u/moesphere)\
**Post date:** [February 20, 2019, 7:08am UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002/1 "2019-02-20T07:08:51Z")

</div>

I am trying to use forward AD in combination with solving an ODE using DifferentialEquations.jl. I want to use an own loss function. However, with the is-in-place version there is a type assertion error. I have come up with following minimal example:

```julia
using DifferentialEquations, Optim, LinearAlgebra

function f(dxdt,x,p,t)
    dxdt[1] = -p[1]*x[1]
    dxdt[2] = p[1]*x[1]
end

tspan = (0.,8.)
tdata = range(0.,stop=8.,length=10)
p0 = [1.]
x0 = [10.,0.]

prob = ODEProblem(f,x0,tspan,p0)

sol = DifferentialEquations.solve(prob; saveat=tdata)

xdata = sol[2,:] + randn(length(sol[2,:])) # in-silico data

function loss_function(x)
    prob_new = remake(prob; p=x)
    sol_new = DifferentialEquations.solve(prob_new; saveat=tdata)
    norm(sol_new[2,:] - xdata)^2
end
optimize(loss_function, [2.0], LBFGS(),autodiff=:forward)

```

Error message is `TypeError: in typeassert, expected Float64, got ForwardDiff.Dual{Nothing,Float64,1}`.

Is there a solution to this? It works fine when the `autodiff` keyword is removed, i.e. using finite differences.

---

<div class="post-metadata">

**Author:** ![pkofod](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pkofod/32/2179_2.png) [@pkofod](https://discourse.julialang.org/u/pkofod)\
**Post date:** [February 20, 2019, 10:41am UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002/2 "2019-02-20T10:41:37Z")

</div>

Why would you rob people of the whole stack trace? 🙂

```julia
julia> optimize(loss_function, [2.0], LBFGS(),autodiff=:forward)
ERROR: MethodError: no method matching Float64(::ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1})
Closest candidates are:
  Float64(::Real, ::RoundingMode) where T<:AbstractFloat at rounding.jl:185
  Float64(::T<:Number) where T<:Number at boot.jl:725
  Float64(::Int8) at float.jl:60
  ...
Stacktrace:
 [1] convert(::Type{Float64}, ::ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1}) at ./number.jl:7
 [2] setindex!(::Array{Float64,1}, ::ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1}, ::Int64) at ./array.jl:769
 [3] f(::Array{Float64,1}, ::Array{Float64,1}, ::Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1}, ::Float64) at ./REPL[15]:2
 [4] (::ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing})(::Array{Float64,1}, ::Array{Float64,1}, ::Vararg{Any,N} where N) at /home/pkofod/.julia/packages/DiffEqBase/tENmC/src/diffeqfunction.jl:107
 [5] initialize!(::OrdinaryDiffEq.ODEIntegrator{CompositeAlgorithm{Tuple{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType}},AutoSwitch{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType},Rational{Int64},Float64}},true,Array{Float64,1},Float64,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1},Float64,Float64,Float64,Array{Array{Float64,1},1},OrdinaryDiffEq.ODECompositeSolution{Float64,2,Array{Array{Float64,1},1},Nothing,Nothing,Array{Float64,1},Array{Array{Array{Float64,1},1},1},ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1},ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Nothing,DiffEqBase.StandardODEProblem},CompositeAlgorithm{Tuple{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType}},AutoSwitch{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType},Rational{Int64},Float64}},OrdinaryDiffEq.CompositeInterpolationData{ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Array{Float64,1},1},Array{Float64,1},Array{Array{Array{Float64,1},1},1},OrdinaryDiffEq.CompositeCache{Tuple{OrdinaryDiffEq.Tsit5Cache{Array{Float64,1},Array{Float64,1},Array{Float64,1},OrdinaryDiffEq.Tsit5ConstantCache{Float64,Float64}},OrdinaryDiffEq.Rosenbrock23Cache{Array{Float64,1},Array{Float64,1},Array{Float64,2},Array{Float64,2},OrdinaryDiffEq.Rosenbrock23ConstantCache{Float64,typeof(identity),typeof(identity)},DiffEqDiffTools.TimeGradientWrapper{ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Float64,1},Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1}},DiffEqDiffTools.UJacobianWrapper{ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1}},LinSolveFactorize{typeof(lu!)},DiffEqDiffTools.JacobianCache{Array{Float64,1},Array{Float64,1},Array{Float64,1},Val{:central},Float64,Val{true}},DiffEqDiffTools.GradientCache{Nothing,Array{Float64,1},Array{Float64,1},Val{:central},Float64,Val{true}}}},AutoSwitch{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType},Rational{Int64},Float64}}}},ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},OrdinaryDiffEq.CompositeCache{Tuple{OrdinaryDiffEq.Tsit5Cache{Array{Float64,1},Array{Float64,1},Array{Float64,1},OrdinaryDiffEq.Tsit5ConstantCache{Float64,Float64}},OrdinaryDiffEq.Rosenbrock23Cache{Array{Float64,1},Array{Float64,1},Array{Float64,2},Array{Float64,2},OrdinaryDiffEq.Rosenbrock23ConstantCache{Float64,typeof(identity),typeof(identity)},DiffEqDiffTools.TimeGradientWrapper{ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Float64,1},Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1}},DiffEqDiffTools.UJacobianWrapper{ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1}},LinSolveFactorize{typeof(lu!)},DiffEqDiffTools.JacobianCache{Array{Float64,1},Array{Float64,1},Array{Float64,1},Val{:central},Float64,Val{true}},DiffEqDiffTools.GradientCache{Nothing,Array{Float64,1},Array{Float64,1},Val{:central},Float64,Val{true}}}},AutoSwitch{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType},Rational{Int64},Float64}},OrdinaryDiffEq.DEOptions{Float64,Float64,Float64,Float64,typeof(DiffEqBase.ODE_DEFAULT_NORM),typeof(opnorm),CallbackSet{Tuple{},Tuple{}},typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN),typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE),typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK),DataStructures.BinaryHeap{Float64,DataStructures.LessThan},DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Nothing,Nothing,Int64,Array{Float64,1},StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}},Array{Float64,1}},Array{Float64,1},Float64}, ::OrdinaryDiffEq.Tsit5Cache{Array{Float64,1},Array{Float64,1},Array{Float64,1},OrdinaryDiffEq.Tsit5ConstantCache{Float64,Float64}}) at /home/pkofod/.julia/packages/OrdinaryDiffEq/HJOah/src/perform_step/low_order_rk_perform_step.jl:580
 [6] initialize!(::OrdinaryDiffEq.ODEIntegrator{CompositeAlgorithm{Tuple{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType}},AutoSwitch{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType},Rational{Int64},Float64}},true,Array{Float64,1},Float64,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1},Float64,Float64,Float64,Array{Array{Float64,1},1},OrdinaryDiffEq.ODECompositeSolution{Float64,2,Array{Array{Float64,1},1},Nothing,Nothing,Array{Float64,1},Array{Array{Array{Float64,1},1},1},ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1},ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Nothing,DiffEqBase.StandardODEProblem},CompositeAlgorithm{Tuple{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType}},AutoSwitch{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType},Rational{Int64},Float64}},OrdinaryDiffEq.CompositeInterpolationData{ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Array{Float64,1},1},Array{Float64,1},Array{Array{Array{Float64,1},1},1},OrdinaryDiffEq.CompositeCache{Tuple{OrdinaryDiffEq.Tsit5Cache{Array{Float64,1},Array{Float64,1},Array{Float64,1},OrdinaryDiffEq.Tsit5ConstantCache{Float64,Float64}},OrdinaryDiffEq.Rosenbrock23Cache{Array{Float64,1},Array{Float64,1},Array{Float64,2},Array{Float64,2},OrdinaryDiffEq.Rosenbrock23ConstantCache{Float64,typeof(identity),typeof(identity)},DiffEqDiffTools.TimeGradientWrapper{ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Float64,1},Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1}},DiffEqDiffTools.UJacobianWrapper{ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1}},LinSolveFactorize{typeof(lu!)},DiffEqDiffTools.JacobianCache{Array{Float64,1},Array{Float64,1},Array{Float64,1},Val{:central},Float64,Val{true}},DiffEqDiffTools.GradientCache{Nothing,Array{Float64,1},Array{Float64,1},Val{:central},Float64,Val{true}}}},AutoSwitch{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType},Rational{Int64},Float64}}}},ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},OrdinaryDiffEq.CompositeCache{Tuple{OrdinaryDiffEq.Tsit5Cache{Array{Float64,1},Array{Float64,1},Array{Float64,1},OrdinaryDiffEq.Tsit5ConstantCache{Float64,Float64}},OrdinaryDiffEq.Rosenbrock23Cache{Array{Float64,1},Array{Float64,1},Array{Float64,2},Array{Float64,2},OrdinaryDiffEq.Rosenbrock23ConstantCache{Float64,typeof(identity),typeof(identity)},DiffEqDiffTools.TimeGradientWrapper{ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Float64,1},Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1}},DiffEqDiffTools.UJacobianWrapper{ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1}},LinSolveFactorize{typeof(lu!)},DiffEqDiffTools.JacobianCache{Array{Float64,1},Array{Float64,1},Array{Float64,1},Val{:central},Float64,Val{true}},DiffEqDiffTools.GradientCache{Nothing,Array{Float64,1},Array{Float64,1},Val{:central},Float64,Val{true}}}},AutoSwitch{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType},Rational{Int64},Float64}},OrdinaryDiffEq.DEOptions{Float64,Float64,Float64,Float64,typeof(DiffEqBase.ODE_DEFAULT_NORM),typeof(opnorm),CallbackSet{Tuple{},Tuple{}},typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN),typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE),typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK),DataStructures.BinaryHeap{Float64,DataStructures.LessThan},DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Nothing,Nothing,Int64,Array{Float64,1},StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}},Array{Float64,1}},Array{Float64,1},Float64}, ::OrdinaryDiffEq.CompositeCache{Tuple{OrdinaryDiffEq.Tsit5Cache{Array{Float64,1},Array{Float64,1},Array{Float64,1},OrdinaryDiffEq.Tsit5ConstantCache{Float64,Float64}},OrdinaryDiffEq.Rosenbrock23Cache{Array{Float64,1},Array{Float64,1},Array{Float64,2},Array{Float64,2},OrdinaryDiffEq.Rosenbrock23ConstantCache{Float64,typeof(identity),typeof(identity)},DiffEqDiffTools.TimeGradientWrapper{ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Float64,1},Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1}},DiffEqDiffTools.UJacobianWrapper{ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1}},LinSolveFactorize{typeof(lu!)},DiffEqDiffTools.JacobianCache{Array{Float64,1},Array{Float64,1},Array{Float64,1},Val{:central},Float64,Val{true}},DiffEqDiffTools.GradientCache{Nothing,Array{Float64,1},Array{Float64,1},Val{:central},Float64,Val{true}}}},AutoSwitch{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType},Rational{Int64},Float64}}) at /home/pkofod/.julia/packages/OrdinaryDiffEq/HJOah/src/perform_step/composite_perform_step.jl:38
 [7] #__init#303(::StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}, ::Array{Float64,1}, ::Array{Float64,1}, ::Nothing, ::Bool, ::Nothing, ::Bool, ::Bool, ::Bool, ::Nothing, ::Bool, ::Bool, ::Float64, ::Bool, ::Rational{Int64}, ::Nothing, ::Nothing, ::Int64, ::Rational{Int64}, ::Int64, ::Int64, ::Rational{Int64}, ::Bool, ::Int64, ::Nothing, ::Nothing, ::Int64, ::Float64, ::Float64, ::typeof(DiffEqBase.ODE_DEFAULT_NORM), ::typeof(opnorm), ::typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN), ::typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK), ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::Int64, ::String, ::typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE), ::Nothing, ::Bool, ::Bool, ::Bool, ::Base.Iterators.Pairs{Symbol,Bool,Tuple{Symbol},NamedTuple{(:default_set,),Tuple{Bool}}}, ::typeof(DiffEqBase.__init), ::ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1},ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Nothing,DiffEqBase.StandardODEProblem}, ::CompositeAlgorithm{Tuple{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType}},AutoSwitch{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType},Rational{Int64},Float64}}, ::Array{Any,1}, ::Array{Any,1}, ::Array{Any,1}, ::Type{Val{true}}) at /home/pkofod/.julia/packages/OrdinaryDiffEq/HJOah/src/solve.jl:273
 [8] (::getfield(DiffEqBase, Symbol("#kw##__init")))(::NamedTuple{(:default_set, :saveat),Tuple{Bool,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}}, ::typeof(DiffEqBase.__init), ::ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1},ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Nothing,DiffEqBase.StandardODEProblem}, ::CompositeAlgorithm{Tuple{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType}},AutoSwitch{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType},Rational{Int64},Float64}}, ::Array{Any,1}, ::Array{Any,1}, ::Array{Any,1}, ::Type{Val{true}}) at ./none:0
 [9] #__solve#302(::Base.Iterators.Pairs{Symbol,Any,Tuple{Symbol,Symbol},NamedTuple{(:default_set, :saveat),Tuple{Bool,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}}}, ::Function, ::ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1},ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Nothing,DiffEqBase.StandardODEProblem}, ::CompositeAlgorithm{Tuple{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType}},AutoSwitch{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType},Rational{Int64},Float64}}, ::Array{Any,1}, ::Array{Any,1}, ::Array{Any,1}, ::Type{Val{true}}) at /home/pkofod/.julia/packages/OrdinaryDiffEq/HJOah/src/solve.jl:6
 [10] (::getfield(DiffEqBase, Symbol("#kw##__solve")))(::NamedTuple{(:default_set, :saveat),Tuple{Bool,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}}, ::typeof(DiffEqBase.__solve), ::ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1},ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Nothing,DiffEqBase.StandardODEProblem}, ::CompositeAlgorithm{Tuple{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType}},AutoSwitch{Tsit5,Rosenbrock23{0,false,LinSolveFactorize{typeof(lu!)},DataType},Rational{Int64},Float64}}) at ./none:0
 [11] #__solve#2(::Bool, ::Base.Iterators.Pairs{Symbol,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}},Tuple{Symbol},NamedTuple{(:saveat,),Tuple{StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}}}, ::Function, ::ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1},ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Nothing,DiffEqBase.StandardODEProblem}, ::Nothing) at /home/pkofod/.julia/packages/DifferentialEquations/UZIcy/src/default_solve.jl:15
 [12] (::getfield(DiffEqBase, Symbol("#kw##solve")))(::NamedTuple{(:saveat,),Tuple{StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}}, ::typeof(solve), ::ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1},ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Nothing,DiffEqBase.StandardODEProblem}) at ./none:0
 [13] loss_function(::Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1}) at ./REPL[23]:3
 [14] vector_mode_gradient!(::DiffResults.MutableDiffResult{1,Float64,Tuple{Array{Float64,1}}}, ::typeof(loss_function), ::Array{Float64,1}, ::ForwardDiff.GradientConfig{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1}}) at /home/pkofod/.julia/packages/ForwardDiff/okZnq/src/apiutils.jl:35
 [15] gradient! at /home/pkofod/.julia/packages/ForwardDiff/okZnq/src/gradient.jl:35 [inlined]
 [16] gradient! at /home/pkofod/.julia/packages/ForwardDiff/okZnq/src/gradient.jl:33 [inlined]
 [17] (::getfield(NLSolversBase, Symbol("##14#18")){Float64,typeof(loss_function),ForwardDiff.GradientConfig{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1,Array{ForwardDiff.Dual{ForwardDiff.Tag{typeof(loss_function),Float64},Float64,1},1}}})(::Array{Float64,1}, ::Array{Float64,1}) at /home/pkofod/.julia/packages/NLSolversBase/LRo6O/src/objective_types/oncedifferentiable.jl:65
 [18] value_gradient!!(::OnceDifferentiable{Float64,Array{Float64,1},Array{Float64,1}}, ::Array{Float64,1}) at /home/pkofod/.julia/packages/NLSolversBase/LRo6O/src/interface.jl:82
 [19] initial_state(::LBFGS{Nothing,LineSearches.InitialStatic{Float64},LineSearches.HagerZhang{Float64,Base.RefValue{Bool}},getfield(Optim, Symbol("##22#24"))}, ::Optim.Options{Float64,Nothing}, ::OnceDifferentiable{Float64,Array{Float64,1},Array{Float64,1}}, ::Array{Float64,1}) at /home/pkofod/.julia/dev/Optim/src/multivariate/solvers/first_order/l_bfgs.jl:158
 [20] optimize(::OnceDifferentiable{Float64,Array{Float64,1},Array{Float64,1}}, ::Array{Float64,1}, ::LBFGS{Nothing,LineSearches.InitialStatic{Float64},LineSearches.HagerZhang{Float64,Base.RefValue{Bool}},getfield(Optim, Symbol("##22#24"))}, ::Optim.Options{Float64,Nothing}) at /home/pkofod/.julia/dev/Optim/src/multivariate/optimize/optimize.jl:33
 [21] #optimize#87 at /home/pkofod/.julia/dev/Optim/src/multivariate/optimize/interface.jl:114 [inlined]
 [22] (::getfield(Optim, Symbol("#kw##optimize")))(::NamedTuple{(:autodiff,),Tuple{Symbol}}, ::typeof(optimize), ::Function, ::Array{Float64,1}, ::LBFGS{Nothing,LineSearches.InitialStatic{Float64},LineSearches.HagerZhang{Float64,Base.RefValue{Bool}},getfield(Optim, Symbol("##22#24"))}, ::Optim.Options{Float64,Nothing}) at ./none:0 (repeats 2 times)
 [23] top-level scope at none:0

```

So without knowing the internals of ODEFunction it’s probably because there is some `convert` (attempt to take the input and make sure it has a specific type for the rest of the code to be type stable) that is based of off some parametric type. This gets initialized in your constructor, so it’s forced to be Float64 at the time of ForwardDiff’ing - causing the error you see. @ChrisRackauckas do you recognize this problem with `ODEProblem`?

edit: oooh… the stack trace (though it has a vertical bar) is much easier to read without line wrapping

---

<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:** [February 20, 2019, 11:17am UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002/3 "2019-02-20T11:17:34Z")

</div>

The automatic differentiation page discusses exactly this issue:

[http://docs.juliadiffeq.org/latest/analysis/sensitivity.html#Examples-using-ForwardDiff.jl-1](http://docs.juliadiffeq.org/latest/analysis/sensitivity.html#Examples-using-ForwardDiff.jl-1)

If you want to propagate derivatives in state, you need to make sure your state variables are dual numbers as well, so you’re missing the type conversion of `u0` to be the same as `p`. Your `remake` should be:

```julia
_prob = remake(prob;u0=convert.(eltype(p),prob.u0),p=p)

```

---

<div class="post-metadata">

**Author:** ![moesphere](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/moesphere/32/3348_2.png) [@moesphere](https://discourse.julialang.org/u/moesphere)\
**Post date:** [February 20, 2019, 11:58am UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002/4 "2019-02-20T11:58:17Z")

</div>

Thanks for the answer!

I should have seen this in the docs.

---

<div class="post-metadata">

**Author:** ![Matthew\_Overlin](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/matthew_overlin/32/10533_2.png) [@Matthew\_Overlin](https://discourse.julialang.org/u/Matthew_Overlin)\
**Post date:** [September 28, 2020, 1:45pm UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002/6 "2020-09-28T13:45:13Z")

</div>

I am getting this issue with delay differential equations problem. I define: dde\_problem = DDEProblem(…); Then, I remake the problem as: tmp\_dde\_problem = remake(dde\_problem; u0=convert.(eltype(my\_parameters), dde\_problem.u0),p=my\_parameters); Then, I try to solve tmp\_dde\_problem and I get an error of: ERROR: LoadError: TypeError: in typeassert, expected Float64, got a value of type ForwardDiff.Dual{Nothing,Float64,9}

I would expect similar behavior when ForwardDiff is used with an ODEProblem and a DDEProblem. I have had success using ForwardDiff with ODEProblem. Could anyone provide insights as to why it might not be working for DDEProblem?  
Thanks.

---

<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:** [September 28, 2020, 1:47pm UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002/7 "2020-09-28T13:47:59Z")

</div>

I’d need more details since this is covered by tests

[https://github.com/SciML/DelayDiffEq.jl/blob/master/test/interface/ad.jl](https://github.com/SciML/DelayDiffEq.jl/blob/master/test/interface/ad.jl)

---

<div class="post-metadata">

**Author:** ![pierre-haessig](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pierre-haessig/32/217129_2.png) [@pierre-haessig](https://discourse.julialang.org/u/pierre-haessig)\
**Post date:** [August 31, 2022, 2:48pm UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002/8 "2022-08-31T14:48:26Z")

</div>

Just searching for infos on the topic of ForwardDiff and type issues, I came across this thread and noticed that the documentation link to [juliadiffeq.org](http://juliadiffeq.org) is broken. Here are updates:

- Archived version of late 2019 [Local Sensitivity Analysis (Automatic Differentiation) · DifferentialEquations.jl](https://web.archive.org/web/20191027231222/http://docs.juliadiffeq.org/latest/analysis/sensitivity.html)
- The closest document on the new sciml.ai website: [Local Sensitivity Analysis (Automatic Differentiation) · DifferentialEquations.jl](https://diffeq.sciml.ai/stable/analysis/sensitivity/)

However, in the new doc, I see no more any usage of `convert` inside the call to `remake`.

---

<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:** [August 31, 2022, 2:50pm UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002/9 "2022-08-31T14:50:54Z")

</div>

> [@pierre-haessig](#):
>
> However, in the new doc, I see no more any usage of `convert` inside the call to `remake`.

It’s handled automatically these days.

---

<div class="post-metadata">

**Author:** ![edljk](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/edljk/32/2429_2.png) [@edljk](https://discourse.julialang.org/u/edljk)\
**Post date:** [November 5, 2023, 8:27pm UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002/10 "2023-11-05T20:27:20Z")

</div>

A related question: as mentioned above and in the [ode faq](https://docs.sciml.ai/DiffEqDocs/stable/basics/faq/), automatic differentiation of the final state with respect to parameters using ForwardDiff works perfectly.

I am interested in computing all the jacobians associated to many final states with respect to parameters but with different initial states (with the same parameters).

I wonder if there is a way to keep the benefit of the multi-threading provided by an `EnsembleProblem` in that context?

Thanks.

---

<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:** [November 5, 2023, 8:34pm UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002/11 "2023-11-05T20:34:31Z")

</div>

> [@edljk](#):
>
> I wonder if there is a way to keep the benefit of the multi-threading provided by an `EnsembleProblem` in that context?

Kind of. It’s a bit easier to just use PolyesterForwardDiff if you just want multithreading though, but yeah if you build the duals yourself you can make use of ensembles.

---

<div class="post-metadata">

**Author:** ![edljk](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/edljk/32/2429_2.png) [@edljk](https://discourse.julialang.org/u/edljk)\
**Post date:** [November 5, 2023, 9:09pm UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002/12 "2023-11-05T21:09:12Z")

</div>

Thanks for the fast answer. Since I am not enough familiar with ForwardDiff and duals I am gonna try PolyesterForwardDiff.jl which looks great!  
Btw, in my setting, PolyesterForwardDiff.jl seems to downgrade many packages like:

⌅ [aae01518] ↓ BandedMatrices v1.1.0 ⇒ v0.17.18  
⌃ [052768ef] ↓ CUDA v5.0.0 ⇒ v4.4.1  
⌃ [13f3f980] ↓ CairoMakie v0.10.12 ⇒ v0.9.4  
⌅ [34da2185] ↓ Compat v4.10.0 ⇒ v3.46.2  
⌃ [a93c6f00] ↓ DataFrames v1.6.1 ⇒ v1.3.6  
⌃ [0c46a032] ↓ DifferentialEquations v7.11.0 ⇒ v7.2.0  
⌃ [e9467ef8] ↓ GLMakie v0.8.12 ⇒ v0.7.4  
⌅ [7ed4a6bd] ↓ LinearSolve v2.15.0 ⇒ v1.34.1  
…

Is it expected?

---

<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:** [November 5, 2023, 9:14pm UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002/13 "2023-11-05T21:14:57Z")

</div>

> <https://github.com/JuliaDiff/PolyesterForwardDiff.jl/pull/12>

---

<div class="post-metadata">

**Author:** ![edljk](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/edljk/32/2429_2.png) [@edljk](https://discourse.julialang.org/u/edljk)\
**Post date:** [November 5, 2023, 11:18pm UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002/14 "2023-11-05T23:18:00Z")

</div>

Probably some difficulties from previous update:

- num\_threads() is not defined anymore
- Defining `num_threads() = Threads.nthreads()` and adapting the example of the [ode faq](https://docs.sciml.ai/DiffEqDocs/stable/basics/faq/):

```julia
using DifferentialEquations
function func(du, u, p, t)
    du[1] = p[1] * u[1] - p[2] * u[1] * u[2]
    du[2] = -3 * u[2] + u[1] * u[2]
end
function f(p)
    prob = ODEProblem(func, eltype(p).([1.0, 1.0]), (0.0, 10.0), p)
    # Lower tolerances to show the methods converge to the same value
    solve(prob, Tsit5(), save_everystep = false, abstol = 1e-12, reltol = 1e-12)[end]
end
J = zeros(2, 2)
using PolyesterForwardDiff
PolyesterForwardDiff.threaded_jacobian!(f, J, [1.5, 1.0], ForwardDiff.Chunk(8))

```

gives the error

```julia
ERROR: Cannot determine ordering of Dual tags ForwardDiff.Tag{DiffEqBase.OrdinaryDiffEqTag, ForwardDiff.Dual{Nothing, Float64, 8}} and Nothing
Stacktrace:
  [1] ≺(a::Type, b::Type)
...

```

---

<div class="post-metadata">

**Author:** ![edljk](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/edljk/32/2429_2.png) [@edljk](https://discourse.julialang.org/u/edljk)\
**Post date:** [November 7, 2023, 9:01am UTC](https://discourse.julialang.org/t/forward-differentiation-and-differential-equations/21002/15 "2023-11-07T09:01:46Z")

</div>

I finally had a try with Dual evaluations using `EnsembleProblem` and following `ForwardDiff.jacobian` implementation.  
The naive implementation below seems to work on this simple example.  
Any comments are welcome.

```julia
using DifferentialEquations, ForwardDiff
using TimerOutputs

to_ode = TimerOutput()
#-------------------------------------------------------------------------------
function func(du, u, p, t)
    du[1] = p[1] * u[1] - p[2] * u[1] * u[2]
    du[2] = -3 * u[2] + u[1] * u[2]
end
function f(u0, λ, εode)
    prob = ODEProblem(func, u0, (0.0, 10.0), λ)
    solve(prob, Tsit5(), save_everystep = false, abstol = εode, 
          reltol = εode)[end]
end
#-------------------------------------------------------------------------------
function jacobiansODE(f, prob, Vu0, λ, solver,
    cfg::ForwardDiff.JacobianConfig{T} = ForwardDiff.JacobianConfig(f, λ);
    εode::Float64 = 1e-4, maxit::Int64 = 1_000) where T

    N = ForwardDiff.chunksize(cfg)
    Vjacobians, Vydual = Array[], Vector{ForwardDiff.Dual}[]

    λdual = cfg.duals
    ForwardDiff.seed!(λdual, λ, cfg.seeds)

    # use multi-threading to evaluate dual variables
    function prob_funcdual(prob, i, repeat)
        remake(prob, u0 = convert.(eltype(λdual), Vu0[i]), p = λdual)
    end
    ∇ensemble_prob = EnsembleProblem(prob, prob_func = prob_funcdual)
    sols = solve(∇ensemble_prob, solver, trajectories = length(Vu0),
                 maxiters = maxit, abstol = εode, reltol = εode,  
                 save_everystep = false)
    for sol ∈ sols 
        push!(Vydual, sol[end])
    end

    for ydual ∈ Vydual
        ydual isa AbstractArray || throw(ForwardDiff.JACOBIAN_ERROR)
        result = similar(ydual, ForwardDiff.valtype(eltype(ydual)), 
                        length(ydual), N)
        ForwardDiff.extract_jacobian!(T, result, ydual, N)
        ForwardDiff.extract_value!(T, result, ydual)
        push!(Vjacobians, [result[k, l].value for k = 1:size(result, 1), 
                                                  l = 1:size(result, 2)])
    end
    return Vjacobians
end
#-------------------------------------------------------------------------------
function testjacobians(nbexp::Int64 = 10)
    εode = 1e-8
    λ = [1.5, 1.0]
    allu0 = [(2 * rand(2) .- 1) / 10 .+ 1 for _ = 1:nbexp]
    prob = ODEProblem(func, allu0[1], (0.0, 10.0), λ)
    jacobians1 = Matrix[]
    for u0 ∈ allu0
        @timeit to_ode "ForwardDiff" push!(jacobians1, 
                                   ForwardDiff.jacobian(p -> f(u0, p, εode), λ))
    end
    @timeit to_ode "jacobiansODE" jacobians2 = jacobiansODE(f, prob, allu0, λ, 
                                                        Tsit5(), εode = εode)
    @show norm(jacobians1 - jacobians2, Inf)

    print(to_ode)
end

# call test function 
testjacobians(10)

```
