2
$\begingroup$

I am calculating numerical solutions of time-dependent Schroedinger equation

$\frac{d\Psi}{dt} = - i H \Psi$

where $\Psi$ is an $N$-element complex vector and $H$ is an $N \times N$ complex matrix, which is ``almost normal''. That is, $H = H_0 + i D$, where $H_0$ is normal (and often Hermitian), $D$ is Hermitian but does not commute with $H_0$, and $||D|| \ll ||H_0||$. The dimension $N$ is betwen 30 and 200.

So far, I have been using a 4-th order Runge-Kutta specialised for linear ODEs (taken from Zingg and Chisholm). However, the accuracy is not as good as I'd like, probably due to strong oscillations arising from the $H_0$ part. Is there a better method out there? I'd prefer one which would not require the calculation of $e^{-iH\tau}$ ($\tau$ - time-step of the ODE solver).

EDIT: more data about the problem:

The eigenvalues of $H_0$ are literally random, as $H_0 = h + diag(dE_i)$, where $dE_i ~ N(0,\sigma)$ or $dE_i \sim L(\alpha,\sigma)$. Matrix $h$ has eigenvalues of the order of unity, and $L(\alpha,\sigma)$ is the symmetric Levy alpha-stable distribution (each diagonal value of $H_0$ is perturbed by an independent distr.) with scale parameter $\alpha \in (0,2)$ and strength $\sigma$. $D$ has eigenvalues of the order of 0.1. I choose the time step to be of the order of $0.01/E_\text{max}$, where $E_\text{max}$ is the maximum absolute value of real or imaginary part of any eigenvalue of $H$.

So in short, $H$ can be nasty :(

Improve this question
$\endgroup$
13
  • $\begingroup$ Are both H_0 and D explicitly time dependent? If H_0 is time independent, then you could diagonalize it once and for all at the start of your calculation and use the interaction picture. $\endgroup$
    – Jonathan
    Jul 21, 2011 at 14:15
  • $\begingroup$ Additionally, do you know anything about the eigenvalues of D? The sign will determine whether the solution grows or decays, which might be relevant to cooking up a better numerical method. $\endgroup$
    – Jonathan
    Jul 21, 2011 at 14:27
  • $\begingroup$ @Jonathan $H_0$ is time-independent. But then my perturbation becomes time-dependent, won't it simply move the oscillations to another place in the calculation? $D$ has only negative eigenvalues, so it should lead to a decay. Funnily, for some forms of $H_0$ the Runge-Kutta algorithm I use now blows up anyway. $\endgroup$
    – Katastrofa
    Jul 21, 2011 at 16:56
  • $\begingroup$ @Jonathan Do you think that using the split-operator approach would help? I mean the approach in pci.uni-heidelberg.de/tc/usr/andreasm/academic/handouthtml/… $\endgroup$
    – Katastrofa
    Jul 21, 2011 at 16:58
  • $\begingroup$ What are the smallest and the largest eigenvalues of $H$ and what are the norm and the smoothness of $D$ and the time interval? $\endgroup$
    – fedja
    Jul 21, 2011 at 18:47

1 Answer 1

Reset to default
4
$\begingroup$

This should really be a comment, but it is too long.

Since Euler is out of question on such a long interval, my next suggestion would be to try to run the 3rd order Runge-Kutta about which I am pretty certain that it is not screwed up anywhere. The recursion step is $$ \begin{aligned} x(t+\tau)&= \cr &x(t) \cr +&\frac \tau 6[H(t)+4H(t+\frac\tau 2)+H(t+\tau)]x(t) \cr +&\frac{\tau^2}6[H(t)^2+2H(t+\tau)H(t+\frac\tau 2)]x(t) \cr +&\frac{\tau^3}{6}H(t)^3x(t) \end{aligned} $$ (of course, you should never multiply matrices, only a matrix by a vector and you should pay attention to the order of matrix multiplication in the third term; I leave the obvious optimizations to you).

Try to compare the results with what you get from your 4th order Runge-Kutta. The maximal step for this 3rd order method would be $10^{-4}$ if you want to get anything meaningful in your setting (if you can afford $10^{-5}$, it would be much better). Let us know what you see. If this method exhibits slower norm growth than the 4th order one, something is screwed up in the programming. If it blows up even faster, it may, indeed, be a precision problem.

Improve this answer
$\endgroup$
5
  • $\begingroup$ Sorry, my $H$ is yours $-iH$, so put all $i$'s accordingly. $\endgroup$
    – fedja
    Jul 22, 2011 at 15:20
  • $\begingroup$ Thanks. My $H$ is time-independent, so I can ignore the order of multiplication. $\endgroup$
    – Katastrofa
    Jul 22, 2011 at 16:31
  • 1
    $\begingroup$ In that case what prevents you from forgetting all that Runge-Kutta stuff altogether and just using the standard Taylor expansion of the exponent with as many terms as you wish (6-7 should be more than enough)? I mean use $x(t+s)=x(t)-isHx(t)-\frac 12s^2H^2x(t)+\frac{i}6s^3H^3x(t)+\dots$. It is not overly computationally intensive (just a 7 matrix by vector multiplications at each step) and if you choose $s\|H\|<0.1$, you should be completely fine. Runge-Kutta is mainly there to accomodate for time dependency of the matrix. $\endgroup$
    – fedja
    Jul 22, 2011 at 18:20
  • $\begingroup$ @fedja Indeed, why didn't I think about it earlier? $\endgroup$
    – Katastrofa
    Jul 22, 2011 at 19:28
  • $\begingroup$ @fedja Would you mind editing the answer to include your last comment? it solved my problem. $\endgroup$
    – Katastrofa
    Jul 27, 2011 at 16:15

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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