Split-Operator-Methode

Split-Operator-Methode

Die Split-Operator-Methode (SOP) ist ein numerisches Verfahren mit dem die zeitabhängige Schrödingergleichung gelöst werden kann. Bei der Methode wird der Hamiltonoperator $ {\hat {H}} $ in einen kinetischen Teil $ {\hat {T}} $ (Impulsteil) und in einen Potentialteil $ {\hat {V}} $ gespalten und einzeln angewendet. Dabei wird von der schnellen Fourier-Transformation (FFT) Gebrauch gemacht, um zwischen Impulsraum und Ortsraum zu unterscheiden.

Die Schrödingergleichung

Die Wellenfunktion $ \psi (x) $ auf einem äquidistanten Gitter dargestellt (Ortsraum)
Die Wellenfunktion $ \psi (k) $ auf einem äquidistanten Gitter dargestellt (Impulsraum)

Die zeitabhängige Schrödingergleichung ist definiert als

$ i\hbar {\frac {\partial \psi (x,t)}{\partial t}}={\hat {H}}(t)\psi (x,t), $

wobei $ \textstyle {\hat {H}}(t)=-{\frac {\hbar ^{2}}{2m}}\Delta +V(t) $ der Hamiltonoperator ist.

Die Wellenfunktion $ \psi (x,t) $ wird im Ortsraum auf einem äquidistanten Gitter dargestellt. Als Startwerte werden die Werte von $ \psi (x,t_{0}) $ zur Zeit $ t_{0} $ an den Gitterpunkten vorgegeben. Durch das Verfahren wird die Wellenfunktion zu einem späteren Zeitpunkt $ t=t_{0}+\Delta t $ berechnet.

Die Wirkung des Hamiltonoperators $ \textstyle {\hat {H}}(t)={\frac {{\hat {\mathbf {p} }}^{2}}{2m}}+{\hat {V}}(t) $ auf eine Wellenfunktion $ {\hat {H}}\psi ={\hat {T}}\psi +{\hat {V}}\psi $ wird mit der schnellen Fourier-Transformation berechnet. Dazu wird neben dem Gitter im Ortsraum auch ein Gitter im Impulsraum benötigt. Die Auflösung im Impulsraum $ \Delta k={\tfrac {2\pi }{L}} $ ist durch die Länge $ L $ des Gitters im Ortsraum festgelegt. Es gilt $ \Delta k\Delta x={\tfrac {2\pi }{N}} $, wobei $ N $ die Anzahl der Gitterpunkte ist.

Anwendung der diskreten Fourier-Transformation

Der Potentialoperator $ {\hat {V}} $ besitzt im Ortsraum eine diagonale Matrixdarstellung und wirkt daher lokal auf jeden Gitterpunkt $ x_{i} $:

$ \left({\hat {V}}(t)\psi \right)(x_{i})=V(x_{i},t)\cdot \psi (x_{i}). $

Genauso wird der kinetische Operator $ {\hat {T}} $ mit seiner diagonalen Darstellung im Impulsraum berechnet. Für jeden Gitterpunkt $ k_{i} $ gilt:

$ \left({\hat {T}}{\tilde {\psi }}\right)(k_{i}))={\frac {\hbar ^{2}}{2m}}{k_{i}}^{2}{\tilde {\psi }}(k_{i}). $

Dabei ist die diskrete Darstellung der Wellenfunktion $ {\tilde {\psi }}(k_{i}) $ im Impulsraum durch die diskrete Fourier-Transformation $ {\hat {Z}}^{\dagger } $ gegeben:

$ {\tilde {\psi }}(k_{i})=\langle k_{i}|\psi \rangle =\sum _{x_{j}}\langle k_{i}|x_{j}\rangle \langle x_{j}|\psi \rangle ={\frac {\Delta x}{\sqrt {2\pi }}}\sum _{x_{j}}e^{-ik_{i}x_{j}}\psi (x_{j}) $

In Vektorschreibweise lautet diese Gleichung

$ {\vec {\tilde {\psi }}}=c^{-1}{\hat {Z}}^{\dagger }{\vec {\psi }} $

mit

$ {\vec {\tilde {\psi }}}:=({\tilde {\psi }}(k_{0}),\dotsm ,{\tilde {\psi }}(k_{N-1}))^{T} $
$ {\vec {\psi }}:=(\psi (x_{0}),\dotsm ,\psi (x_{N-1}))^{T} $
$ Z_{ij}:=N^{-{\frac {1}{2}}}e^{+ik_{i}x_{j}} $
$ c:={\tfrac {\sqrt {2\pi N}}{L}}. $

Entsprechend erhält man für die Rücktransformation in den Ortsraum

$ \psi (x_{j})={\frac {\Delta k}{\sqrt {2\pi }}}\sum _{k_{i}}e^{+ik_{i}x_{j}}{\tilde {\psi }}(k_{i}) $

beziehungsweise

$ {\vec {\psi }}=c{\hat {Z}}{\vec {\tilde {\psi }}} $

mit den Gitterschrittweiten $ \Delta x={\tfrac {L}{N}} $ bzw. $ \Delta k={\tfrac {2\pi }{L}} $. Hierbei ist $ L $ die Länge des Gitters im Ortsraum und $ N $ die Zahl der Punkte im Orts- und Impulsraum. Die Konstante $ c $ wird nur benötigt, wenn die richtige Normierung der Funktion $ {\tilde {\psi }} $ gewünscht wird. Die Fourier-Transformation erhält die Norm der Vektoren $ {\vec {\tilde {\psi }}} $ und $ {\vec {\psi }} $.

Split-Operator-Methode

Die Berechnung der $ e $-Funktion eines Operators wird in der Diagonaldarstellung des Operators besonders einfach. Die Split-Operator-Methode verwendet eine Zerlegung des Hamiltonoperators in die Operatoren für kinetische Energie $ {\hat {T}} $ und für potentielle Energie $ {\hat {V}} $, welche im Impuls- bzw. Ortsraum Diagonalform annehmen.

Der durch die Nicht-Vertauschbarkeit von $ {\hat {T}} $ und $ {\hat {V}} $ entstehende Fehler kann durch die symmetrische Aufspaltung

$ e^{-{\frac {i}{\hbar }}{\hat {H}}\Delta t}\approx e^{-{\frac {i}{\hbar }}{\frac {\hat {T}}{2}}\Delta t}\cdot e^{-{\frac {i}{\hbar }}{\hat {V}}\Delta t}\cdot e^{-{\frac {i}{\hbar }}{\frac {\hat {T}}{2}}\Delta t} $

auf Terme der Größenordnung $ {\Delta t}^{3} $ reduziert werden: Mit $ {\hat {X}}:=-{\tfrac {i}{\hbar }}{\hat {T}}\Delta t $ und $ {\hat {Y}}:=-{\tfrac {i}{\hbar }}{\hat {V}}\Delta t $ erhält man für die rechte Seite

$ \exp \left({\frac {\hat {X}}{2}}\right)\exp \left({\hat {Y}}\right)\exp \left({\frac {\hat {X}}{2}}\right)=\exp \left({\frac {\hat {X}}{2}}+{\hat {Y}}+{\frac {\hat {X}}{2}}+\underbrace {{\frac {1}{2}}\left[{\frac {\hat {X}}{2}},{\hat {Y}}\right]+{\frac {1}{2}}\left[{\hat {Y}},{\frac {\hat {X}}{2}}\right]} _{0}+{\frac {1}{12}}\left[\left[{\frac {\hat {X}}{2}},{\hat {Y}}\right],{\hat {X}}+2{\hat {Y}}\right]+\dotsm \right). $

Der führende Fehlerterm ist somit proportional zu $ {\Delta t}^{3}\left[\left[{\hat {T}},{\hat {V}}\right],{\hat {T}}+2{\hat {V}}\right] $.

Diagonalform

Eine Koordinatentransformation $ {\hat {Z}}^{\dagger } $ vom Orts- in den Impulsraum ermöglicht eine einfache Berechnung von

$ e^{-{\frac {i}{\hbar }}{\frac {\hat {T}}{2}}\Delta t}\psi ={\hat {Z}}e^{-{\frac {i}{\hbar }}{\frac {\hat {T}}{2}}\Delta t}{\hat {Z}}^{\dagger }\psi . $

Mit der diagonalen Darstellung des Operators der kinetischen Energie

$ {\hat {T}}={\hat {Z}}^{\dagger }{\hat {T}}{\hat {Z}}={\begin{pmatrix}{\frac {\hbar ^{2}}{2m}}{k_{0}}^{2}&0&\cdots &0\\0&\ddots &0&0\\\vdots &\vdots &\ddots &\vdots \\0&0&\cdots &{\frac {\hbar ^{2}}{2m}}k_{N-1}^{2}\\\end{pmatrix}} $

erhält man

$ e^{-{\frac {i}{\hbar }}{\frac {\hat {T}}{2}}\Delta t}={\begin{pmatrix}\ddots &0&\cdots &0\\0&e^{-{\frac {i}{\hbar }}{\frac {\Delta t}{2}}{\frac {\hbar ^{2}}{2m}}k_{i}^{2}}&0&0\\\vdots &\vdots &\ddots &\vdots \\0&0&\cdots &\ddots \\\end{pmatrix}}. $

Die Koordinatentransformation erfolgt auf dem $ N $-Punkt-Gitter $ x_{0},\dotsm ,x_{N-1} $ mit Hilfe der diskreten Fourier-Transformation:

$ {\tilde {\psi }}(k_{i})=N^{-{\frac {1}{2}}}\sum _{j=0}^{N-1}\psi (x_{j})e^{-ik_{i}x_{j}} $     für   $ i=0,\dotsm ,N-1 $

oder $ {\vec {\tilde {\psi }}}={\hat {Z}}^{\dagger }{\vec {\psi }} $.

Numerischer Algorithmus

Durch Zusammenfassen der aufeinanderfolgenden Terme $ {\hat {Z}}e^{-{\frac {i}{\hbar }}{\frac {\hat {T}}{2}}\Delta t}{\hat {Z}}^{\dagger } $ zweier Zeitschritte lässt sich die Zahl der Fourier-Transformationen, d. h. der numerische Aufwand, reduzieren: $ {\hat {Z}}^{\dagger }{\hat {Z}}=1 $, und die beiden $ e $-Funktionen mit $ {\frac {\hat {T}}{2}} $ ergeben $ e^{-{\frac {i}{\hbar }}{\hat {T}}\Delta t} $.

Die Wellenfunktion nach $ n $ Zeitschritten erhält man also durch:

  • Fourier-Transformation von $ \psi _{0} $
  • Multiplikation mit den Diagonalelementen $ e^{-{\frac {i}{\hbar }}{\frac {\Delta t}{2}}{\frac {\hbar ^{2}}{2m}}k_{i}^{2}} $ (halber Zeitschritt)
  • Rücktransformation
  • Multiplikation mit den Diagonalelementen $ e^{-{\frac {i}{\hbar }}V_{i}\Delta t} $
  • Fourier-Transformation
  • Multiplikation mit den Diagonalelementen $ e^{-{\frac {i}{\hbar }}\Delta t{\frac {\hbar ^{2}}{2m}}k_{i}^{2}} $ (ganzer Zeitschritt)
  • usw., bis beim letzten Schritt noch einmal eine Multiplikation mit halben Zeitschritt wie in der zweiten Zeile notwendig wird.

Literatur

  • I. N. Bronstein, K. A. Semendjajew, G. Musiol, H. Muehlig: Taschenbuch der Mathematik. Deutsch Harri GmbH, 2008.
  • T. Fließbach: Quantenmechanik: Lehrbuch zur Theoretischen Physik III. 5. Auflage, Spektrum Akademischer Verlag, 2008, ISBN 978-3-8274-2020-6.
  • Herbert Sager: Fourier-Transformation. vdf Hochschulverlag, Zürich 2012, ISBN 978-3-7281-3393-9.
  • A. Askar, A. S. Cakmak: Explicit integration method for the time‐dependent Schrodinger equation for collision problems. In: Journal of Chemical Physics. Band 68, Nr. 6, 1978, S. 2794–2798, doi:10.1063/1.436072.
  • J. B. Delos: Theory of Electronic Transitions in Slow Atomic Collisions. In: Physical Review. Band 176, Nr. 1, 1968, S. 141–150, doi:10.1103/PhysRev.176.141.
  • Juha Javanainen, Janne Ruostekoski: Symbolic calculation in development of algorithms: split-step methods for the Gross–Pitaevskii equation. In: Journal of Physics A. Band 39, 2006, S. L179–L184, doi:10.1088/0305-4470/39/12/L0.
  • Michael Hintenender: Propagation von Wellenpaketen. In: MPQ-Berichte. MPQ163. Garching 1992 (online).