L=30.;
N=200;
Nt=400;
x=linspace(-L,L,N+1);
dx=2*L/N;
dt=0.5*dx;
v=.4;
g=1/sqrt(1-v^2);
x0=L*.4;
phi=zeros(Nt,N+1);
phi(1,:)=tanh(g*(x-x0)/sqrt(2))-tanh(g*(x+x0)/sqrt(2))+1;
phi(2,:)=tanh(g*(x+v*dt-x0)/sqrt(2))-tanh(g*(x-v*dt+x0)/sqrt(2))+1;
for k=1:Nt-2
    phi(k+2,:)=2*phi(k+1,:)-phi(k,:)+dt^2/dx^2*(phi(k+1,:)*shl(N)...
        -2*phi(k+1,:)+phi(k+1,:)*shr(N))+dt^2*(phi(k+1,:)-phi(k+1,:).^3);
end
mesh(-phi)