I'm trying to use the Kutta-Merson to get the same results as in the book Solitons, Nonlinear Evolution Equations and Inverse Scattering - M. J. Ablowitz - pg 140
The author propose using the Kutta-Merson in the following scheme:
$$ i\frac{du_j}{dt} + (u_{j-1} - 2u_j + u_{j+1})/h^2 + |u_j|^2(u_{j-1} + u_{j+1}) = 0 $$
given condition of periodicity and the following initial condition
$$ u(x,0) = 0.5 + 0.05 \cos(\mu x) + 10^{-5} i \sin(\mu x) $$ where $L = 2\pi\sqrt{2}$ is the length of the interval and $\mu = 2\pi/L$
I'm not very familiar with this method for PDE, but here is what i try to do, using scilab:
I define the laplacian and the nonlinear using scilab functions, and put the scheme in a form of Kutta's method
$$\frac{du_j}{dt} = F(u_j)$$
I use the condition of periodicity to define F. Here is the Laplacian:
function l =L(u)
A = zeros(M+1,M+1)
// Condição periódica imposta ao laplaciano
A(1,1) = -2; A(1,2) = 1; A(1,M+1) = 1
for j=2:M
A(j,j-1) = 1; A(j,j) = -2; A(j,j+1) = 1
end
A(M+1,1) = 1; A(M+1,M) = 1; A(M+1,M+1) = -2
A = h^(-2)*A
l = A*u
endfunction
And here is the nonlinear therm
function n = N(u)
n(1) = abs(u(1))^2 * (u(M+1) + u(2))
for j=2:M
n(j) = abs(u(j))^2 *(u(j-1) + u(j+1))
end
n(M+1) = abs(u(M+1))^2 * (u(M) + u(1))
endfunction
Then F becomes:
function f = F(u)
f = %i*(L(u) + N(u))
endfunction
Finally I use this to calculate the solution
for n=1:Nt
k1 = k*F(U(:,n))
k2 = k*F(U(:,n) + 1/3*k1)
k3 = k*F(U(:,n) + 1/6*k1 + 1/6*k2)
k4 = k*F(U(:,n) + 1/8*k1 + 3/8*k3)
k5 = k*F(U(:,n) + 1/2*k1 - 3/2*k3 + 2*k4)
U(:,n+1) = U(:,n) + 1/6*k1 + 2/3*k4 + 1/6*k5
end
I'm using 5000 points at time grid and 20 for x-grid. I really don't know nothing about so if I'm doing some weird thing let me know and recommend some material about it please