How to create SecondOrderODEProblem from ODESystem

Specific Domains
1

I want to simulate a single body orbiting (0, 0) by VelocityVerlet() and using a system defined through ModelingToolkit. Functional code without modeling toolkit is as follows:

using DifferentialEquations
function ddu_kepler(du, u, p,t)  # In this form, u is r (position vector)
    rx, ry = u
    pos_mag³ = hypot(rx, ry)^3
    ddu = [-rx/pos_mag³, -ry/pos_mag³]
    return ddu
end

r0 = [1.0, 0.0]
v0 = [0, 1/√r0[1]]
tspan = (0.0, 100.0)
prob = SecondOrderODEProblem(ddu_kepler, v0, r0, tspan)
sol = solve(prob_partitioned, VelocityVerlet(), dt=0.1)

My best attempt with ModelingToolkit is as follows:

using ModelingToolkit
@variables t rx(t) ry(t)
∂ₜ = Differential(t)
∂ₜ² = ∂ₜ^2

@named kepler_eqns = ODESystem([∂ₜ²(rx) ~ -rx/hypot(rx, ry)^3,
                                ∂ₜ²(ry) ~ -ry/hypot(rx, ry)^3]
)
kepler_eqns_simplified = structural_simplify(kepler_eqns)
u0 = [rx=>1, ry=>0, ∂ₜ(rx)=>0, ∂ₜ(ry)=>1]
tspan = (0.0, 10.0)

Which can be solved by non-second order solvers:

prob = ODEProblem(kepler_eqns_simplified, u0, tspan)
sol = solve(prob, Tsit5())

However, I am unable to define a SecondOrderODEProblem from the ODESystem:

prob = SecondOrderODEProblem(kepler_eqns_simplified, u0, tspan)
prob = SecondOrderODEProblem(kepler_eqns_simplified, u0, tspan, [])
prob = SecondOrderODEProblem(kepler_eqns, u0, tspan)
prob = SecondOrderODEProblem(kepler_eqns, u0, tspan, [])

results in the following errors (it seems strange that error is different when no parameters are passed explicitly):

image
image1806×446 102 KB

Is this not implemented, or am I doing it wrong?

2

You cannot. There is information lost if you have already lowered to a first order ode. Mtk would need to allow for building it directly. Worth an issue.

1 Like
3

Thanks for the quick reply.