function sol = prob4
% prob4, prob4f of tutorial.

b = 2;
c = 1;
e = 1 - c/b;   % Equilibrium point.

sol = dde23(@prob4f,7,0.8,[0, 60],[],b,c);

plot(sol.x,sol.y,[0, 60],[e, e],'r')
legend('y(t)','1 - c/b')
title(['Problem 4. Cooke Epidemic Model with' ...
       ' b = ',num2str(b),' , c = ',num2str(c),'.'])
xlabel('time t')

%=============================================================
function yp = prob4f(t,y,Z,b,c)
yp = b * Z * (1 - y) - c * y;