Doing repeated Euler-steps in a simulation, Differentialequations.jl

Specific Domains Modelling & Simulations
1

Hi!

I am researcher working with problems related to optimal stopping/switching problems, and I have an issue where I need to, in a for loop, advance one time step (size dt) of an SDE for lots of previous states. And I am doing this recursively in a dynamic programming way, so I cannot directly generate complete trajectories here. In the end I want to do further computations on the results of the below, repeatedly in the for loop.

So far I have been doing it using a function like this:

function cpu_euler_maryama_alt(X_prev,time,dt,StepProc,p)
    d,M = size(X_prev)
    dX = [col .+ StepProc.drift(col,p,time) .* dt .+ StepProc.dispersion(col,p,time) * randn(Float32,d)*sqrt(dt)  for col in eachcol(X_prev)]
    return reduce(hcat,dX)
end

where StepProc is a struct containing the drift and dispersion of the SDE This works fine and is pretty efficient. I am trying to see if I can use Differentialequations.jl to accomplish the same thing, but I end up having to do a lot of allocations. Here is an example of what I have tried

function one_step(SDEprob,v0,dt,t)
    function start_values(prob,i,repeat)
        remake(prob,u0=v0[:,i],tspan=(t,t+dt))
    end
    d,M = size(v0)
    ensembleprob = EnsembleProblem(SDEprob,prob_func=start_values)
    sol = solve(ensembleprob,EnsembleThreads(),trajectories=M,dt=dt,saveat=[t+dt])

end

and above SDEprob is an instance of SDEProblem.

Do you know if it seems likely to do this efficiently with Differentialequations.jl?

2

https://docs.sciml.ai/DiffEqDocs/stable/basics/integrator/

Or use discrete callbacks

3

Thank you for the suggestions, a few questions though.

So far I have been trying to solve this using an EnsembleProblem as above, but it does not seem that init(prob) works if prob is an EnsembleProblem?

Further, I am doing this backwards in time, so I start at time N, perform an Euler-step for a set spatial points, do some computations (including solving a regression problem), then go back one time step to N-1, and repeat the process to N-2, N-3, …
The specific set of spatial points I have been using are coming from pre-generated paths of the same SDEprocess, but like I said, my problem right now is going backwards in time and doing these computations. Is this compatible with integrator or callbacks?

4

It currently does not.

Just flip the tspan, i.e. instead of (0,n) do (n,0).

5

Again, thank you for your comments and suggestions!

OK, just to be clear. Saying I am going backwards in time is a bit ambiguous, and I think I was not so clear. So this is what is happening:

I generate a path from time 0 to N, as usual, saving paths at timepoints 0,1,2, …, N-1, N. Then I want to go “backwards”, starting at time N, doing a FORWARD Euler-step (or similar) to time N+1 and do calculations. Then I want to go to time N-1, do a forward step to time N and do computations, and so on, going back to N-2, forward step to N-1 and then calculations.

And the points I use (so the spatial points that I input to the Euler step) to go forward with each each Euler step are the ones I began by generating in the usual way.

I have a hard time wrapping my head around about what happens if we flip the tspan around. Does this try to define some backward stochastic differential equation based on the normal SDE provided, or does it just do a normal Euler step with negative dt?

6

For the you also want to reverse the noise process. See https://github.com/SciML/DiffEqNoiseProcess.jl/blob/master/test/reversal_test.jl