I solved Schrödinger equation for a following tow-level time-dependent Hamiltonian numerically in two ways:

import numpy as np

def H(t):
    return np.array([[t,0.5],[0.5,-t]],dtype="complex64")

First, the state is treated as wave vector:

from scipy.integrate import solve_ivp
def schrodinger(t, X):
    dXdt = -1j * H(t).dot(X)
    return dXdt

t_list = np.arange(-T, T, 2*T / 1000.0)
init_state= np.array([1.,0.],dtype="complex64")
solution = solve_ivp(fun=lambda t, X: schrodinger(t, X), t_span=[-T,T], y0=init_state, t_eval=t_list, method="RK45", vectorized=True)

Second, the state is treated as (vectorized) density matrix:

def schrodinger_rho(t, X):
    unit = np.eye(2, 2, dtype="complex64")
    dXdt = -1j * (np.kron(unit, H(t)) - np.kron(H(t).T, unit)).dot(X)
    return dXdt
t_list = np.arange(-T, T, 2*T / 1000.0)
init_state= np.array([1.,0.,0.,0.],dtype="complex64")
solution_rho = solve_ivp(fun=lambda t, X: schrodinger_rho(t, X), t_span=[-T,T], y0=init_state, t_eval=t_list, method="RK45", vectorized=True)

The numerical results are as follows:

wave vector density matrix

Here is a question : why cannot the state preserve its unitarity, say $P_0+P_1=1$, when the state is treated as wave vector? The larger the value of T, the stronger this tendency of violation becomes. Moreover, the numerical results become worse when I use method="BDF".

I want to treat the state as wave vector because of the numerical cost. Is there any way to improve this phenomenon?

  • 1
    $\begingroup$ Explicit Runge Kutta methods rarely conserves integrals of motion. Did you try with an implicit Gauss Runge Kutta method? $\endgroup$ Mar 3, 2021 at 22:50
  • $\begingroup$ @G. Fougeron But isn't it true that for a smaller integration time step the integrals of motion are conserved better, even if the method is not conservative? So for a small enough time step you'll get all conservation laws that the underlying PDE has? $\endgroup$ Mar 4, 2021 at 6:28
  • $\begingroup$ This is true for the short time behavior, but not necessarily for long-time behavior. For more precise info on the subject, I recommend the book by Hairer : Geometric Numerical Integration. Structure-Preserving Algorithms for Ordinary Differential Equations. If you can't find it, then maybe chackout his list of preprints unige.ch/~hairer/preprints.html $\endgroup$ Mar 4, 2021 at 10:21
  • $\begingroup$ @G.Fougeron I did not try it because there is no choice of implicit Runge Kutta method in complex domain in scipy. Instead, I found odeintw where LSODA can be applied in complex domain. This method works. $\endgroup$
    – wayna
    Mar 4, 2021 at 13:38
  • 1
    $\begingroup$ There you go. Then it makes perfect sense. If you check out the resources I linked, you will find proofs / explanations for this 😊 $\endgroup$ Mar 5, 2021 at 8:01


Your Answer

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