%
%  lect24.m:  VN analysis for leapfrog
%             (uncommented script)
%

N = 16;
L = pi;  c = 1;

dx = 2*L/N;

xj = 0:dx:2*L;

figure(1)
clf;subplot(1,2,1)
hold on
for nn = -N/2:1:N/2;

alpha = nn * (pi/L);
plot([0 2*pi],[nn nn],'k')
plot(xj,nn + real(0.3*exp(i*alpha*xj)),'o')
plot(xj,nn + real(0.3*exp(i*alpha*xj)),'r')

end

axis([0 2*L -N/2-1 N/2+1])
title('\bf VN modes for N=16 (real part)')
ylabel('\bf u_j(\alpha) - shifted vertically')
xlabel('\bf x_j')

subplot(1,2,2)
hold on
for nn = -N/2:1:N/2;

alpha = nn * (pi/L);
plot([0 2*pi],[nn nn],'k')
plot(xj,nn + imag(0.3*exp(i*alpha*xj)),'o')
plot(xj,nn + imag(0.3*exp(i*alpha*xj)),'r')

end

axis([0 2*L -N/2-1 N/2+1])
title('\bf VN modes for N=16 (imag part)')
ylabel('\bf u_j(\alpha) - shifted vertically')
xlabel('\bf x_j')

figure(2)
clf;subplot(1,2,1)
hold on
for nn = 0:1:N;

alpha = nn * (pi/L);
plot([0 2*pi],[nn nn],'k')
plot(xj,nn + 0.3*real(exp(i*alpha*xj)),'o')
plot(xj,nn + 0.3*real(exp(i*alpha*xj)),'r')

end

axis([0 2*L -1 N+1])
title('\bf aliasing of VN modes (real)')
ylabel('\bf u_j(\alpha) - shifted vertically')
xlabel('\bf x_j')

subplot(1,2,2)
hold on
for nn = 0:1:N;

alpha = nn * (pi/L);
plot([0 2*pi],[nn nn],'k')
plot(xj,nn + imag(0.3*exp(i*alpha*xj)),'o')
plot(xj,nn + imag(0.3*exp(i*alpha*xj)),'r')

end

axis([0 2*L -1 N+1])
title('\bf aliasing of VN modes (imag)')
ylabel('\bf u_j(\alpha) - shifted vertically')
xlabel('\bf x_j')

     
nn = -N/2:1:N/2-1;
alpha = nn * (pi/L);


figure(3);clf
subplot(2,1,1)
plot(alpha,c*alpha,'x')
hold on
dt = 1.0*dx;  lambda = c*dt/dx;
plot(alpha,real(asin(lambda*sin(alpha*dx)))/dt,'ro')
%plot(alpha,imag(asin(lambda*sin(alpha*dx)))/dt,'r.')
title('\bf \omega_1(\alpha) = sin^{-1}(\lambda sin(\alpha x))')
ylabel('\bf \omega_1')
subplot(2,1,2)
plot(alpha,c*alpha,'x')
hold on
dt = 0.75*dx;  lambda = c*dt/dx;
omega = asin(lambda*sin(alpha*dx))/dt;
plot(alpha,real(asin(lambda*sin(alpha*dx)))/dt,'ro')
%plot(alpha,imag(asin(lambda*sin(alpha*dx)))/dt,'r.')
xlabel('\bf \alpha')
ylabel('\bf \omega_1')