Solve delay differential equations (DDEs) with general delays
sol = ddesd(ddefun,delays,history,tspan)
sol = ddesd(ddefun,delays,history,tspan,options)
Function handle that evaluates the right side of the differential equations y′(t) = f(t,y(t),y(d(1),...,y(d(k))). The function must have the form
dydt = ddefun(t,y,Z)
where t corresponds to the current t, y is a column vector that approximates y(t), and Z(:,j) approximates y(d(j)) for delay d(j) given as component j of delays(t,y). The output is a column vector corresponding to f(t,y(t),y(d(1),...,y(d(k))).
Function handle that returns a column vector of delays d(j). The delays can depend on both t and y(t). ddesd imposes the requirement that d(j) ≤ t by using min(d(j),t).
If all the delay functions have the form d(j) = t – τj, you can set the argument delays to a constant vector delays(j) = τj. With delay functions of this form, ddesd is used exactly like dde23.
Specify history in one of three ways:
Interval of integration from t0=tspan(1) to tf=tspan(end) with t0 < tf.
Optional integration argument. A structure you create using the ddeset function. See ddeset for details.
on the interval [t0,tf], where delays d(j) can depend on both t and y(t), and t0 < tf. Inputs ddefun and delays are function handles. See the function_handle reference page for more information.
Parameterizing Functions explains how to provide additional parameters to the functions ddefun, delays, and history, if necessary.
ddesd returns the solution as a structure sol. Use the auxiliary function deval and the output sol to evaluate the solution at specific points tint in the interval tspan = [t0,tf].
yint = deval(sol,tint)
The structure sol returned by ddesd has the following fields.
Mesh selected by ddesd
Approximation to y(x) at the mesh points in sol.x.
Approximation to y′(x) at the mesh points in sol.x
Solver name, 'ddesd'
sol = ddesd(ddefun,delays,history,tspan,options) solves as above with default integration properties replaced by values in options, an argument created with ddeset. See ddeset and Types of DDEs in the MATLAB® documentation for details.
Commonly used options are scalar relative error tolerance 'RelTol' (1e-3 by default) and vector of absolute error tolerances 'AbsTol' (all components are 1e-6 by default).
Use the 'Events' option to specify a function that ddesd calls to find where functions g(t,y(t),y(d(1)),...,y(d(k))) vanish. This function must be of the form
[value,isterminal,direction] = events(t,y,Z)
and contain an event function for each event to be tested. For the kth event function in events:
value(k) is the value of the kth event function.
isterminal(k) = 1 if you want the integration to terminate at a zero of this event function and 0 otherwise.
direction(k) = 0 if you want ddesd to compute all zeros of this event function, +1 if only zeros where the event function increases, and -1 if only zeros where the event function decreases.
Row vector of locations of all events, i.e., times when an event function vanished
Matrix whose columns are the solution values corresponding to times in sol.xe
Vector containing indices that specify which event occurred at the corresponding time in sol.xe
sol = ddesd(@ddex1de,@ddex1delays,@ddex1hist,[0,5]);
solves a DDE on the interval [0,5] with delays specified by the function ddex1delays and differential equations computed by ddex1de. The history is evaluated for t ≤ 0 by the function ddex1hist. The solution is evaluated at 100 equally spaced points in [0,5]:
tint = linspace(0,5); yint = deval(sol,tint);
and plotted with
This problem involves constant delays. The delay function has the form
function d = ddex1delays(t,y) %DDEX1DELAYS Delays for using with DDEX1DE. d = [ t - 1 t - 0.2];
The problem can also be solved with the syntax corresponding to constant delays
delays = [1, 0.2]; sol = ddesd(@ddex1de,delays,@ddex1hist,[0, 5]);
or using dde23:
sol = dde23(@ddex1de,delays,@ddex1hist,[0, 5]);
 Shampine, L.F., "Solving ODEs and DDEs with Residual Control," Applied Numerical Mathematics, Vol. 52, 2005, pp. 113-127.