I'm using 4th order RK to solve the schroedinger equation in atomic units. Say I want to simulate 400fs in intervals of h=10fs, then in atomic units this is h=413a.u and 400fs=16500a.u. 4RK involves repeated multiplication by the timestep, so with h=413 everything blows up. Could someone explain what I'm missing?

Edit: to be clear, I'm solving

\begin{align} \frac{d}{dt} \psi = -i H \psi \end{align} Where $\psi$ is a vector and $H$ a hermitian matrix.

  • $\begingroup$ Runge Kutta methods are explicit schemes in time. There is a relationship between the time step and the spatial discretisation to fulfill in order to get rid of instabilities. I supose that the Hamiltonian operator has been discretised. $\endgroup$ – HBR Feb 1 '18 at 16:48
  • $\begingroup$ Could you elaborate please? The system is a lattice with lattice constant 8 (so space is discrete) but the Hamiltonian is in a basis of Fock states. $\endgroup$ – Insert_Username Feb 1 '18 at 17:11
  • $\begingroup$ Currently I am with mobile phone... and I am lazy to write here. Search on google about heat diffusion equation temporal scheme. In wikipedia you will find the restriction I was talking about. Regarding to the basis of Fock states... it is too specific to know what is that... $\endgroup$ – HBR Feb 1 '18 at 17:26

This phenomena of an explicit Runge-Kutta method is known as stability. For a given equation, the eigenvalues of your system (here, the Hamiltonian) determine a maximum time step for which the numerical integration scheme will be stable (i.e. not blow up). Here's a good resource which mentions the stability region for RK4. As these resources show, there are special methods such as implicit methods which can be used in order to take larger time steps.

But even explicit Runge-Kutta methods can be enhanced. For example, using adaptive Runge-Kutta methods with step-stabilization (PI-adaptive step controllers) will allow the method to only take smaller time steps when they are needed for stability, and this then allows the steps to be (usually) larger on average. Standard ODE solver software like those in MATLAB, Fortran, or Julia make use of this. So I would recommend using whatever ODE solvers you have available in your language, since not only will they have (potentially PI-controlled) adaptive methods, but also integrators for stiff systems. These will be much more efficient than RK4 which is really a method just for teaching and not for production use (other explicit RK methods have larger stability regions and/or less error for the same amount of work).


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.