Error while using DynamicalSystems.jl to solve for a complex nonlinear system

Hi,

I was having quite a bit of trouble with DynamicalSystems.jl. In particular, I would like to implement an elastic double pendulum (so, 3 degrees of freedoms). The problem with this system is that it is nonlinear, second order and its quite troublesome and tedious to separate the different second order derivatives by hand.

I have solved this system numerically once already, separating the derivatives by using,

E₁₁, E₂₂, E₃₃ = Symbolics.solve_for([E₁ ~ 0, E₂ ~ 0, E₃ ~ 0], [d(d(θ)), d(d(ϕ)), d(d(r))])

(The E’s inside the parentheses are the equations of motion I have)

The problem is that this doesn’t seem to work with DynamicalSystems.jl (I am quite new to this package).

My code looks like this:

using Plots, ModelingToolkit, DifferentialEquations, LaTeXStrings, DynamicalSystems
@variables t θ(t) ϕ(t) r(t)
@parameters m₁, m₂, l₁, l₂, k, g
d = Differential(t)
M = m₁+m₂
ω₃² = k/m₂
μ = m₂/m₁
l₀ = g/(ω₃²)
l = l₂ + l₀ ##Defining some constants.
E₁ = …
E₂ = …
E₃ = … //These equations contain both the variables and their first-order
derivatives
E₁₁, E₂₂, E₃₃ = Symbolics.solve_for([E₁ ~ 0, E₂ ~ 0, E₃ ~ 0], [d(d(θ)), d(d(ϕ)), d(d(r))])
function truncated(du, u, p, t)
m₁ = p[1] ; m₂ = p[2]; l₁ = p[3]; l₂ = p[4]; k = p[5]; g = p[6]
u[1] = θ
u[2] = dθ
u[3] = ϕ
u[4] = dϕ
u[5] = r
u[6] = dr
du[1] = u[2]
du[3] = u[4]
du[5] = u[6]
du[2] = E₁₁
du[4] = E₂₂
du[6] = E₃₃
end

I then define the initial state, the parameters and the dynamical system as:

u_0 = [0, 1e-7, 0, 1e-7, 0.05, 0]
p_0 = [0.4, 0.8, 0.786, 0.7, 225.02, 9.8]
DS = ContinuousDynamicalSystem(truncated, u_0, p_0)

This gives me the following error:

MethodError: no method matching Float64(::Num)
Closest candidates are:
(::Type{T})(::Real, ::RoundingMode) where T<:AbstractFloat at rounding.jl:200
(::Type{T})(::T) where T<:Number at boot.jl:772
(::Type{T})(::AbstractChar) where T<:Union{AbstractChar, Number} at char.jl:50

Stacktrace:
[1] convert(#unused#::Type{Float64}, x::Num)
@ Base .\number.jl:7
[2] setindex!(A::Vector{Float64}, x::Num, i1::Int64)
@ Base .\array.jl:966
[3] truncated(du::Vector{Float64}, u::Vector{Float64}, p::Vector{Float64}, t::Float64)
@ Main .\In[29]:19
[4] ODEFunction
@ C:\Users\shiha.julia\packages\SciMLBase\nzvrv\src\scimlfunctions.jl:2096 [inlined]
[5] initialize!(integrator::OrdinaryDiffEq.ODEIntegrator{Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}, true, Vector{Float64}, Nothing, Float64, Vector{Float64}, Float64, Float64, Float64, Float64, Vector{Vector{Float64}}, ODESolution{Float64, 2, Vector{Vector{Float64}}, Nothing, Nothing, Vector{Float64}, Vector{Vector{Vector{Float64}}}, ODEProblem{Vector{Float64}, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, SciMLBase.AutoSpecialize, typeof(truncated), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}, OrdinaryDiffEq.InterpolationData{ODEFunction{true, SciMLBase.AutoSpecialize, typeof(truncated), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Vector{Vector{Float64}}, Vector{Float64}, Vector{Vector{Vector{Float64}}}, OrdinaryDiffEq.Tsit5Cache{Vector{Float64}, Vector{Float64}, Vector{Float64}, OrdinaryDiffEq.Tsit5ConstantCache{Float64, Float64}, typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}}, DiffEqBase.DEStats, Nothing}, ODEFunction{true, SciMLBase.AutoSpecialize, typeof(truncated), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, OrdinaryDiffEq.Tsit5Cache{Vector{Float64}, Vector{Float64}, Vector{Float64}, OrdinaryDiffEq.Tsit5ConstantCache{Float64, Float64}, typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}, OrdinaryDiffEq.DEOptions{Float64, Float64, Float64, Float64, PIController{Rational{Int64}}, typeof(DiffEqBase.ODE_DEFAULT_NORM), typeof(LinearAlgebra.opnorm), Bool, CallbackSet{Tuple{}, Tuple{}}, typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN), typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE), typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK), DataStructures.BinaryHeap{Float64, DataStructures.FasterForward}, DataStructures.BinaryHeap{Float64, DataStructures.FasterForward}, Nothing, Nothing, Float64, Tuple{}, Tuple{}, Tuple{}}, Vector{Float64}, Float64, Nothing, OrdinaryDiffEq.DefaultInit}, cache::OrdinaryDiffEq.Tsit5Cache{Vector{Float64}, Vector{Float64}, Vector{Float64}, OrdinaryDiffEq.Tsit5ConstantCache{Float64, Float64}, typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False})
@ OrdinaryDiffEq C:\Users\shiha.julia\packages\OrdinaryDiffEq\GxEzW\src\perform_step\low_order_rk_perform_step.jl:736
[6] __init(prob::ODEProblem{Vector{Float64}, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, SciMLBase.AutoSpecialize, typeof(truncated), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, alg::Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}, timeseries_init::Tuple{}, ts_init::Tuple{}, ks_init::Tuple{}, recompile::Type{Val{true}}; saveat::Tuple{}, tstops::Tuple{}, d_discontinuities::Tuple{}, save_idxs::Nothing, save_everystep::Bool, save_on::Bool, save_start::Bool, save_end::Bool, callback::Nothing, dense::Bool, calck::Bool, dt::Float64, dtmin::Nothing, dtmax::Float64, force_dtmin::Bool, adaptive::Bool, gamma::Rational{Int64}, abstol::Float64, reltol::Float64, qmin::Rational{Int64}, qmax::Int64, qsteady_min::Int64, qsteady_max::Int64, beta1::Nothing, beta2::Nothing, qoldinit::Rational{Int64}, controller::Nothing, fullnormalize::Bool, failfactor::Int64, maxiters::Float64, internalnorm::typeof(DiffEqBase.ODE_DEFAULT_NORM), internalopnorm::typeof(LinearAlgebra.opnorm), isoutofdomain::typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN), unstable_check::typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK), verbose::Bool, timeseries_errors::Bool, dense_errors::Bool, advance_to_tstop::Bool, stop_at_next_tstop::Bool, initialize_save::Bool, progress::Bool, progress_steps::Int64, progress_name::String, progress_message::typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE), userdata::Nothing, allow_extrapolation::Bool, initialize_integrator::Bool, alias_u0::Bool, alias_du0::Bool, initializealg::OrdinaryDiffEq.DefaultInit, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
@ OrdinaryDiffEq C:\Users\shiha.julia\packages\OrdinaryDiffEq\GxEzW\src\solve.jl:493
[7] CoupledODEs(prob::ODEProblem{Vector{Float64}, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, SciMLBase.AutoSpecialize, typeof(truncated), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, diffeq::NamedTuple{(:alg, :abstol, :reltol), Tuple{Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}, Float64, Float64}}; special_kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
@ DynamicalSystemsBase C:\Users\shiha.julia\packages\DynamicalSystemsBase\byMNf\src\core_systems\continuous_time_ode.jl:94
[8] CoupledODEs(prob::ODEProblem{Vector{Float64}, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, SciMLBase.AutoSpecialize, typeof(truncated), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, diffeq::NamedTuple{(:alg, :abstol, :reltol), Tuple{Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}, Float64, Float64}})
@ DynamicalSystemsBase C:\Users\shiha.julia\packages\DynamicalSystemsBase\byMNf\src\core_systems\continuous_time_ode.jl:89
[9] CoupledODEs(f::Function, u0::Vector{Float64}, p::Vector{Float64}; t0::Int64, diffeq::NamedTuple{(:alg, :abstol, :reltol), Tuple{Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}, Float64, Float64}})
@ DynamicalSystemsBase C:\Users\shiha.julia\packages\DynamicalSystemsBase\byMNf\src\core_systems\continuous_time_ode.jl:85
[10] CoupledODEs(f::Function, u0::Vector{Float64}, p::Vector{Float64})
@ DynamicalSystemsBase C:\Users\shiha.julia\packages\DynamicalSystemsBase\byMNf\src\core_systems\continuous_time_ode.jl:80
[11] top-level scope
@ In[31]:1
[12] eval
@ .\boot.jl:368 [inlined]
[13] include_string(mapexpr::typeof(REPL.softscope), mod::Module, code::String, filename::String)
@ Base .\loading.jl:1428

I expect I am making an elementary error since I am very new to this package. I have mixed up the syntax from quite a few different packages here since I am not too sure how to proceed without separating the second order derivatives by hand (I am using Plots, ModelingToolkit, DifferentialEquations, LaTeXStrings, DynamicalSystems; apart from the latter, the rest I have used earlier to numerically solve this very system).
Thanks in advance!