Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Is there an idiomatic way to terminate integration after n callbacks in DifferentialEquations.jl

First of all, I am using the DifferentialEquations.jl library, which is fantastic! Anyway, my question is as follows:

Say for example, I have the following differential equation:

function f(du, u, t)
    du[1] = u[3]
    du[2] = u[4]
    du[3] = -u[1] - 2 * u[1] * u[2]
    du[4] = -u[2] - u[1]^2 + u[2]^2
end

and I have a callback which is triggered every time the trajectory crosses the y axis:

function condition(u, t, integrator)
    u[2]
end

However, I need the integration to terminate after exactly 3 crossings. I am aware that the integration can be terminated by using the effect:

function affect!(integrator)
    terminate!(integrator)
end

but what is the proper way to allow for a counting of the number of callbacks until the termination criterion is met. Furthermore, is there a way to extend this methodology to n events with n different counts?

In my research I often need to look at Poincare maps and the first, second, third, etc. return to the map so I am in need of a framework that allows me to perform this counting termination. I am still new to Julia and so am trying to reinforce good idiomatic code early on. Any help is appreciated and please feel free to ask for clarification.

like image 975
Rolfe Power Avatar asked Apr 09 '19 04:04

Rolfe Power


1 Answers

There is a userdata keyword argument to solve which can be useful for this. It allows you to pass objects to the integrator. These objects can be used in creative ways by the callback functions.

If you pass userdata = Dict(:my_key=>:my_value) to solve, then you can access this from integrator.opts.userdata[:my_key].

Here is a minimal example which controls how many times the callback is triggered before it actually terminates the simulation:

function f(du, u, t)
    du[1] = sin(t)
end
function condition(u, t, integrator)
    u[1] 
end

function affect!(integrator)
    integrator.opts.userdata[:callback_count] +=1
    if integrator.opts.userdata[:callback_count] == integrator.opts.userdata[:max_count]
        terminate!(integrator)
    end
end


callback = ContinuousCallback(condition, affect!)

u0 = [-1.]
tspan = (0., 100.)

prob = ODEProblem(f, u0, tspan)
sol = solve(prob; callback=callback, userdata=Dict(:callback_count=>0, :max_count=>3))
like image 193
Korsbo Avatar answered Oct 15 '22 12:10

Korsbo