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).