Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

FunctionCallingCallback not correctly synchronized with events #126

Open
MartinOtter opened this issue Oct 18, 2020 · 3 comments
Open

FunctionCallingCallback not correctly synchronized with events #126

MartinOtter opened this issue Oct 18, 2020 · 3 comments

Comments

@MartinOtter
Copy link

The FunctionCallingCallback seems to be not correctly synchronized with events as shown with the following test model:

module Test_bouncingBall_withFunctionCallBack

using DifferentialEquations

function f(du,u,p,t) #u[1]: height, u[2]:velocity
    du[1] = u[2]
    du[2] = -9.81
    return nothing
end


# Own result storage (emulated by print statements)
function outputs!(x, t, integrator)::Nothing
    println("... Store result at time = $t: h = ", integrator.u[1], ", v = ", integrator.u[2])
    return nothing
end


# State event
function condition(out,u,t,integrator)
    out[1] = u[1] 
    return nothing
end

function affect!(integrator, event_index)
    println("... State event at ", integrator.t)
    
    # Store result before event
    outputs!(integrator.u, integrator.t, integrator)
    
    # Process event
    integrator.u[2] = -0.8*integrator.u[2]
    auto_dt_reset!(integrator)
    set_proposed_dt!(integrator, integrator.dt) 
    
    # Store result after event
    outputs!(integrator.u, integrator.t, integrator)   
    return nothing    
end

u0     = [1.0,0.0]
tspan  = (0.0,0.6)
tspan2 = 0.0:0.1:0.6 
prob = ODEProblem(f,u0,tspan)

cb1 = FunctionCallingCallback(outputs!, funcat=tspan2)
cb2 = VectorContinuousCallback(condition,affect!,1)
callbackSet = CallbackSet(cb1,cb2)

sol = solve(prob, 
            Tsit5(), 
            save_everystep=false, 
            save_start=false, 
            save_end=false,
            adaptive=true,
            dt=0.1,
            callback=callbackSet)
end

This results in the following output:

... Store result at time = 0.0: h = 1.0, v = 0.0
... Store result at time = 0.1: h = 0.9509500000000003, v = -0.9809999999999998
... State event at 0.4515236409857323
... Store result at time = 0.4515236409857323: h = 1.6653345369377348e-15, v = -4.429446918070022
... Store result at time = 0.4515236409857323: h = 1.6653345369377348e-15, v = 3.543557534456018
... Store result at time = 0.2: h = 1.6653345369377348e-15, v = 3.543557534456018
... Store result at time = 0.3: h = 1.6653345369377348e-15, v = 3.543557534456018
... Store result at time = 0.4: h = 1.6653345369377348e-15, v = 3.543557534456018
... Store result at time = 0.5: h = 0.418002671515621, v = 2.087004452526053
... Store result at time = 0.6: h = 0.418002671515621, v = 2.087004452526053

This means the affect-function is called at time = 0.45, but the FunctionCallingCallback for time = 0.2, 0.3, 0.4 is called after the affect-function is called at time = 0.45.

@MartinOtter
Copy link
Author

This above result was generated with DifferentialEquations.jl version v6.15.0.

@MartinOtter
Copy link
Author

Unfortunately, the issue is worse than expected:

It does not matter whether results are stored at an event point (I thought this is a work-around): The FunctionCallingCallback is still called for time instants before the event, but using the state information returned from affect!(..). This means that the result storage of saveat points before the event is completely wrong.

@MartinOtter
Copy link
Author

Any progress here?

I have the following (bad) work-around:
At an event instant, I call outputs!(..) myself for all time instants stored in integrator.sol.t/sol.x from the previous event until the current time instant and store the result in my result data structure and only afterwards process the event.

I have now an application where my current work-around is not working: I am using IDA() to avoid solving large equation systems in the model. Before calling outputs!(..), I need to compute the derivatives by interpolation in the integrator (integrator(du, t, Val{1})). However, interpolation seems to be only possible within the current step:

# Error message:
[IDAS ERROR]  IDAGetDky
  Illegal value for t. t = 0.2 is not between tcur - hu = 0.67939 and tcur = 0.77939.  # event is at 0.75

┌ Warning: IDAGetDky failed with error code =
│   flag = -26
└ @ Sundials ...\packages\Sundials\k9hc3\src\simple.jl:20

and I need it for previous steps since outputs!(..) is called for previous steps.

It seems my only option is currently to call IDA(..) directly (with all the logic for event handling and communication points), and not using DifferentialEquations.

@ChrisRackauckas ChrisRackauckas transferred this issue from SciML/DifferentialEquations.jl Sep 9, 2022
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

2 participants