Solving stiff SDDE problems

[anirudhpammi]

I tried using this package to solve a stiff problem, I am only able to use the EM or LambaEM solvers which result in some instability. I want to use stiff solvers like ImplicitEM or ImplicitRKMil. Can you help me out? Thanks!

[ChrisRackauckas]

Minimal example?

[anirudhpammi]

using BenchmarkTools, PyPlot, DifferentialEquations, StochasticDelayDiffEq

function yamada_norm!(du,u,p,t)
#no noise ODE
b1,b2,mu1,mu2,s=p;
du[1] = (b1mu1) - ( (b1u[1])(1+u[3]) ) ;
du[2] = (b2mu2) - ( (b2u[2])( 1 + (s*u[3]) ) );
du[3] = ( u[1] - u[2] - 1 )*u[3] ;

end

function yamada_dde(du,u,h,p,t)

b1, b2, mu1, mu2, s, kappa, tau = p;
histI = h(p,t-tau)[3];
du[1] = (b1*mu1) - ( (b1*u[1])*(1+u[3]) ) ;
du[2] = (b2*mu2) - (  (b2*u[2])*( 1 + (s*u[3]) )  );
du[3] = ( u[1] - u[2] - 1 )*u[3] + kappa*(histI);

end

function yamada_sigDDE(du,u,h,p,t)
b1, b2, mu1, mu2, s, kappa, tau = p;
du[1] = b1*0.01;
du[2] = 0;
du[3] = 0;
end

function yamada_sigODE(du,u,p,t)
# noise component of ODE
b1, b2, mu1, mu2, s = p;
du[1] = b1*0.01;
du[2] = 0;
du[3] = 0;

end

b1 = 0.01;
b2 = 0.02;
mu1 = 2.4;
mu2 = 2.2;
s = 5;
kappa = 0.05;
tau = 1100.0;
RT_max = 100;
tspan_hist = (0.0,tau);
tspan_dde = (tau,tau*RT_max)
lags = [tau];
u0_hist = [mu1; mu2;3];
u0_dde = [mu1; mu2;0];
p_hist = [b1,b2,mu1,mu2,s];
p_dde = [b1,b2,mu1,mu2,s,kappa,tau];

u0_hist = [mu1; mu2;3];

p_hist = [b1,b2,mu1,mu2,s];
p_dde = [b1,b2,mu1,mu2,s,kappa,tau];

prob_SDE = SDEProblem(yamada_norm!,yamada_sigODE,u0_hist,tspan_hist,p_hist);
sol_sde = solve(prob_SDE,ImplicitEM());

u0_dde = [mu1; mu2;0];
h(p,t) = sol_sde(t);

prob = SDDEProblem(yamada_dde, yamada_sigDDE, u0_dde, h, tspan_dde, p_dde; constant_lags=tau);

sol = solve(prob,ImplicitEM());

[HTSykora]

Hi,

SROCK methods are working, and they are also stabilized methods. The necessary pull request that make the implicit methods work are no yet included, but they are already done. So until this package is not integrated into the package, I recommend, that you use SROCK methods.

If you desperately need implicit methods, hit me up on the Julia slack channel, and I'll help you make the implicit method work, but unfortunately I will be available only from Monday.

[anirudhpammi]

Hi,
Thanks, the S-ROCK methods work and I don't seem to have the instabilities now. While I don't need the implicit method desperately, could you please tell me if I will have any significant difference/speedup if I use ImplicitEM() v/s SROCK1()? Thanks!

[ChrisRackauckas]

The SROCK are better if you only have large real eigenvalues because they are stable (for real eigenvalues) while being explicit. You only need to go implicit when you have large complex eigenvalues.

[ChrisRackauckas]

Looks like this is all good now.