DifferentialEquations.jl: Problems with event detection

Specific Domains Numerics
,
1

Hi folks,

I am having problems getting event detection to work properly in DifferentialEquations.jl. I have reduced the code as far as possible:

using DifferentialEquations
using PyPlot

function evaluate!(t, y, δy)
    rm = norm(y[1:3])
    δy[1:3] = y[4:6]
    δy[4:6] = -mu .* y[1:3] ./ rm^3
end

function condition(t, y, integrator)
    r = y[1:3]
    v = y[4:6]
    rm = norm(r)
    vm = norm(v)
    ecc = norm(((vm^2 - mu/rm)*r - v * (r ⋅ v))/mu)
    sma = -mu/(2*(vm^2/2 - mu/rm))
    ese = (r ⋅ v) / sqrt(mu * sma)
    ece = rm * vm^2 / mu - 1.0
    E = atan2(ese, ece)
    ano = 2*atan2(sqrt(1 + ecc)*sin(E/2), sqrt(1 - ecc)*cos(E/2))
    push!(ts, t)
    push!(as, ano)
    ano
end

function affect!(integrator)
    push!(te, integrator.t)
end

const mu = 398600.4418

ts = Vector()
as = Vector()
te = Vector()

t1 = 5562.095583019713

r = [
6068279.27,
-1692843.94,
-2516619.18,
]/1000

v = [
-660.415582,
5495.938726,
-5303.093233,
]/1000

y0 = [r; v]

prob = ODEProblem(
    evaluate!,
    y0, (0.0, t1),
    callback=ContinuousCallback(condition, affect!),
)
res = solve(prob, Vern9())

idx = sortperm(ts)
t = ts[idx]
ano = as[idx]
plot(t, ano)
for t in te
    plot([te, te], [-pi, pi], "r")
end
grid()

This is the output:

The blue line is the output of the condition function. I would expect the callback to only be triggered at the upcrossing between 4000 and 5000 seconds. But somehow it also gets triggered at the downcrossing. Am I doing something wrong?

2

This means it will trigger at upcrossings and downcrossings.

cb = ContinuousCallback(condition,affect!) # Both upcrossings and downcrossings, the default
cb = ContinuousCallback(condition,affect!,nothing) # Only upcrossings
cb = ContinuousCallback(condition,nothing,affect!) # Only downcrossings

It is mentioned further down the page:

Notice that passing only one affect! is the same as ContinuousCallback(condition,affect!,affect!), i.e. both upcrossings and downcrossings will activate the event. Using ContinuousCallback(condition,affect!,nothing)will thus be the same as above because the first event is an upcrossing.

but reading the top, this seems to be a bit more confusing than I remember it being. Would you mind doing a PR for the docs to make it clearer?

3

Thanks!

Here is the PR: https://github.com/JuliaDiffEq/DiffEqDocs.jl/pull/31

There is a sarcastic German proverb that roughly translates to:
“The person who is able to read always has the upper hand.” :joy:

1 Like