This Matlab code illustrates the dependence on paramenter analysis for a discrete dymanical system via repeated simulations.
clear all, close all, clc endtime = []; winner = []; values = [1 1.5 2 3 5 6]; % choose a set of values for lambda
Plot the direction field (no rescaling!)
for lambda = values a2=0.05; b2=0.005; a1=lambda*a2; b1 =lambda*b2; figure M=10; % number of samples points x1min=0; x1max=5; % domain specification x2min=0; x2max=2; [X1,X2]=meshgrid(x1min:(x1max-x1min)/M:x1max,x2min:(x2max-x2min)/M:x2max); dX1=-a1*X2-b1*X1.*X2; % x1-component dX2=-a2*X1-b2*X1.*X2; % x2-component quiver(X1,X2,dX1,dX2,0); % matlab routine axis([x1min x1max x2min x2max]); hold on xlabel('Red Army'); ylabel('Blue Army'); x=[5;2]; n=1; sprintf(strcat('Battle evolution for lambda =',num2str(lambda),'\n')) fprintf('n Red Army Blue Army \n') fprintf('%2.0f %5.2f %5.2f\n', 0, x(1), x(2)) while x(1)>0 & x(2)>0 n=n+1; rhs = [ -a1*x(2)-b1*x(1)*x(2); -a2*x(1)-b2*x(1)*x(2)]; xnew = x + rhs; plot([x(1),xnew(1)], [x(2),xnew(2)],'--ro',... 'MarkerFaceColor','k','MarkerSize',4) x = xnew; fprintf('%2.0f %5.2f %5.2f\n', n, x(1), x(2)) end if x(1)>0 winner = [winner,1]; else winner = [winner,2]; end title(strcat('Battle outcome for \lambda =',num2str(lambda))) xlabel('Red Army'); ylabel('Blue Army'); endtime = [endtime n]; end %endtime=endtime(2:end);
ans = Battle evolution for lambda =1 n Red Army Blue Army 0 5.00 2.00 2 4.85 1.70 3 4.72 1.42 4 4.62 1.15 5 4.54 0.89 6 4.47 0.64 7 4.42 0.40 8 4.40 0.17 9 4.38 -0.05 ans = Battle evolution for lambda =1.5 n Red Army Blue Army 0 5.00 2.00 2 4.78 1.70 3 4.59 1.42 4 4.43 1.16 5 4.31 0.91 6 4.21 0.68 7 4.14 0.45 8 4.09 0.24 9 4.06 0.03 10 4.06 -0.18 ans = Battle evolution for lambda =2 n Red Army Blue Army 0 5.00 2.00 2 4.70 1.70 3 4.45 1.43 4 4.24 1.17 5 4.08 0.93 6 3.95 0.71 7 3.85 0.50 8 3.78 0.30 9 3.74 0.10 10 3.72 -0.09 ans = Battle evolution for lambda =3 n Red Army Blue Army 0 5.00 2.00 2 4.55 1.70 3 4.18 1.43 4 3.87 1.19 5 3.63 0.98 6 3.43 0.78 7 3.27 0.59 8 3.15 0.42 9 3.07 0.26 10 3.02 0.10 11 3.00 -0.05 ans = Battle evolution for lambda =5 n Red Army Blue Army 0 5.00 2.00 2 4.25 1.70 3 3.64 1.45 4 3.15 1.24 5 2.74 1.07 6 2.40 0.91 7 2.12 0.78 8 1.88 0.67 9 1.68 0.57 10 1.52 0.48 11 1.38 0.40 12 1.26 0.33 13 1.17 0.26 14 1.10 0.20 15 1.04 0.15 16 1.00 0.09 17 0.98 0.04 18 0.96 -0.01 ans = Battle evolution for lambda =6 n Red Army Blue Army 0 5.00 2.00 2 4.10 1.70 3 3.38 1.46 4 2.79 1.27 5 2.31 1.11 6 1.90 0.98 7 1.55 0.88 8 1.25 0.79 9 0.98 0.73 10 0.74 0.67 11 0.52 0.63 12 0.32 0.61 13 0.14 0.59 14 -0.04 0.58
disp('lambda battle length winner')
disp([values' endtime' winner'])
lambda battle length winner 1.0000 9.0000 1.0000 1.5000 10.0000 1.0000 2.0000 10.0000 1.0000 3.0000 11.0000 1.0000 5.0000 18.0000 1.0000 6.0000 14.0000 2.0000