%
%  code9Wa.m -- 30 oct 02 -- djm
%
%  p385, ex1 plotting code
%

%  clear workspace & graphics ("help clear", etc)

clear;  close all

%  set vector of t-points

dt = 0.05;
t  = 0:dt:12;

%  calculate solutions

x1 =  exp(-t/2).*cos(t);
x2 = -exp(-t/2).*sin(t);

%  plot line graphs

figure(1);  
subplot(2,1,1);  hold on
plot(t,x1)

%  label plot ("\bf"-prefix makes text boldface)

title(['\bf x_1(t) for example 1, p385'])
xlabel('\bf t-axis')
ylabel('\bf x_1(t)-axis')

subplot(2,1,2);  hold on
plot(t,x2)

%  label plot ("\bf"-prefix makes text boldface)

title(['\bf x_2(t) for example 1, p385'])
xlabel('\bf t-axis')
ylabel('\bf x_2(t)-axis')

pause

%  plot phase plane

figure(2);  hold on
plot(x1,x2,'b')
for j=1:20:101
	plot(x1(j),x2(j),'rx')
end
 
axis equal
axis([-0.4 1 -0.8 0.6])

%  label plot ("\bf"-prefix makes text boldface)

title(['\bf phase plane for example 1, p385'])
xlabel('\bf x_1(t)-axis')
ylabel('\bf x_2(t)-axis')

pause

xt1 = -0.4:0.05:1.0;
xt2 = -0.8:0.05:0.6;

[xg1,xg2] = meshgrid(xt1,xt2);

xp1 = -0.5*xg1 + xg2;
xp2 = - xg1 - 0.5*xg2;

r = sqrt(xp1.^2 + xp2.^2);

quiver(xt1,xt2, xp1./r, xp2./r,0.25,'k')
quiver(xt1,xt2,-xp1./r,-xp2./r,0.25,'k.')