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 
 Function handle that returns a column vector of delays d(j).
The delays can depend on both t and y(t). If
all the delay functions have the form d(j) = t – τ_{j}, you
can set the argument 
 Specify

 Interval of integration from 
 Optional integration argument. A structure you create
using the 
sol = ddesd(ddefun,delays,history,tspan)
integrates
the system of DDEs
$${y}^{\prime}(t)=f\left(t,y(t),y(d(1)),\mathrm{...},y(d(k))\right)$$
on the interval [t_{0},t_{f}],
where delays d(j) can depend
on both t and y(t),
and t_{0} < t_{f}.
Inputs ddefun
and delays
are
function handles. See Create Function Handle 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 
 Approximation to y(x)
at the mesh points in 
 Approximation to y′(x)
at the mesh points in 
 Solver name, 
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 for
more information.
Commonly used options are scalar relative error tolerance 'RelTol'
(1e3
by
default) and vector of absolute error tolerances 'AbsTol'
(all
components are 1e6
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 k
th event function in events
:
value(k)
is the value of the k
th
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.
If you specify the 'Events'
option
and events are detected, the output structure sol
also
includes fields:
 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 
 Vector containing indices that specify which event occurred
at the corresponding time in 
The equation
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
plot(tint,yint);
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]);
For more examples of solving delay differential equations see ddex2
and ddex3
.
[1] Shampine, L.F., “Solving ODEs and DDEs with Residual Control,” Applied Numerical Mathematics, Vol. 52, 2005, pp. 113127.