Hale [7] cites predatorprey models obtained by
introducing a resource limitation on the prey and assuming the
birth rate of predators responds to changes in the magnitude of
the population y_{1} of prey and the population y_{2} of
predators only after a time delay t. Starting with the system
of ODEs [15]


Recall that you solve ODEs with dde23 by setting lags to []. When this is done, the argument Z that dde23 supplies to the functions it calls is the empty array. You can use this to code the evaluation of both sets of equations in the same function by testing isempty(Z) to find out which set to evaluate. A more straightforward approach is to use two functions for the two sets of equations. Solve the DDE with t = 1. Plot in one figure y_{2}(t) against y_{1}(t) for the two solutions. This phase plane plot of the solution of the ODEs should be a closed curve corresponding to a limit cycle. To achieve this you will need to tighten the error tolerances with a command like
opts = ddeset('RelTol',1e5,'AbsTol',1e8);The figure makes the point that introducing a delay into an ODE model can have a profound effect on the solution. If you experiment with t, you will find this to be true even for small delays. It is also interesting to remove the resource term 1 y_{1}(t)/m and see how the orbits change as t is changed.