%% Example 6.5
% Differentiation Using the Fourier Spectral Method and Second-Order Central Difference Formula

close all
clear example f fd x N dX X Y k YY YYd Ydl Y2 A Yd2

%% Choice of example a or b : exact solution
example = 'b';

switch example
    case 'a'
        f   = @(x) cos(3*x);
        fd  = @(x) -3*sin(3*x);
        
    case 'b'
        f   = @(x) 2*pi*x - x.^2;
        fd  = @(x) 2*pi - 2*x;
end;


%% Initialization :
x   = [0:0.01:2*pi];

N   = 8;
dX  = 2*pi/N;
X   = [0:dX:2*pi-dX]';
Y   = f(X);



%% Fourier Spectral Method :
k   = [[0:N/2-1]';0;[-N/2+1:-1]'];

YY  = fft(Y)/N;
YYd = i*k.*YY;
Yd1 = ifft(YYd)*N;


%% 2nd order Finite Difference :

 %extension of domain (assuming periodicity)
  Y2  = [f(2*pi-dX);Y;f(2*pi)];

 %differentiation operator for interior points
A   = diag(-ones(N+1,1),-1) + diag(ones(N+1,1),1);
A   = A/(2*dX);
Yd2 = A*Y2;
Yd2 = Yd2(2:N+1);



%% Plot
switch example
    case 'a'
           %exact
        plot(x/pi,fd(x),'k-','LineWidth',2);
        hold on;

           %spectral
        plot(X/pi,Yd1,'ko','MarkerFaceColor','k');

           %FD (for N=16 and N=8)
        plot(X/pi,Yd2,'k:^','MarkerFaceColor','k','LineWidth',2);
        %plot(X/pi,Yd2,'k-.v','MarkerFaceColor','k','LineWidth',2);

        axis([-0.1 2.2 -4 4]);
        
    case 'b'
           %exact
        plot(x/pi,fd(x),'k-','LineWidth',2);
        hold on;
        
           %spectral
        plot(X/pi,Yd1,'k:o','MarkerFaceColor','k','LineWidth',2);
        
           %FD
        plot(X/pi,Yd2,'k-.^','MarkerFaceColor','k','LineWidth',2);
        
        axis([-0.1 2.0 -10 10]);

end;

grid on;
