Jump Diffusion Process. When does a jump occur?

[SSkor]

I continue study SDE with Julia and I straggle with integrator

I’d like to get simulations like this

And I get something like that with this code:

rate(u,p,t) = 2
affect!(integrator) = (integrator.u[1] = 1)
jump = VariableRateJump(rate,affect!)
jump_prob = JumpProblem(prob,Direct(),jump)
monte_prob = MonteCarloProblem(jump_prob)
sol = solve(monte_prob,SRIW1(),num_monte=3,parallel_type=:threads,dt=0.001,
adaptive=false) 

But when I change (integrator.u[1] = 1) to (integrator.u[1] = integrator.u[1]/2) the 1st jump occurs in the same time for all trajectories.

What does an integrator mean for the jump process?

[ChrisRackauckas]

Does this happen when parallel_type=:none? Need to double check the thread-safety of the RNG in this case.

[SSkor]

Yes, this happens.
sol = solve(monte_prob,SRIW1(),num_monte=n,parallel_type=:none,dt=dt,adaptive=false)

[ChrisRackauckas]

Interesting. Can you get a simplified example which recreates this problem and post it in an issue to DiffEqJump.jl? Issues · SciML/JumpProcesses.jl · GitHub . It looks like the first jump is re-using the same random number each time. We have tests for that when doing looping but might have omitted the test when using a MonteCarloProblem, but that is obviously a good think to fix here.

[SSkor]

If you don’t mind, at first I post the code here. Maybe I do something wrong when I convert the solution into array.

using DifferentialEquations, Plots

rate(u,p,t) = 2
affect!(integrator) = (integrator.u[1] = integrator.u[1]/2)
jump = VariableRateJump(rate,affect!)
jump_prob = JumpProblem(prob,Direct(),jump)
monte_prob = MonteCarloProblem(jump_prob)
sol = solve(monte_prob,SRIW1(),num_monte=3,parallel_type=:none,dt=0.001,adaptive=false)

#I need the trajectories of equal lengths, so I do that
A = Array{Float64}(undef, 1000, 3)
for i in 1:3
    A[:,i] = Array(sol[i])[1,1:1000]
end

pyplot() 
p = plot(A, legend = false, title="$proc_name $n")
display(p)

[ChrisRackauckas]

Oh, the jumps add a time point in there themselves since they introduce a discontinuity. If you don’t want that save, you do save_positions=(false,false) in VariableRateJump. Then Array(sol) will work and line up with the time points you’re expecting. I am not sure if the first jump is caused actually the same or if this rescaling you did made them look the same, so give this a shot first.

[SSkor]

I’ve changed this line
jump = VariableRateJump(rate, affect!, save_positions=(false,false))

I plot the solution without rescaling. Nothing changes.

[ChrisRackauckas]

Yup okay, this requires an issue. Somehow the callback init which is supposed to re-sample that random number isn’t being called. I think it’s confined to only variable rate jumps though.

[SSkor]

Okay, I’m posting an issue!

What the sense in this line for the jumps? affect!(integrator) = (integrator.u[1] = 1) or integrator.u[1] = integrator.u[1]
As I understand, it affects on the size of jumps at time points, doesn’t it? But in what time points?

[ChrisRackauckas]

I’m confused why that would matter too, hence the need to investigate