DifferentialEquations.jl saving data without redundant calculation of control inputs

General Usage
1

I’m trying to save some data, which is also calculated within given ODE.

For example, with a linear system, I usually write a function

A = Matrix(I, n, n)
B = Matrix(I, n, m)
function dynamics!(dx, x, p, t; u)
    dx .= A*x + B*u
end

and apply helper functions provided in ComponentArrays.jl as

prob = ODEProblem(apply_inputs(dynamics!; u=(x, p, t) -> -K*x), x0)

with given gain matrix K.
Then, when saving input data, I usually use SavingCallback in DifferentialEquations.jl with affect! like

function affect!(x, t, integrator)
    u = -K*x
    (; x=x, u=u)
end

I think that in this case I have to calculate u twice at each instant when saving data.

  1. Is it correct?
  2. If so, how can I avoid redundant calculation?
2

Check out SimulationLogs.jl. You could change your dynamics! function to:

function dynamics!(dx, x, p, t; u)
    @log u
    dx .= A*x + B*u
end

Then access the control vector at each time point from the solution as

solution_log = get_log(sol)
solution_log.u

You can also interpolate at different time points with

solution_log_interpolated = get_log(sol, time_points)

Note that “logging” actually works by going back and re-calculating values, so there is some redundant calculation there. But it’s generally a pretty small price to pay compared to the running time of the simulation.

1 Like
3

Thanks for sharing a nice package!

Can I export some dynamical equations with @log from a custom package?
For example,

using MyPackage: predefined_dynamics!  # with @log k = 1
using SimulationLogs


prob = ODEProblem(predefined_dynamics!, x0, tspan)
sol = solve(prob)
out = get_log(sol)
@show out.k  # should be an array filled with 1

Sorry for a weird question if I understood incorrectly (I’m not familiar with that package :slight_smile: ).

4

Yep, that should work as long as MyPackage is using SimulationLogs and has the @log k = 1 defined somewhere where it will be called by predefined_dynamics!. For example:

function predefined_dynamics!(dx, x, p, t)
    @log k = 1
    @log u = -k*x
    dx .= A*x + B*u
end

or

function unity_gain(x, p, t)
    @log k = 1
    return k
end

function predefined_dynamics!(dx, x, p, t)
    k = unity_gain(x, p, t)
    @log u = -k*x
    dx .= A*x + B*u
end

It wouldn’t work, though, if k is defined statically somewhere in MyPackage with a @log in front of it, because then it wouldn’t be called again by get_log.

5

Awesome!
But as it re-calculate values after simulation, some stochastic behaviour would be hard to be logged.

As a remedy, it would be possible to update (stochastic) parameters and log the parameters.
To reproduce the same result, rng of the stochastic parameter update callback should be reset.

Am I right? Or, do you have a good idea for this?

6

Ah, I see. Unfortunately, since your dynamics function is being called by the solver multiple times per simulation time step, the random values generated while logging will get out of sync with those generated during the initial simulation. I’m not sure what the right solution is here.

1 Like