function sol = prob_4_5(G) % % Written by: % -- % John L. Weatherwax 2005-04-27 % % email: wax@alum.mit.edu % % Please send comments and especially bug reports to the % above email address. % %----- global c D n c = 0.5; D = 5; G = 1; n = 100; % x = y(1), y = y(2), lambda = y(3), % I_1 = y(4), I_2 = y(5), I_3 = y(6). y0 = [0.8*n; 0.2*n; 0; 0; 0; 0]; sol = dde23(@odes,[],y0,[0, D],[],G); sol = dde23(@ddes,D,sol,[D, 4*D],[],G); plot(sol.x,[sol.y(1:2,:); 100*sol.y(3,:)]); legend('x(t)','y(t)','100 \lambda(t)',0) disp( 'Equilibrium solution #1:' ); disp( [ n 0 0 ] ); disp( 'Equilibrium solution #2:' ); disp( [ G*n/(c*D) n*(1-G/(c*D)) c*D-G ] ); disp( 'The long time solution limit''s is:' ); sol.y(1:3,end) %============================================================ function dydt = odes(t,y,Z,G) global c D n dydt = zeros(6,1); dydt(1) = - y(1)*y(3) + G*y(2); dydt(2) = - dydt(1); dydt(4) = exp(G*t)*y(2); dydt(5) = t*exp(G*t)*y(1)*y(3); dydt(6) = exp(G*t)*y(1)*y(3); dydt(3) = (c/n)*exp(-G*t)*((dydt(4)+dydt(5))-G*(y(4)+y(5))); function dydt = ddes(t,y,Z,G) global c D n dydt = zeros(6,1); dydt(1) = - y(1)*y(3) + G*y(2); dydt(2) = - dydt(1); dydt(4) = exp(G*t)*y(2) - exp(G*(t - D))*Z(2); dydt(5) = D*exp(G*t)*y(1)*y(3) - y(6); dydt(6) = exp(G*t)*y(1)*y(3) - exp(G*(t - D))*Z(1)*Z(3); dydt(3) = (c/n)*exp(-G*t)*((dydt(4)+dydt(5))-G*(y(4)+y(5)));