%
%  code02.m -- airy & asymptotics plot (djm - 12 sept 2001)
%
%- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

zmin = -25;  zmax = 10;  dz = 0.1;
z = zmin:dz:zmax;

%  basic airy plot

clf;  subplot(2,1,1)
plot(z,real(airy(z)),'k')
hold on
plot([zmin zmax],[0 0])
axis([zmin zmax -0.6 0.6])

%  labels

title('\bf airy function & large argument asymptotic forms')
xlabel('\bf z-axis')
ylabel('\bf Ai(z)')

disp('return for more')
pause

%
%  asymptotic plots (holmes, p298)
%

%  negative real (with envelope)

z_n    = zmin:dz:-dz;
zeta_n = (2/3)*(-z_n).^(3/2); 
eta_n  = (5/72)./zeta_n;
a_env  = (1/sqrt(pi)).*(-z_n).^(-1/4);
airy_n = a_env .* (cos(zeta_n -pi/4) + eta_n.*sin(zeta_n -pi/4));

plot(z_n,a_env,'r:');  plot(z_n,-a_env,'r:')
plot(z_n,airy_n,'r')

%  positive real

z_p    = dz:dz:zmax;
zeta_p = (2/3)*(z_p).^(3/2); 
eta_p  = (5/72)./zeta_p;
airy_p = (1/2/sqrt(pi)).*(z_p).^(-1/4).* ...
         exp(-zeta_p) .* (1-eta_p);

plot(z_p,airy_p,'r');  plot(z,real(airy(z)),'k')

disp('return for more')
pause

%
%  error plots
%

subplot(2,1,2)
plot(log(-z_n),log(abs(airy(z_n)-airy_n)),'m')
hold on
plot(log( z_p),log(abs(airy(z_p)-airy_p).*exp(zeta_p)),'g')
axis([0 3 -12 0])

plot([0 3],-2+[0 3*(-0.25-3)],'k--')
text(1.5,-6.2,'\bf slope = -3.25')

%  labels

title('\bf error plots')
xlabel('\bf log(|z|)')
ylabel('\bf log(error norm)')
text(0.5,-2,'\bf green (z>0) ; magenta (z<0)')