When plotting the numerical solution at step 4, initial conditions of the equation was taken from analytical solution (by using un = u.copy() at the beginning of the FD loop, array u is the analytical solution of function at t=0 in our code). I think we need to use initial conditions below in the picture.

When plotting the numerical solution at step 4, initial conditions of the equation was taken from analytical solution (by using un = u.copy() at the beginning of the FD loop, array u is the analytical solution of function at t=0 in our code). I think we need to use initial conditions below in the picture.