clear all; close all; clc;

L=20;
n=50;
x2=linspace(-L/2,L/2,n+1);
x=x2(1:n);
dx=x(2)-x(1);
dt=.01;
c=4;
CFL=dt/dx^2;
Time=4;
t=0:dt:Time;
time_steps=Time/dt;

u0=sech(x)';
u1=sech(x)';
usol(:,1)=u0;
usol(:,2)=u1;

e1=ones(n,1);
A=spdiags([e1 -2*e1 e1],[-1 0 1],n,n);
A(n,1)=1; A(1,n)=1;

for j=1:time_steps-1
    u2=u1-i*CFL*A*u1-i*2*dt*(conj(u1).*u1).*u1;  % fwd euler
    u1=u2;              % fwd euler
%   u2=u0 -i*CFL*A*u1-i*2*dt*(conj(u1).*u1).*u1;    % leap frog
%    u0=u1; u1=u2;       % leap frog
    usol(:,j+2)=u2;
end

waterfall(x,t,abs(usol)');