clear all
close all
%-----------------------------------------------------------------------

%
%------------------------------------------------------------------------

% parameters
L = 4*pi; % domain length 
N=1024;   % number of spatial sampling points
h=0.001;  % dt (time stepping interval)           
nu=0.02;  % physical constant
T=30000;  % end time
tt=100;   % output times
x=(2*pi/N)*(-N/2:N/2-1)'; % initialize spatial grid
% 
k=((2*pi/L)*[0:N/2 -N/2+1:-1]'); % array of k values corresponding to the fft of the spatial array
%

u=sin(x); %initial condition

% plotting routine initialisation
plot(x,u,'r','LineWidth',2)  
xlabel('x','FontName','Arial','FontSize',16)
ylabel('u(x,t)','FontName','Times New Roman','FontSize',16)
axis([-L/2 L/2 -1 1])
title(['Burgers Equation, t=' num2str(0)],'Color','k')
hold off
pause(0.1)
%
uhat=fft(u); % fft of initial condition

% two different ways of calculating the rhs in fourier space
%rhs = @(nu,k,u) -nu*k.^2.*u+fft(ifft(-i*k.*u).*ifft(u)); % RHS function
rhs = @(nu,k,u) -nu*k.^2.*u-0.5*i*k.*fft((ifft(u)).^2); % RHS function
%

% RK 4 routine
for n=1:T
    k1 =rhs(nu,k,uhat);
    k2 =rhs(nu,k,uhat+0.5*h*k1);
    k3 =rhs(nu,k,uhat+0.5*h*k2);
    k4 =rhs(nu,k,uhat+h*k3);
    uhat = uhat + (h/6)*(k1 + 2*k2 + 2*k3 + k4);
    %plotting routine
    if mod(n,tt)==0
        u=real(ifft(uhat));
        figure(1)
        plot(x,u,'r','Linewidth', 2)
        xlabel('x','FontName','Arial','FontSize',16)
        ylabel('u(x,t)','FontName','Times New Roman','FontSize',16)
        axis([-L/2 L/2 -1 1])
        title(['Burgers Equation, t=' num2str(n)],'Color','k')
        hold off
        pause(0.1)
    end;
end;