DiffEqs: Hybrid (continuous+discrete) system / periodic callback

Just wanted to try the MWE also with the “clean” solution of subtyping DEDataVector and letting DiffEq do everything (saving the time trajectories, interpolation, etc.). Besides the problems I described here (What does "user_cache: method has not been implemented for the integrator" mean? - #3 by asprionj) and the necessity to loop over the cache, there’s another issue with the interpolation of the “discrete” user fields. Here’s the code for the MWE:

mutable struct CtrlSimTypeScalar{T} <: DEDataVector{T}
    x::Array{T,1}
    ctrl_x::T # controller state
    ctrl_u::T # controller output
end

# NOTE: non-in-place formulation does not work; should mention this in docs
function f(dx,x,p,t)
    dx[1] = x.ctrl_u - x[1]
end

ctrl_f(e, x) = 0.85(e + 1.2x)

x0 = CtrlSimTypeScalar([0.0], 0.0, 0.0)
prob = ODEProblem(f, x0, (0.0, 8.0))

function ctrl_fun(int)
    e = 1 - int.u[1] # control error

    # pre-calculate values to avoid overhead in iteration over cache
    x_new = int.u.ctrl_x + e
    u_new = ctrl_f(e, x_new)

    # iterate over cache...
    for c in full_cache(int)
        c.ctrl_x = x_new
        c.ctrl_u = u_new
    end
end

integrator = init(prob, Tsit5()) # only Tsit5 seems to work

# take 8 steps with a sampling time of 1s
Ts = 1.0
for i in 1:8
    ctrl_fun(integrator)
    step!(integrator,Ts,true)
end

Then, for illustrating what I mean, do some plotting:

# extract solution, plot...
sol = integrator.sol
plot(sol, label="system")

# actual full steps, with steppre line-style ("correct")
plot!(sol.t, [u.ctrl_u for u in sol.u], linetype=:steppre, color=:black, marker=true, label="u[.].ctrl_u")

# interpolated control signal ("wrong")
t_ctrl = 0:0.01:sol.t[end]
plot!(t_ctrl, [sol(t).ctrl_u for t in t_ctrl], color=:red, label="u(t).ctrl_u")

newplot(1)

It seems that:

  1. The cache does not contain the current stop, because changing the cache only affects future points. Is this always the case? Is there a closer description of the role and inner workings of the cache?

  2. The interpolation method used by DiffEq for the user fields is “last value”, which in this context is wrong. Can this behaviour be changed by some kwarg or similar? If not, would adding this be a simple addition or would we have to change the internals a lot?

Why bother? Isn’t this just cosmetics? Nope: When calculating intermediate quantities from the overall solution (that is, continuous state and the discrete signals), the correct value for the discrete signals at any given time has to be used. And I’m always against implementing a “core functionality” such as interpolating the solution again “on the outside”.

1 Like