Defining a history function for Delay Differential Equations using past data points

I was trying to figure out how to define delay differential equations using DifferentialEquations.jl, but even after reading through the documentation, I am somewhat confused on how to define the history function. I am sure this is simpler than what it looks like, but I am confused mostly because the examples in the documentation only use historyfunction = 1 for all time < initial time. I am confused at the conceptual level, so I’m not yet providing a MWE here; I might do so if I understand it better and want to rephrase this question.

Assume we have some array with datapoints like pastdata = rand(10).
In this case, we can assume the datapoints are equidistant and sampled at times t=(0,1,2,…,9).
Let’s also assume a delay of tau=3.

At this stage, we would write our equation, in my case it’s an SDE with delay, but I believe this doesn’t matter for my question. So we are told we need to define a history function h that gets passed into the problem definition, with a signature similar to h(p,t) (for the out-of-place version).

I assume that internally the solver will only access this history function as long as the time t at which we are solving the equation allows for t-tau to be defined in the history function; this would mean that, if we solve the equation starting at t=10, the history function should stop being used when we reach t=13, because then the needed u(t=13-tau=10) would already have been simulated by the solver.

If so, how can we define our history function to make it so the solver reads the corresponding value from our pastdata array? We are told the interpolation happens internally, so defining the history function as an interpolation of our array seems to conflict with this. How can we ensure the solver knows that the value for t=1 here would be pastdata[2], for example? In the extreme case, I would imagine that if I defined tspan=(0,20), the solver would just give me back interpolated values of pastdata for all the t in [0,9], and simulate the leftover t’s until we reach t=20.

Sorry for the confusing write-up, but it’s confusing because I myself am confused.

You would just use an interpolation like DataInterpolations.jl for the h(p,t) that you give to the solver. That will then by used for any t < tspan[1].

All right, thanks for that! I will give a try. Does it need to be defined using a type from DataInterpolations.jl, or are the types from Interpolations.jl also ok?

That’s fine.