% Showing the two-soliton solutions to the KdV equation

x=-6:0.06:6;
t=-0.3:0.003:0.3;
mx=min(x);
Mx=max(x);
mt=min(t);
Mt=max(t);
y=[mx Mx];

[X,T]=meshgrid(x,t);

set(gca,'FontSize',[16]);

Z12=12*(3+4*cosh(2*X-8*T)+cosh(4*X-64*T))./((3*cosh(X-28*T)+cosh(3*X-36*T)).^2);
v1=16;
v2=4;
Z1 = 8*sech(2*(X-v1*T)+0.5*log(3)).^2;
Z2 = 2*sech(X-v2*T-0.5*log(3)).^2;

MZ=max(max(Z1+Z2));

figure(1)
clf
subplot(2,1,1)
surfl(X,T,Z12);
shading interp
colormap(gray)
axis tight
ylabel('t')
xlabel('x')
zlabel('u(x,t)')
view([-35 65])
title('two-soliton solution')
subplot(2,1,2)
surfl(X,T,Z1+Z2);
shading interp
colormap(gray)
axis tight
ylabel('t')
xlabel('x')
zlabel('u(x,t)')
view([-35 65])
title('(one+one)-soliton solution')

figure(2)
clf
subplot(2,2,1)
surf(X,T,(Z12));
MZ12=max(max(Z12));
shading interp
axis tight
colormap('default')
ylabel('t')
xlabel('x')
zlabel('u(x,t)')
view([0 90])
colorbar
hold on
plot3(y-0.25*log(3),y/v1,[1 1]*MZ12,'k--','LineWidth',2);
plot3(y+0.5*log(3),y/v2,[1 1]*MZ12,'k--','LineWidth',2);
title('two-soliton solution')
axis([mx Mx mt-0.0*(Mt-mt) Mt+0.0*(Mt-mt)]);

subplot(2,2,3)
surf(X,T,(Z1+Z2));
MZ12=max(max(Z1+Z2));
shading interp
colormap('default')
ylabel('t')
xlabel('x')
zlabel('u(x,t)')
view([0 90])
colorbar
hold on
plot3(y-0.25*log(3),y/v1,[1 1]*MZ12,'k--','LineWidth',2);
plot3(y+0.5*log(3),y/v2,[1 1]*MZ12,'k--','LineWidth',2);
title('(one+one)-soliton solution')
axis([mx Mx mt-0.0*(Mt-mt) Mt+0.0*(Mt-mt)]);

subplot(2,2,2)
surf(X,T,log(Z12+1));
MZ12=max(max(log(Z12+1)));
shading interp
axis tight
colormap('default')
ylabel('t')
xlabel('x')
zlabel('u(x,t)')
view([0 90])
colorbar
hold on
plot3(y-0.25*log(3),y/v1,[1 1]*MZ12,'k--','LineWidth',2);
plot3(y+0.5*log(3),y/v2,[1 1]*MZ12,'k--','LineWidth',2);
title('two-soliton solution [log(z+1)]')
axis([mx Mx mt-0.0*(Mt-mt) Mt+0.0*(Mt-mt)]);
%axis tight

subplot(2,2,4)
surf(X,T,log(Z1+Z2+1));
MZ12=max(max(log(Z1+Z2+1)));
shading interp
colormap('default')
ylabel('t')
xlabel('x')
zlabel('u(x,t)')
view([0 90])
colorbar
hold on
plot3(y-0.25*log(3),y/v1,[1 1]*MZ12,'k--','LineWidth',2);
plot3(y+0.5*log(3),y/v2,[1 1]*MZ12,'k--','LineWidth',2);
title('(one+one)-soliton solution [log(z+1)]')
axis([mx Mx mt-0.0*(Mt-mt) Mt+0.0*(Mt-mt)]);