The equations of the Marchuk immunology model discussed in
[6] are


yplot = sol.y; yplot(1,:) = 1e4*yplot(1,:);and so forth. To solve this problem accurately over the whole interval, you will need to reduce the tolerances to, say, a relative tolerance of 10^{5} and an absolute tolerance of 10^{8}.
dde23 is sufficiently robust that it can solve this problem in a straightforward way. However, you can be much more confident of the results if you ensure that the solver is always working with equations that have smooth coefficients. There are two kinds of discontinuities in this problem. As in Example 5, use the 'Jumps' option to tell the solver about the discontinuity at t = 10^{6}. As in Example 7, terminate the integration when the event function y_{4}(t)  0.1 vanishes. Use a parameter state with value +1 if y_{4}(t) £ 0.1 and 1 otherwise. The problem is to be solved with y_{4}(0) = 0, so initialize state to +1. Thereafter, each time that the solver returns, check whether you have reached the end of the interval. If sol.x(end) < 60, change the sign of state and call dde23 again with the previous solution as history. In the function for evaluating the DDEs, set x(y_{4}) = 1 if state is +1 and x(y_{4}) = (1 y_{4}) 10/9 otherwise.