Finite difference and Runge-Kutta for PDEs

[EmpaDoc]

I made a small program to simulate the time development of a 1D wavepacket obeying the Schrodinger equation, mostly in order to learn a new programming language - so in order to not have to invoke big numerical methods packages, I opted for the simplest solution: The standard three-point finite-difference discretization for the Laplacian (kinetic term), and the 4-order Runge-Kutta method to propagate the wavepacket forward.

I notice that sooner or later, my wavepacket will start to develop oscillations on the smallest wavelengths, followed by a blowup. If I include a step potential to reflect against, this happens once the wavepacket encounters the potential. Otherwise it seems to happen when the wavepacket scatters against itself across the periodic boundary (as it expands).

Changing my three-point finite-difference to a five-point stencil only makes the problem worse.

I know that there are much better methods for PDE solving (e.g. Crank-Nicolson). But that's not really the point here. Now I am curious: What is the cause of this behavior? And-can I tweak the FD+RK so I can get it to work?

 

[AlephZero]

The first thing I would try is reducing your time step (e.g. halve it). Possibly (by happenstance) your time step is close to the limit for stability, so it appears to be working OK until the instability has grown from 10^-6 or (10^-16) to be bigger than the true solution.

Using a solution scheme that has a higher order of accuracy than the order of the differential equation almost always results in conditional stability. That's why a 4th order method like RK is rarely "better" than a well-chosen lower order method like Crank Nicholson - there are very few 4th order differential equations in physics!

 

[EmpaDoc]

Forget it! I just hadn't debugged the program properly...
Sorry. :(
At least I know now that FD+RK4 works fine with the Schrodinger equation.

 

[Sumert]

Can you explain what was the problem in your code ? I am also trying such a way but norm of wave function isn't constant and differs from analytical solution (Im trying to propagate free gaussian wave packet).

 

[LudusRex]

Thank you for sharing your experience with using finite difference and Runge-Kutta methods for solving PDEs. It is commendable that you took the initiative to learn a new programming language and explore these numerical methods on your own.

The behavior you are observing, where the wavepacket develops oscillations and eventually blows up, is a common issue when using finite difference and Runge-Kutta methods for solving PDEs. This is known as numerical instability and can be caused by a variety of factors.

One possible cause is the choice of grid spacing in your finite difference discretization. If the grid spacing is too large, it can lead to inaccurate solutions and numerical instability. You can try decreasing the grid spacing to see if it improves the stability of your simulation.

Another factor to consider is the time step used in the Runge-Kutta method. If the time step is too large, it can also lead to numerical instability. You can experiment with different time steps to see if it helps improve the stability of your simulation.

Additionally, the boundary conditions you are using can also play a role in the stability of your simulation. If the boundary conditions are not properly handled, it can lead to numerical instability. You can review and adjust your boundary conditions to see if it helps improve the stability of your simulation.

It is important to note that while finite difference and Runge-Kutta methods are simple to implement, they may not always be the most efficient or accurate methods for solving PDEs. As you mentioned, there are other methods such as Crank-Nicolson that may be more suitable for your problem. It would be beneficial to explore and compare different methods to find the most appropriate one for your specific PDE.

In summary, the cause of the behavior you are observing is likely due to numerical instability, and you can try adjusting the grid spacing, time step, and boundary conditions to improve the stability of your simulation. However, it is also important to consider other methods that may better suit your problem. Overall, your curiosity and willingness to experiment with different methods is commendable and will ultimately help you gain a deeper understanding of PDEs and numerical methods.