Solving DDE system using derivative at a previous time

Hi. I’m trying to figure out to use DifferentialEquations.jl to solve a system of DDEs that depend not only on the values at previous times, but on the time derivatives of those values at previous times. For example, suppose that instead of solving this system (the Lotka-Volterra system)

x'(t) = 0.5x(t)(1 - y(t-\tau)) \\ y'(t) = -0.5y(t)(1 - x(t-\tau))

we were trying to solve this system.

x'(t) = 0.5x(t)(1 - y(t-\tau)) \\ y'(t) = -0.5y(t)(1 - x'(t-\tau))

(note the x' in the second equation).

For the first system, we could set it up using DifferentialEquations like this:

function lv_model(du,u,h,p,t)
  tau = p
  hist1 = h(p, t-tau)[1]
  hist2 = h(p, t-tau)[2]
  du[1] = 0.5*u[1]*(1 - hist2)
  du[2] = -0.5*u[2]*(1 - hist1)
end

and that works fine when building the DDEProblem. How could the second system be constructed? Is there a way to get the derivatives at a previous time from h(p, t-tau)?

Thanks in advance for any help. I’m still relatively new to Julia.

h(p, t-tau, Val{1}) gives the first derivative

1 Like

Thanks so much for the quick response! Given how well-thought out everything in Julia seems to be, I assumed there’d be an easy answer and it looks like there is.

Hi Chris. I’m having some trouble when I drop in the Val{1}. It results in a lengthy error trace, but the relevant part at the beginning is

ERROR: MethodError: no method matching h(::Float64, ::Float64, ::Type{Val{1}})
Closest candidates are:
  h(::Any, ::Any) at <other location in code>
Stacktrace:
  [1] (::DelayDiffEq.Hi.......

Seems like the definition of h doesn’t expect the third argument. Do I need do do something different?

Here’s my full script.

using DifferentialEquations
function lv_model(du,u,h,p,t)
  tau = p
  hist1 = h(p, t-tau)[1]
  hist2 = h(p, t-tau)[2]
  dhist1 = h(p, t-tau, Val{1})[1]
  dhist2 = h(p, t-tau, Val{1})[2]

  du[1] = 0.5*u[1]*(1 - hist2)
  du[2] = -0.5*u[2]*(1 - hist1)
end

h(p, t) = [1.0,2.0]

tau = 0.2
lags = [tau]

p = (tau)
tspan = (0.0,30.0)
u0 = [1.0,2.0]

prob = DDEProblem(lv_model,u0,h,tspan,p; constant_lags=lags)

sol = solve(prob)

using Plots
plot(sol, idxs=(0,1))
plot!(sol, idxs=(0,2))

You need to implement it!
This says for all times -tau…0.0 and all parameter values hist1 is 1.0 and his2 is 2.0.

h(p, t) = [1.0,2.0]

You also need to write a function that specifies the first derivative with respect to time for all history and parameters.

h(p, t, ::Type{Val{1}}) = [0.0, 0.0] # h' is 0 since h was constant
1 Like

Thanks @contradict. That got me unstuck. This is also helping me to better understand the differences between Julia and Python.

For anyone who might benefit from this later, here’s a complete working example.

using DifferentialEquations
function lv_model(du,u,h,p,t)
  tau = p
  hist1 = h(p, t-tau)[1]
  hist2 = h(p, t-tau)[2]
  dhist1 = h(p, t-tau, Val{1})[1]
  dhist2 = h(p, t-tau, Val{1})[2]

  du[1] = 0.5*u[1]*(1 - dhist2)
  du[2] = -0.5*u[2]*(1 - hist1)
end

h(p, t) = [1.0,2.0]
h(p, t, ::Type{Val{1}}) = [0.0, 0.0]

tau = 0.2
lags = [tau]

p = (tau)
tspan = (0.0,30.0)
u0 = [1.0,2.0]

prob = DDEProblem(lv_model,u0,h,tspan,p; constant_lags=lags)

sol = solve(prob)

using Plots
# plot(sol, idxs=(1,2))
plot(sol, idxs=(0,1))
plot!(sol, idxs=(0,2))
1 Like

and note it is mentioned in the docs:

https://docs.sciml.ai/DiffEqDocs/stable/types/dde_types/#SciMLBase.DDEProblem

though briefly.

1 Like

Yeah, and I did actually find it after the fact, although I still probably would have needed to ask about how to use it. :slight_smile:

Thanks for making this great package, @ChrisRackauckas . Blows the Python DDE packages out of the water.