+ All Categories
Home > Documents > Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der...

Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der...

Date post: 29-Oct-2019
Category:
Upload: others
View: 8 times
Download: 0 times
Share this document with a friend
109
Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls Universität Tübingen
Transcript
Page 1: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Theorie und Numerik derBorn-Oppenheimer-Approximation

zweiter Ordnung

Daniel Weber

Fachbereich MathematikEberhard Karls Universität Tübingen

Page 2: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls
Page 3: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Theorie und Numerik derBorn-Oppenheimer-Approximation

zweiter Ordnung

Daniel Weber

Oktober 2018

Masterarbeit

Fachbereich MathematikEberhard Karls Universität Tübingen

Betreut von

Prof. Dr. S. Teufel

Page 4: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls
Page 5: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Abstract

Im Verlauf dieser Arbeit wollen wir die Born-Oppenheimer-Approximation qualitativmit Hilfe von asymptotischer Analysis studieren und den quantitativen Mehrwert vonKorrekturen höherer Ordnung numerisch untersuchen. Wir erarbeiten uns zunächstdie mathematischen Grundlagen der zeit- und raumadiabatischen Störungstheorie, da-bei formulieren wir insbesondere das zeit- sowie das raumadiabatische Theorem. An-schließend stellen wir die Born-Oppenheimer-Approximation vor und präsentieren einsystematisches Schema für Korrekturen höherer Ordnungen. Dabei leiten wir die zeit-abhängige Born-Oppenheimer-Approximation zweiter Ordnung in expliziter Form her.In einem abschließenden Beispielmodell verdeutlichen wir den Mehrwert der korrektenKorrekturterme anhand von numerischen Simulationen.

Page 6: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls
Page 7: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Inhaltsverzeichnis

1 Einführung 1

2 Adiabatische Störungstheorie für die Schrödinger Gleichung 72.1 Das zeitadiabatische Theorem der Quantenmechanik . . . . . . . . . . . 7

2.1.1 Mathematische Grundlagen . . . . . . . . . . . . . . . . . . . . . 82.1.2 Das zeitadiabatische Theorem . . . . . . . . . . . . . . . . . . . . 10

2.2 Raumadiabatische Störungstheorie . . . . . . . . . . . . . . . . . . . . . 172.2.1 Rahmenbedingungen . . . . . . . . . . . . . . . . . . . . . . . . . 172.2.2 Das raumadiabatische Theorem . . . . . . . . . . . . . . . . . . . 212.2.3 Effektive Dynamik . . . . . . . . . . . . . . . . . . . . . . . . . . 22

3 Die zeitabhängige Born-Oppenheimer-Approximation 253.1 Die konventionelle Born-Oppenheimer-Approximation . . . . . . . . . . 25

3.1.1 Anwendung auf molekulare Dynamik . . . . . . . . . . . . . . . . 253.1.2 Der effektive Hamiltonoperator . . . . . . . . . . . . . . . . . . . 28

3.2 Korrekturen für die Approximation höherer Ordnungen . . . . . . . . . 29

4 Numerische Simulation 414.1 Umsetzung der Anwendung . . . . . . . . . . . . . . . . . . . . . . . . . 41

4.1.1 Rahmenbedingungen . . . . . . . . . . . . . . . . . . . . . . . . . 414.1.2 Numerische Methoden . . . . . . . . . . . . . . . . . . . . . . . . 494.1.3 Implementierung . . . . . . . . . . . . . . . . . . . . . . . . . . . 504.1.4 Konvergenztests . . . . . . . . . . . . . . . . . . . . . . . . . . . 51

4.2 Durchführung der Approximation . . . . . . . . . . . . . . . . . . . . . . 554.2.1 Vergleich der Fehler . . . . . . . . . . . . . . . . . . . . . . . . . 554.2.2 Vergleich der Laufzeiten . . . . . . . . . . . . . . . . . . . . . . . 584.2.3 Verhalten für längere Zeiten . . . . . . . . . . . . . . . . . . . . . 59

5 Fazit und kurzer Ausblick 63

A Appendix 65A.1 Hermitesche Polynome . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65A.2 Quelltext . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67

A.2.1 run.py . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67A.2.2 parameter.py . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71A.2.3 function_data_base.py . . . . . . . . . . . . . . . . . . . . . . . 73

v

Page 8: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Symbole und Zeichen 93

Abbildungsverzeichnis 95

Tabellenverzeichnis 97

Literaturverzeichnis 99

vi

Page 9: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

1 Einführung

Die Born-Oppenheimer-Approximation findet in vielen Bereichen der Naturwissenschaf-ten wie der Physik und Chemie Anwendung und ist eine der wichtigsten und grundle-gendsten Approximationen in der Quantenmechanik. Sie ist benannt nach Max Bornund J. Robert Oppenheimer und bezeichnet eine Näherung zur Vereinfachung der Schrö-dingergleichung von aus mehreren Teilchen bestehenden Systemen. Dieser Näherungliegt die adiabatische Störungstheorie zu Grunde. Im Allgemeinen versteht man unter„Störungstheorie“ physikalische Verfahren, bei denen ein kompliziertes Problem zu-nächst idealisiert wird, das heißt solange durch Ignorieren erwartungsgemäß kleinerEinflüsse vereinfacht wird, bis es auf ein Problem mit bekannter Lösung reduziert ist.Anschließend werden die zuvor ignorierten Einflüsse dem System als kleine Störungenwieder hinzugefügt und eine approximative Lösung berechnet. Für das Verständnis desdynamischen Verhaltens von komplexen Systemen spielt eine Trennung von Zeitskaleneine wichtige Rolle und hilft eine Lösung für die Dynamik zu finden. Liegt solch eineTrennung vor, ist es meist möglich, einfache dynamische Gesetzmäßigkeiten für gewisselangsame Freiheitsgrade aus der zugrunde liegenden Struktur der schnellen herzuleiten.In dieser Arbeit betrachten wir quantenmechanische Systeme, die besagtes Verhaltenaufweisen.

Die Dynamik eines quantenmechanischen Systems, das heißt dessen zeitliche Ent-wicklung, wird im allgemeinen durch die Schrödingergleichung

i~ ∂∂t

Ψ = − ~2

2m∆xΨ + VΨ =: HΨ (1.1)

beschrieben. Startwert sei die aus einem Hilbertraum H stammende WellenfunktionΨ0 := Ψ

∣∣t=0 ∈ H. Mit der Naturkonstante ~ bezeichnen wir das Plancksche Wirkungs-

quantum, m sei die Masse der Teilchen, ∆x der Laplace-Operator und V das Potential.Der Hamiltonoperator H ist die Summe der kinetischen Energien aller Teilchen plus derpotentiellen Energie der zu dem System gehörenden Teilchen und beschreibt damit dieGesamtenergie des Systems. Trotz der prinzipiell einfachen mathematischen Strukturder Schrödingergleichung ist eine analytische Lösung häufig sehr schwierig oder nichtmöglich. Abhängig von der Dimension des zugrundeliegenden Konfigurationsraums istselbst eine numerische Lösung mit der aktuellen Rechenleistung heutiger Computerundenkbar. Folgerichtig ist es von großem Interesse, Konzepte zu finden, die eine best-mögliche approximative Lösung von (1.1) liefern. Eines dieser Konzepte ist es, mittelsder angesprochenen Trennung der Skalen die Anzahl der relevanten Freiheitsgrade zureduzieren. Die effektive Gleichung

i~ ∂∂tψ(t, x) = Heffψ(t, x), ψ0 := ψ

∣∣t=0 ∈ Heff , (1.2)

Page 10: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 1. Einführung

lebt auf einem Konfigurationsraum mit weniger Dimensionen und beschreibt mit einemeinfacheren effektiven HamiltonoperatorHε

eff dennoch die für uns interessante Dynamik.Der effektive Hamiltonoperator ist nicht nur einfacher in dem Sinne, dass er auf einemRaum mit weniger Dimensionen lebt, sondern liefert zudem eine bessere physikalischeIntuition. Durch die Reduzierung der Dimension ist (1.2) einfacher zu lösen.

Die adiabatische Störungstheorie eröffnet uns einen Weg, eine effektive Gleichung inquantenmechanischen Systemen zu finden. Sie nutzt aus, dass schnelle Freiheitsgradevon langsamen Freiheitsgraden abhängen, welche sich wiederum selbständig mit derZeit entwickeln. Dies wird als adiabatische Entkopplung bezeichnet.

Abhängig von der Struktur des zugehörigen Hamiltonoperators unterscheiden wirzwischen zwei Arten von adiabatischen Systemen. In der zeitadiabatischen Störungs-theorie ist der Hamiltonoperator zeitabhängig. Teil dieser Theorie ist das zeitadiabati-sche Theorem, welches häufig auch nur als das adiabatische Theorem bezeichnet wird.Der Begriff „adiabatisch“ stammt von dem griechischen „adiabatos“, was übersetzt„nicht hindurchtretend“ heißt. Er wird traditionell auch in der Thermodynamik ver-wendet, wo er jedoch eine andere Bedeutung hat. In dem Zusammenhang wie wir denBegriff verwenden, ist die Eigenschaft eines Prozesses gemeint, sich unter allmählich än-dernden Bedingungen anzupassen. Beginnt das System in einem Eigenzustand des ur-sprünglichen Hamiltonoperators, endet es im entsprechenden Eigenzustand des letztenHamiltonoperators. Dies ist sogleich ein Teil der Aussage des adiabatischen Theorems.

In der raumadiabatischen Theorie hängt der Hamiltonoperator nicht explizit vonder Zeit ab. Sie bildet die Basis für die Born-Oppenheimer-Approximation, auf diewir im Verlauf dieser Arbeit hinarbeiten, indem das klassische adiabatische Theoremauf ein raumadiabatisches Theorem verallgemeinert wird. Ortsvektoren übernehmendabei die Rolle der Zeit aus der zeitadiabatischen Theorie. Der Konfigurationsraumdes betrachteten Quantensystems wird in einen Unterraum der schnellen und in einenUnterraum der langsamen Freiheitsgrade geteilt. Ein klassisches Beispiel für solch einSystem sind Moleküle mit schweren Kernen und leichten Elektronen, also bestehendaus zwei Arten von Teilchen mit stark unterschiedlicher Masse. Aufgrund dieses Mas-senunterschieds bewegen sich die Elektronen um ein vielfaches schneller als die Kerne.Eine Unterteilung in schnelle und langsame Freiheitsgrade ist möglich. Die effektiveDynamik für die langsamen Freiheitsgrade, genauer gesagt für die Kerne, ist als Born-Oppenheimer-Näherung bekannt und für das Verständnis der Bewegung von Molekülenvon außerordentlicher Bedeutung. Vereinfacht gesagt entwickeln sich die Kerne in derBorn-Oppenheimer-Näherung in einem effektiven Potential, das von einem Energieni-veau der Elektronen, also einem Eigenwert des effektiven Hamiltonoperators, erzeugtwird. Indessen wird der Zustand der Elektronen sofort zu einem Eigenzustand, der dermomentanen Konfiguration der Kerne entspricht.

Die Approximation der Schrödingergleichung des Gesamtsystems (1.1) durch dieeffektive Gleichung (1.2) verursacht einen Fehler. Diesen gilt es zu minimieren. DasStudium der adiabatischen Störungstheorie hat ergeben, dass der Fehler kleiner als dasProdukt ε · C ist, wobei 0 < ε � 1 adiabatischer Parameter genannt wird und dieKonstante C von diesem unabhängig ist. Es sei angemerkt, dass die in dieser Arbeitangegebenen Konstanten aus Notationsgründen generisch sind und von Gleichung zuGleichung unterschiedlich sein können. Der adiabatische Parameter steuert die Tren-nung der Skalen. Je kleiner ε ist, desto besser sind langsame und schnelle Freiheitsgrade

2

Page 11: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

voneinander getrennt und desto kleiner ist der durch die Approximation bedingte Feh-ler. Da dieser Parameter in der Praxis durch das physikalische Problem festgelegt istund nur begrenzt verkleinert werden kann, ist es von Interesse, den Fehler auf andereArt undWeise zu optimieren. Dies führt zur Born-Oppenheimer-Approximation höhererOrdnungen. Ziel dieser Arbeit ist es, ein systematisches Schema für solche Korrektu-ren auszuarbeiten und die Born-Oppenheimer-Approximation zweiter Ordnung konkretherzuleiten. In unserem Fall wird damit die obere Schranke des Fehlers von ε · C aufε2 · C verkleinert.

Ein sich in dem Konfigurationsraum Rd+q entwickelndes Quantensystem könnenwir mit dem adiabatischen Parameter durch die zeitabhängige Schrödingergleichungmit einem zeitunabhängigen Potential beschreiben. In Abschnitt 2.2 werden wir (1.1)für x ∈ Rd und y ∈ Rq zu

i~ε ∂∂tΨ ε(t, x, y) = HεΨ ε(t, x, y)

mit dem selbstadjungierten Hamiltonoperator

Hε = − ~2

2mε2∆x −~2

2m∆y + V (x, y)

umformen. Das Trennen des Systems in langsame und schnelle Freiheitsgrade entsprichtder Aufteilung des Hilbertraums H = L2(Rd) ⊗Hf . Das bedeutet, dass die Variable xdie langsamen und y die schnellen Freiheitsgrade repräsentiert.

Auf die in Kapitel 3 thematisierte molekulare Dynamik übertragen sind die Orts-vektoren x = {x1, . . . , xK} und y = {y1, . . . , yl} die Positionen der K Kerne und l Elek-tronen, das heißt das System entwickelt sich in Rd+q = R3(K+l). Fordern wir, dass alleKerne die gleiche Masse mn haben und wählen für die Masse der Elektronen me = 1 so-wie das Plancksche Wirkungsquantum ~ = 1, können wir den adiabatischen Parameterdefinieren durch

ε =√memn

.

Auf diese Weise ist im molekularen Fall der Hamiltonoperator Hε von der Form

Hε = −12ε

2∆x −12∆y + V (x, y).

Den hinteren Teil des Terms definieren wir durch

Hf(x) := −12∆y + V (x, y).

Dieser selbstadjungierte Hamiltonoperator hängt von x ∈ Rd ab und wirkt auf den Hil-bertraum der schnellen Freiheitsgrade, den wir mit Hf bezeichnen. Das Spektrum vonHf muss dafür der sogenannten Lückenbedingung genügen. Es muss einen isolierten Teildes Spektrums geben, der durch eine Lücke vom restlichen Teil des Spektrums getrenntist. Die Lückenbedingung ist das Fundament der adiabatischen Störungstheorie in dieserArbeit und darüber hinaus von wesentlicher Bedeutung für viele Ergebnisse in diesemGebiet. Im molekularen Fall besteht der isolierte Teil in der Regel nur aus einer einzel-nen Energieoberfläche Ej(x). Zustände, deren Elektronen sich im zu Ej(x) korrespon-dierenden Eigenzustand befinden, sind von der Form Ψ ε(t, x, y) = ψε(t, x)φj(x, y) und

3

Page 12: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 1. Einführung

zusammengesetzt aus einer nur von x abhängigen Funktion ψε(t, x) und der Eigenfunk-tion φj(x, y) zum Eigenwert Ej(x). Lösungsverfahren zur Ermittlung der Eigenwertewerden wir nicht weiter thematisieren, da diese sehr komplex sein können. Stattdessensehen wir die entsprechenden Eigenwertprobleme als bereits gelöst und die Eigenwerteals gegeben an. In Abschnitt 2.2.3 werden wir sehen, dass eine effektive Schrödinger-gleichung, die die Dynamik der langsamen Freiheitsgrade, also im molekularen Fall dieder Kerne, beschreibt, gegeben ist durch

i~ε ∂∂tψε(t, x) = Hε

effψε(t, x), ψε(t0) = ψε0 ∈ L2(Rd).

Den effektiven Hamiltonoperator bestimmen wir am Ende von Abschnitt 3.1 und er-halten

Hεeff = −1

2ε2∆x + Ej(x) + 1

2ε2VBH(x). (1.3)

Der Ausdruck VBH bezeichnet das sogenannte Born-Huang-Potential, welches zusammenmit dem Eigenwert der schnellen Freiheitsgrade Ej(x) die effektive potentielle Energiebildet. Ist die Lückenbedingung erfüllt, haben die schnellen Freiheitsgrade nur überdieses effektive Potential Einfluss auf die Dynamik der langsamen Freiheitsgrade.

Der Fehler, der bei der Verwendung von (1.3) entsteht, ist von der Ordnung O(ε).Aus diesem Grund wird (1.3) auch konventioneller Born-Oppenheimer-Hamiltonopera-tor erster Ordnung genannt. Mithilfe von Korrekturtermen kann er verbessert werden,sodass der Fehler der Born-Oppenheimer-Approximation n-ter Ordnung von der Ord-nung O(εn) ist. Die allgemeine Konstruktion höherer Ordnungen ist sehr technisch,weshalb in Kapitel 3 dieser Arbeit nur das Grundkonzept skizziert wird. Den Born-Oppenheimer-Hamiltonoperator zweiter Ordnung werden wir hingegen explizit herlei-ten. Er ist in etwa von der Form

HεBO = −1

2ε2∆x + Ej(x) + 1

2ε2VBH(x)− ε2M(x, p),

wobei wir in dieser Darstellung den später auftauchenden Berry-Zusammenhangskoef-fizienten vernachlässigt haben. Bei dem letzten Term handelt es sich um die Korrektur,auf die wir später genauer eingehen werden.

Die theoretische Verbesserung des Fehlers durch den Schritt von der Born-Oppen-heimer-Approximation erster Ordnung zur zweiten Ordnung werden wir anhand einesexpliziten Beispielmodells überprüfen. Dazu studieren wir die Dynamik einer zweidi-mensionalen Schrödingergleichung, deren effektive Dynamik nur noch durch eine eindi-mensionale Gleichung beschrieben wird. Es wird deutlich werden, dass die Berechnungder effektiven Dynamik um ein vielfaches schneller vonstatten geht als die direkte Be-rechnung im Gesamtsystem. Als überzeugendes Ergebnis wird sich hervorheben, dassder Fehler der Approximation zweiter Ordnung selbst oder gerade über lange Zeitenhinweg signifikant kleiner ist als der der ersten Ordnung. Bemerkenswert in unseremModell ist zudem, dass sich bei beiden Ordnungen der Fehler ungefähr so verhält wiedie Schranke ε ·C beziehungsweise ε2 ·C und sich ihr nicht nur asymptotisch annähert.

Gliederung der Arbeit

Als Grundlage dieser Arbeit dienten hauptsächlich [Teu03], [Böl09] und [PST07]. Zu-nächst soll in Kapitel 2 in die Thematik der adiabatischen Störungstheorie eingeführt

4

Page 13: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

werden, beginnend mit der zeitadiabatischen Theorie. Nachdem wir alle nötigen ma-thematischen Voraussetzungen und Begrifflichkeiten eingeführt haben, sind wir in derLage, das klassische zeitadiabatische Theorem 2.3 zu formulieren und im Anschlusszu beweisen. Dieses werden wir dann innerhalb der raumadiabatischen Theorie zu demraumadiabatischen Theorem 2.6 verallgemeinern, bevor wir uns der effektiven Dynamikwidmen. In Kapitel 3 stellen wir die zeitabhängige Born-Oppenheimer-Approximationvor. Zu diesem Zweck erklären wir die konventionelle Vorgehensweise, indem wir diemathematischen Ergebnisse aus Kapitel 2 auf molekulare Dynamik übertragen. In die-sem Zuge berechnen wir dann auch den effektiven Hamiltonoperator und präsentierenim darauf folgenden Abschnitt ein systematisches Schema für Korrekturen höherer Ord-nungen. Die zeitabhängige Born-Oppenheimer-Approximation zweiter Ordnung leitenwir dabei explizit her. In einem abschließenden Beispielmodell in Kapitel 4 soll derMehrwert der korrekten Korrekturterme verdeutlicht werden. Nachdem wir das zwei-dimensionale Modell vorgestellt haben, gehen wir auf die Umsetzung wie zum Beispieldie verwendeten numerischen Methoden und die Implementierung ein. Zum Schlusspräsentieren wir die Ergebnisse der Durchführung und ziehen ein abschließendes Fazit.

Danksagung

An dieser Stelle will ich die Gelegenheit nutzen, den Personen zu danken, die mich beider Entstehung dieser Arbeit unterstützt und begleitet haben.

Zuallererst geht mein Dank an Prof. Stefan Teufel, der mir nicht nur die Möglichkeitgab, meine Masterarbeit in seiner Arbeitsgruppe zu schreiben, sondern darüber hinausein herausragender Betreuer für mich war. Mit seiner Hilfsbereitschaft und seinem En-gagement hat er maßgeblich zum Gelingen dieser Arbeit beigetragen.

Ein außerordentlicher Dank geht zudem an Mario Laux, der mit seinem breitenmathematischen und physikalischen Wissen eine große Unterstützung für mich war.Er hatte immer ein offenes Ohr für mich und nahm sich stets Zeit, mir bei meinenProblemen zu helfen. Das war keinesfalls selbstverständlich und ich weiß das sehr zuschätzen.

Ich danke allen meinen Freunden, die ich in Tübingen gefunden habe und die meinStudium zu dem gemacht haben, was es war. Insbesondere geht mein Dank an ValentinBolz und Sebastian Künkele, die mich durch das gesamte Studium begleitet und es zueiner Zeit gemacht haben, an die ich immer gerne zurück denken werde.

Meinen Eltern Brunhilde und Waldemar sowie meiner Schwester Marina und meinerGroßmutter Elfriede gebührt ganz besonderer Dank. Durch die finanzielle und vor allememotionale Unterstützung haben sie mir mein Studium erst ermöglicht und gaben mirRückhalt, um es erfolgreich zu absolvieren. Danke für eure fürsorgliche, warme undherzliche Art.

Zu guter Letzt danke ich meiner Freundin Carmen Müller, die selbst in den Mo-menten an mich geglaubt hat, als ich es nicht getan habe. Du bist eine Stütze für mich,auf die ich mich immer verlassen kann. Dafür danke ich dir von ganzem Herzen.

5

Page 14: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls
Page 15: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

2 Adiabatische Störungstheorie fürdie Schrödinger Gleichung

Ziel dieses Kapitels ist es, die adiabatische Störungstheorie vorzustellen und die für dieBorn-Oppenheimer-Approximation nötigen Grundlagen zu schaffen. Am Anfang gebenwir eine kurze Übersicht über die Entstehungsgeschichte der adiabatischen Theorie inder Quantenmechanik, bevor wir die mathematischen Voraussetzungen für die zeitadia-batische Theorie angehen. Nachdem wir in diesem Zuge das zeitadiabatische Theoremder Quantenmechanik studiert haben, gehen wir zur raumadiabatischen Störungstheorieüber. Das darin präsentierte raumadiabatische Theorem ist von wesentlicher Bedeutungfür das Verständnis der Born-Oppenheimer-Approximation in Kapitel 3. Die Basis fürdieses Kapitels bildet [Teu03] in Verbindung mit Teilen aus [Böl09].

2.1 Das zeitadiabatische Theorem der Quantenmechanik

Das klassische zeitadiabatische Theorem ist das einfachste seiner Art und gleichzeitigdie Grundlage der adiabatischen Theoreme in der Quantenmechanik. Das Präfix zeitwird des öfteren weggelassen, sodass mit dem adiabatischen Theorem meist das zeita-diabatische Theorem gemeint ist. Es geht zurück auf Born und Fock [BF28], welche 1928das Theorem von Ehrenfest [Ehr16] von 1916 erweiterten. Eine elegante mathematischeFormulierung gelang allerdings erst Tosio Kato [Kat50] im Jahre 1950.

Die Theorie um das zeitadiabatische Theorem beschäftigt sich mit Quantensyste-men, die von einem langsam von der Zeit abhängigen Hamiltonoperator beschriebenwerden. Die Zeitabhängigkeit des Hamiltonoperators geht oft von der Zeitabhängigkeiteines externen Parameters aus. Das Theorem entspricht in der Form, wie wir es prä-sentieren, nicht unbedingt der üblichen Darstellung, aber erlaubt uns die bestmöglicheVerallgemeinerung auf das später folgende raumadiabatische Theorem. Das klassischezeitadiabatische Theorem sagt aus, dass wenn sich das System zur Zeit t = t0 in ei-nem Eigenzustand des zugrundeliegenden Hamiltonoperators H(t0) befindet und zu-dem gewisse Voraussetzungen erfüllt sind, das System zum Zeitpunkt t = T in einenEigenzustand des Hamiltonoperators H(T ) übergegangen ist. Wichtig hierfür sind un-ter anderem, dass die Veränderung des Systems hinreichend langsam, also adiabatisch,ist und eine spektrale Lücke existiert, die ein Teil vom Spektrum von dem Rest isoliert.

Im Nachfolgenden werden wir die mathematischen Grundlagen einführen, die fürdas Verständnis dieser Theorie nötig sind. Weitere Details können zudem in [Teu03]nachgeschlagen werden.

Page 16: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 2. Adiabatische Störungstheorie für die Schrödinger Gleichung

2.1.1 Mathematische Grundlagen

Im weiteren Verlauf der Arbeit sei ‖·‖H die von einem Skalarprodukt 〈·, ·〉H induzierteNorm und dist(·, ·) die induzierte Metrik in einem Hilbertraum H. Sei A : H → H einlinearer Operator auf dem Hilbertraum H, dann bezeichnet

‖A‖ := sup {‖Av‖H | v ∈ H mit ‖v‖H ≤ 1}

die Operatornorm von A. In den meisten Fällen wird für uns H der Raum L2(Λ) derquadratintegrierbaren Funktionen auf dem Raum Λ = Rd sein. Dann ist

‖ϕ‖H =(∫

Λϕ(x)ϕ(x) dx

) 12

=(∫

Λ|ϕ(x)|2 dx

) 12,

für ϕ ∈ L2(Λ) und ϕ die komplexe Konjugation von ϕ. Für zwei Operatoren A und Bsei [A,B] = AB −BA deren Kommutator.

Sei also H der Hilbertraum L2(Rd) mit zugrunde liegendem Konfigurationsraum Rd.Wir betrachten ein Quantensystem, welches durch eine Wellenfunktion ψε : Rt → Hbeschrieben wird. Der Zustand zum Zeitpunkt t ∈ R von ψε wird bestimmt durch diezeitabhängige Schrödingergleichung. Diese lässt sich als Anfangswertproblem auffassen:

i~ ddtψ

ε(t) = Hψε(t), ψε(t0) = ψε0 ∈ DH . (2.1)

Hierbei ist H der selbstadjungierte Hamiltonoperator mit dichter Domäne DH ⊆ H, derals Folgerung aus dem Spektralsatz eine stark stetige einparametrige unitäre GruppeU(t) = e−i/~Ht ∈ L(H) mit t ∈ R erzeugt. Die Domäne kann formalisiert werden durch

DH = {ϕ ∈ H | t→ U(t)ϕ ist differenzierbar} .

Beachte, dass der exponentielle Ausdruck nur symbolisch zu verstehen ist, da H nichtzwangsweise beschränkt ist. Die eindeutige globale Lösung von (2.1) ist gegeben durch

ψε(t) = e−i~Htψε0.

Wir bezeichnen mit P : H → ker(H − E) ⊆ H die spektrale Projektion auf den Eigen-raum zu einem Eigenwert E von H. Es ist leicht zu sehen, dass P mit H kommutiert,denn es ist

PHϕ−HPϕ = PHϕ− EPϕ = P (Hϕ− Eϕ) ∀ϕ ∈ DH .

Außerdem gilt für alle ϕ′ ∈ ker(H − E), dass 〈ϕ, (H − E)ϕ′〉 = 0 und mit der Selbst-adjungiertheit von H, also insbesondere E ∈ R, auch 〈Hϕ − Eϕ,ϕ′〉 = 0. Somit istHϕ− Eϕ ∈ ker(H − E)⊥ = ker(P ).

Folglich ist also PH ⊂ H ein H-invarianter Unterraum von H und für einen Eigen-zustand ψε0 ∈ PH, liegt ebenso ψε(t) = e−i/~Etψε0 für alle t ∈ R im Unterraum PH. DieProjektion P kommutiert also auch mit U(t) = e−i/~Ht.

Was ändert sich nun, wenn sowohl H als auch P und E von der Zeit abhängen? Seidazu H(s), s ∈ R eine Familie selbstadjungierter Operatoren auf einem Hilbertraum H.Wir wollen (2.1) entsprechend anpassen und benötigen dazu folgende Definition.

8

Page 17: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

2.1. Das zeitadiabatische Theorem der Quantenmechanik

Definition 2.1 (unitärer Propagator). Ein unitärer Propagator U : R2 → L(H) ist eine(stark stetige) Familie unitärer Operatoren (U(s, t))s,t∈R mit

(i) U(s, r)U(r, t) = U(s, t) für alle s, r, t ∈ R,(ii) U(s, s) = IdH für alle s ∈ R,(iii) U ist in beiden Variablen gleichzeitig stark stetig, das heißt

limtn→tsn→s

U(tn, sn)ϕ = U(t, s)ϕ für alle ϕ ∈ H.

Wir sind nun interessiert an der Lösung des Anfangswertproblems

i~ ddsU

εf (s, s0) = H(εs)U εf (s, s0), U εf (s0, s0) = IdH. (2.2)

Denn löst der zu H(t) korrespondierende unitäre Propagator U εf (s, s0) die Gleichung(2.2), dann ist ψε(s) = U εf (s, s0)ψε0 eine Lösung von

i~ ddsψ

ε(s) = H(εs)ψε(s), ψε(s0) = ψε0 ∈ DH .

Bevor wir fortfahren, sei auf die oben verwendeten Indizes „ε“ und „f“ eingegan-gen. Der Index ε verdeutlicht die Abhängigkeit eines Operators beziehungsweise derLösung einer Gleichung von dem Parameter ε mit 0 < ε � 1. Er wird adiabatischerParameter genannt und nimmt eine große Rolle in der adiabatischen Theorie ein, daer die Trennung der Zeitskalen, auf der sich das System verändert, kontrolliert. Dielangsame beziehungsweise makroskopische Zeitskala bezeichnen wir mit t = εs und mits sei die schnelle Skala respektive mikroskopische Zeiteinheiten gemeint. Je kleiner ε,desto langsamer geht die Veränderung des Hamiltonoperators also vonstatten. Die Be-zeichnung mikroskopisch und makroskopisch wirkt momentan noch etwas aus der Luftgegriffen, wird allerdings in Abschnitt 2.2.1 klarer. Operatoren und Funktionen, die aufder schnellen Zeitskala wirken, kennzeichnen wir mit dem Index „f“, stellvertretend fürdas englische Wort „fast“.

Die langsame Zeitskala ist für uns von größerer Bedeutung, da wir die Dynamikfür ausreichend lange mikroskopische Zeiten beobachten müssen, um nicht nur trivialeBewegungen festzustellen. Dementsprechend müssen wir Anpassungen vornehmen. DieAbbildung τ : s 7→ s/ε und damit einhergehend eine Änderung der Variablen s = t/εtransformiert die obige Gleichung (2.2) auf die langsame Zeitskala:

i~ε ddtU

ε(t, t0) = H(t)U ε(t, t0), U ε(t0, t0) = IdH. (2.3)

H(t) und P (t), also die spektrale Projektion auf den Eigenraum zum Eigenwert E(t) vonH(t), sind nun unabhängig von ε, allerdings bleibt die Lösung U ε(t, t0) = U εf (t/ε, t0/ε)weiterhin von ε abhängig, was wieder mit dem Index verdeutlicht wird.

Die Existenz einer Lösung der Schrödingergleichung wird mittels folgender Propo-sition gesichert, welche auch in [Teu03] nachzulesen ist. Diese ist eine direkte Folgerungaus dem allgemeinen Ergebnis für Kontraktionshalbgruppen, zu finden in [RS75] Theo-rem X.70.

9

Page 18: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 2. Adiabatische Störungstheorie für die Schrödinger Gleichung

Proposition 2.2. Sei J ⊆ R ein offenes Intervall und H(t), t ∈ R eine Familie vonselbstadjungierten Operatoren auf einem Hilbertraum H mit einer gemeinsamen dichtenDomäne D ⊂ H, ausgestattet mit der Graphennorm von H(t) für ein t ∈ J , sodassH(·) ∈ C1

b(J,L(D,H)).Dann existiert ein unitärer Propagator U ε, sodass für t, t0 ∈ J und ψε0 ∈ D eine

Lösung der zeitabhängigen Schrödingergleichung

i~ε ddtψ

ε(t) = H(t)ψε(t), ψε(t0) = ψε0 (2.4)

gegeben ist durch ψε(t) = U ε(t, t0)ψε0.

Es sei angemerkt, dass wie im einparametrigen autonomen Fall auch, U(s, t) eineneindeutig korrespondierenden Erzeuger H(t) hat und dass U(s, t) die Domäne D vonH(t) invariant lässt, das heißt U(s, t)D = D ist.

Auch für den nicht autonomen Fall gilt [P (t), H(t)] = 0. Jedoch kommutiert P (t)nun nicht mehr mit dem unitären Propagator U ε(t, t0). Für t 6= t0 ist

P (t)U ε(t, t0) 6= U ε(t, t0)P (t0) bzw. U ε(t0, t)P (t)U ε(t, t0) 6= P (t0).

Dies wird bei beidseitiger Betrachtung der Ableitungen nach t deutlich.

ddt(U ε(t0, t)P (t)U ε(t, t0)

)(2.5)

= i~εU ε(t0, t)[H(t), P (t)]U ε(t, t0) + U ε(t0, t)P (t)U ε(t, t0)

= U ε(t0, t)P (t)U ε(t, t0) 6= 0.

Hierbei wurde (2.3) und daraus resultierend

ddtU

ε(t0, t) = ddt(U ε(t, t0)∗

)= i

~εU ε(t0, t)H(t)

genutzt. Für ψε0 ∈ P (t0)H ist U ε(t, t0)ψε0 also nicht notwendigerweise in dem spektralenUnterraum P (t)H. Dies gilt nur approximativ, wenn P (t) klein ist. Das bereits erwähn-te klassische zeitadiabatische Theorem macht eine Aussage über diese Eigenschaft. DerPropagator U ε(t, t0) transportiert approximativ den zeitabhängigen spektralen Unter-raum von H(t0), also dessen Eigenzustände, wenn diese sich hinreichend langsam mitder Zeit verändern.

2.1.2 Das zeitadiabatische Theorem

Als Folgerung aus den vorangegangenen Resultaten betrachten wir spektrale Unterräu-me, in denen Teile des Spektrums durch eine Lücke vom Rest getrennt sind. Genau-er gesagt sei σ(H(t)) = {E ∈ C | (H(t)− E · Id) ist nicht beschränkt invertierbar} dasSpektrum von H(t), welches eine Teilmenge σ∗(t) ⊂ σ(H(t)) enthält, sodass es zwei be-schränkte stetige Funktionen f± ∈ Cb(R,R) gibt, die ein Intervall I(t) = [f−(t), f+(t)]definieren mit

σ∗(t) ⊂ I(t) und inft∈R

dist(I(t), σ(H(t))\σ∗(t)

)=: g > 0. (2.6)

10

Page 19: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

2.1. Das zeitadiabatische Theorem der Quantenmechanik

f+(t)

σ∗(t)

f−(t)

σ(H(t))

t

Abbildung 2.1: Durch eine Lücke (weiß) lokal getrenntes Spektrum σ∗(t).

Abbildung 2.1 veranschaulicht diese Eigenschaft.Sei P∗(t) die spektrale Projektion auf den zu σ∗(t) korrespondierenden Eigenraum

von H(t). Das zeitadiabatische Theorem der Quantenmechanik in seiner einfachstenForm sagt nun aus, dass bei hinreichender Regularität von H(t) eine Konstante C <∞unabhängig von ε existiert, sodass∥∥∥P⊥∗ (t)U ε(t, t0)P∗(t0)

∥∥∥L(H)

≤ Cε(1 + |t− t0|) (2.7)

mit der zu P∗(t) orthogonalen Projektion P⊥∗ (t) = Id − P∗(t). In Worten heißt dies:Gilt für den Anfangszustand ψε0 ∈ P∗(t0)H, dann bleibt der Zustand zu einem späterenZeitpunkt ψ(t) gegeben durch die Lösung von (2.4) bis auf einen Fehler der OrdnungO(ε(1 + |t− t0|)‖ψε0‖) in dem Unterraum P∗(t)H. Die analoge Aussage gilt, wenn wirin dem orthogonalen Komplement von P∗(t)H beginnen. Man nennt dies adiabatischeInvarianz getrennter spektraler Unterräume. Der Vorgang, bei dem ein Raum in spek-trale Unterräume getrennt wird, die in irgendeiner Weise langsam von einem Parameterabhängen und approximativ invariant unter der quantenmechanischen Zeitentwicklungsind, heißt adiabatische Entkopplung.

Das adiabatische Theorem wird häufig in obiger Form (2.7) angegeben, allerdingsist die Aussage aus dem Beweis von Kato [Kat50] stärker. In seiner Arbeit definiertKato den adiabatischen Hamiltonoperator

Ha(t) := H(t) + i~ε[P∗(t), P∗(t)

](2.8)

zusammen mit dem adiabatischen Propagator U εa (t, t0), der als Lösung der Gleichung

i~ε ddtU

εa (t, t0) = Ha(t)U εa (t, t0), U εa (t0, t0) = IdH (2.9)

gegeben ist. Der adiabatische Propagator löst (2.9) auf einer gemeinsamen dichten Do-mäne D, vergleiche Proposition 2.2. Er ist so definiert, dass er die spektralen Teilräume

11

Page 20: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 2. Adiabatische Störungstheorie für die Schrödinger Gleichung

zu unterschiedlichen Zeiten exakt ineinander überführt, das heißt

P∗(t)U εa (t, t0) = U εa (t, t0)P∗(t0) für alle t, t0 ∈ R. (2.10)

Infolgedessen lässt sich analog zu (2.5) mit dem Korrekturterm

M(t) := i~ε[P∗(t), P∗(t)]

berechnen, dass

ddt(U εa (t0, t)P∗(t)U εa (t, t0)

)= i

~εU εa (t0, t)[Ha(t), P∗(t)]U εa (t, t0) + U εa (t0, t)P∗(t)U εa (t, t0)

= i~εU εa (t0, t)[H(t), P∗(t)]U εa (t, t0) + i

~εU εa (t0, t)[M(t), P∗(t)]U εa (t, t0)

+ U εa (t0, t)P∗(t)U εa (t, t0)

= i~εU εa (t0, t)[M(t), P∗(t)]U εa (t, t0) + U εa (t0, t)P∗(t)U εa (t, t0)

= U εa (t0, t)(P∗(t)− P∗(t)

)U εa (t, t0) = 0.

Dabei wurde genutzt, dass

i~ε

[M(t), P∗(t)] = i~ε

[i~ε[P∗(t), P∗(t)

], P∗(t)

]= −P∗(t).

Dies wird bei genauerer Betrachtung der ersten Ableitung von P∗(t) nach t klar:

P∗(t) = ddt(P∗(t)

)2 = P∗(t)P∗(t) + P∗(t)P∗(t).

Bemerke außerdem, dass dies

P∗(t)P∗(t)P∗(t) = P∗(t)P∗(t)P∗(t) + P∗(t)P∗(t)P∗(t) = 2P∗(t)P∗(t)P∗(t)

und daher P∗(t)P∗(t)P∗(t) = 0 impliziert. Die analoge Aussage gilt für P⊥∗ (t). Letzt-endlich ist somit

P∗(t) =(P∗(t) + P⊥∗ (t)

)P∗(t)

(P∗(t) + P⊥∗ (t)

)= P∗(t)P∗(t)P∗(t) + P∗(t)P∗(t)P⊥∗ (t) + P⊥∗ (t)P∗(t)P∗(t) + P∗(t)⊥P∗(t)P⊥∗ (t)= P∗(t)P∗(t)P⊥∗ (t) + P⊥∗ (t)P∗(t)P∗(t)

Durch Multiplikation mit P∗(t) und P⊥∗ (t) erkennen wir die Ergänzungen

P∗(t)P∗(t) = P∗(t)P∗(t)P⊥∗ (t) und P∗(t)P∗(t) = P⊥∗ (t)P∗(t)P∗(t) (2.11)

ebenso wie

P⊥∗ (t)P∗(t) = P⊥∗ (t)P∗(t)P∗(t) und P∗(t)P⊥∗ (t) = P∗(t)P∗(t)P⊥∗ (t). (2.12)

Wir können nun die starke Version des zeitadiabatischen Theorems formu-lieren.

12

Page 21: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

2.1. Das zeitadiabatische Theorem der Quantenmechanik

Theorem 2.3 (Kato 1950). Sei H(·) ∈ C2b(J,L(D,H)) eine Familie selbstadjungierter

Operatoren, die die Voraussetzungen aus Proposition 2.2 erfüllt. Des Weiteren genügeσ∗(t) ⊂ σ(H(t)) der Lückenbedingung (2.6) auf J ⊂ R.

Dann ist P∗(·) ∈ C2b(J,L(H)) und es existiert eine Konstante C < ∞, sodass für

t, t0 ∈ J‖U ε(t, t0)− U εa (t, t0)‖L(H) ≤ Cε(1 + |t− t0|). (2.13)

Die Aussage (2.13) impliziert aufgrund von (2.10) offensichtlich (2.7), da wir Ha sokonstruiert haben, dass ‖P⊥∗ (t)U εa (t, t0)P∗(t0)‖ = 0 ist und somit∥∥∥P⊥∗ (t)U ε(t, t0)P∗(t0)

∥∥∥L(H)

≤∥∥∥P⊥∗ (t)U εa (t, t0)P∗(t0)

∥∥∥L(H)

+∥∥∥P⊥∗ (t)

(U ε(t, t0)− U εa (t, t0)

)P∗(t0)

∥∥∥L(H)

≤ ‖U ε(t, t0)− U εa (t, t0)‖L(H).

Es ergeben sich somit nicht nur eine approximative Invarianz des spektralen Unter-raums, sondern zusätzlich Informationen über die effektive Zeitentwicklung innerhalbdes entkoppelten Teilraums.

Wir wollen das zeitadiabatische Theorem nun beweisen. Die Beweisidee geht zurückauf [Kat50] und die nachfolgend verwendete Darstellung entspricht der aus [Teu03].

Beweis von Theorem 2.3. Zunächst zeigen wir, dass die Regularität von H(t) als Funk-tion von t mit der Lückenbedingung (2.6) auf die Regularität der spektralen ProjektionP∗(t) übertragen wird. Aufgrund der Lückenbedingung existiert eine Familie von positivorientierten geschlossenen Kurven (Γ (t))t∈J ⊂ C, die jeweils σ∗(t) umkreisen, sodass

inft∈J

dist(Γ (t), σ(H(t)

)= g

2 .

Wir nutzen Rieszs’ Projektionsformel (vgl. [HS96]),

P∗(t) = i2π

∮Γ (t)

R(ζ, t) dζ ,

wobei R(ζ, t) := (H(t) − ζId)−1. Aufgrund der Stetigkeit der Funktionen f± in (2.6)existiert für jedes ξ ∈ J eine Umgebung I(ξ) von ξ, sodass

inft∈I(ξ)

dist(Γ (ξ), σ(H(t)

)≥ g

4 .

Für t ∈ I(ξ) erhalten wir also die Darstellung

P∗(t) = i2π

∮Γ (ξ)

R(ζ, t) dζ . (2.14)

Leiten wir (2.14) für t ∈ I(ξ) nach t ab, ergibt sich

dn

dtnP∗(t) = i2π

∮Γ (ξ)

dn

dtnR(ζ, t) dζ , (2.15)

falls R(ζ, ·) ∈ Cnb (J,L(H)). Dies wiederum sichert uns folgendes Lemma:

13

Page 22: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 2. Adiabatische Störungstheorie für die Schrödinger Gleichung

Lemma 2.4. Sei H(·) ∈ Cnb (I,L(D,H)) für ein offenes Intervall I ⊂ J . Dann giltR(ζ, ·) ∈ Cnb (I,L(H,D)) für alle ζ ∈ C\∪t∈Iσ(H(t)).

Beweis. Mit

H(t) = H(t)− ζId + ζId ⇔ H(t)R(ζ, t) = Id + ζR(ζ, t) (2.16)⇔ R(ζ, t)H(t) = Id +R(ζ, t)ζ,

erhalten wir

R(ζ, t)−R(ζ, s) = R(ζ, t)(Id + ζR(ζ, s)

)−(Id +R(ζ, t)ζ

)R(ζ, s)

= R(ζ, t)H(s)R(ζ, s)−R(ζ, t)H(t)R(ζ, s)= R(ζ, t)(H(s)−H(t))R(ζ, s).

Es folgt

‖R(ζ, t)−R(ζ, s)‖ ≤ ‖R(ζ, t)(H(s)−H(t))R(ζ, s)‖≤ ‖R(ζ, t)‖ · ‖(H(s)−H(t))‖ · ‖R(ζ, s)‖.

Mit R(ζ, t) ∈ L(H,D) für t ∈ I und H(t) ∈ L(D,H) gilt also R(ζ, ·) ∈ Cb(I,L(H,D)).Für die Differenzierbarkeit betrachten wir die Ableitung der Identität

Id =(H(t)− ζ

)R(ζ, t)

nach t:R(ζ, t) = −R(ζ, t)H(t)R(ζ, t). (2.17)

Nach Voraussetzung ist H(·) ∈ Cnb (I,L(D,H)) und mit R(ζ, ·) ∈ Cb(I,L(H,D)) wirdR(ζ, ·) ∈ Cnb (I,L(H,D)) impliziert.

Wenden wir Lemma 2.4 für n = 2 auf die Gleichung (2.15) an, erhalten wir diegewünschte Regularität von P∗(t), genauer P∗(·) ∈ C2

b(J,L(H)). Das bedeutet ins-besondere, dass Ha(t) den Voraussetzungen aus Proposition 2.2 genügt und dadurchdefiniert U εa einen Propagator auf D.

Wir wollen nun Gleichung (2.13) zeigen. Dazu drücken wir die Differenz der Propa-gatoren mithilfe von (2.3), (2.8) und (2.9) als Differenz ihrer Erzeuger aus.

U ε(t, t0)− U εa (t, t0) = −U ε(t, t0)∫ t

t0

ddt′(U ε(t′, t0)U εa (t′, t0)

)dt′

= − i~εU ε(t, t0)

∫ t

t0U ε(t′, t0)

(H(t′)−Ha(t′)

)U εa (t′, t0) dt′

= −U ε(t, t0)∫ t

t0U ε(t′, t0)

[P∗(t′), P∗(t′)

]U εa (t′, t0) dt′ . (2.18)

Wir haben außerdem benutzt, dass U εa (t, t0)D = D für alle t, t0 ∈ J . Eine erste naiveVermutung ergibt, dass der Integrand in (2.18) von der Ordnung O(1) ist. Die Kernideeist allerdings, dass dieser oszilliert und als Zeitableitung einer Funktion der Ordnung

14

Page 23: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

2.1. Das zeitadiabatische Theorem der Quantenmechanik

O(ε) dargestellt werden kann. Dies wird beispielsweise anhand der Funktion eit/ε klar.Sie ist von der Ordnung O(1), jedoch erhalten wir durch Integration

∥∥∥∥∫ t

t0eit′/ε

∥∥∥∥ =∥∥∥−iε

(eit/ε − eit0/ε

)∥∥∥ = O(ε).

Um die folgende Definition zu motivieren, nehmen wir für den Moment den Spezi-alfall σ∗(t) = {E(t)} an, also dass σ∗(t) nur aus einem einzigen Eigenwert E(t) besteht.Setzen wir

F (t) = −(R(E(t), t)P∗(t)P∗(t) + P∗(t)P∗(t)R(E(t), t)

)und benutzen (

H(t)− E(t))P∗(t) = P∗(t)

(H(t)− E(t)

)= 0,

sehen wir

[H(t), F (t)] = [H(t)− E(t), F (t)]= −

(H(t)− E(t)

)R(E(t), t

)P∗(t)P∗(t)−

(H(t)− E(t)

)P∗(t)P∗(t)R

(E(t), t

)+R(E(t), t)P∗(t)P∗(t)(H(t)− E(t)) + P∗(t)P∗(t)R(E(t), t)(H(t)− E(t))

= −P∗(t)P∗(t)−(H(t)− E(t)

)P∗(t)P∗(t)R

(E(t), t

)+R

(E(t), t

)P∗(t)P∗(t)

(H(t)− E(t)

)+ P∗(t)P∗(t)

=[P∗(t), P∗(t)

]. (2.19)

Wir müssen beachten, dass die Resolvente eigentlich nicht für Elemente aus dem Spek-trum definiert ist. Unter Berücksichtigung von (2.12) operiert die reduzierte ResolventeR(E(t), t) allerdings nur auf P⊥∗ (t)H und ist aufgrund der Lückenbedingung gleichmäßigbeschränkt,

∥∥∥R(E(t), t)P⊥∗ (t)

∥∥∥L(H,D)

≤ 1dist

(E(t), σ

(H(t)|P⊥∗ (t)

)) = 1g. (2.20)

Ein Teil des Integranden aus (2.18) kann also nach (2.19) als Kommutator von H(t)und F (t) geschrieben werden. Dies wollen wir nun auch für den allgemeinen Fall einesbeliebigen σ∗(t) zeigen. Wir definieren dazu

F (t) := 12πi

∮Γ(t)

P⊥∗ (t)R(ζ, t)R(ζ, t) dζ + adj.,

wobei hier und im Folgenden „adj.“ die Adjungierte aller links davon stehenden Termebezeichnet. Mit (2.15), (2.16) und R(ζ, t)(H(t)−ζ) = −R(ζ, t)H(t) als Umformung von

15

Page 24: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 2. Adiabatische Störungstheorie für die Schrödinger Gleichung

(2.17) berechnen wir nun für den allgemeinen Fall

[H(t), F (t)] = 12πi

∮Γ(t)

[H(t), P⊥∗ (t)R(ζ, t)R(ζ, t)

]dζ − adj.

= 12πi

∮Γ(t)

P⊥∗ (t)(H(t)R(ζ, t)R(ζ, t)−R(ζ, t)R(ζ, t)H(t)

)dζ − adj.

= 12πi

∮Γ(t)

P⊥∗ (t)((

Id + ζR(ζ, t))R(ζ, t)−R(ζ, t)R(ζ, t)H(t)

)dζ − adj.

= 12πi

∮Γ(t)

P⊥∗ (t)R(ζ, t) + P⊥∗ (t)R(ζ, t)R(ζ, t)(ζ −H(t)

)dζ − adj.

= P⊥∗ (t)P∗(t)−1

2πi

∮Γ(t)

P⊥∗ (t)R(ζ, t)R(ζ, t)(H(t)− ζ

)dζ − adj.

= P⊥∗ (t)P∗(t) + 12πi

∮Γ(t)

P⊥∗ (t)R2(ζ, t)H(t) dζ − adj.

= P⊥∗ (t)P∗(t)− adj.

=[P∗(t), P∗(t)

].

Bei der vorletzten Gleichheit haben wir ausgenutzt, dass ζ → P⊥∗ (t)R2(ζ, t) holomorphinnerhalb von Γ(t) ist und im letzten Schritt haben wir (2.11) und (2.12) verwendet.

Wie bereits erwähnt, wollen wir die Gleichheit der beiden Kommutatoren [H(t), F (t)]und

[P∗(t), P∗(t)

]in dem Integrand in (2.18) ausnutzen. Sei dazu

A(t) := −i~εU ε(t0, t)F (t)U ε(t, t0).

Zusammen mit P∗(·) ∈ C2b(R,L(H)) und der Definition von F (t) ist F (·) ∈ C1

b(R,L(H))und wir berechnen analog zu (2.5), dass

ddtA(t) = U ε(t0, t)[H(t), F (t)]U ε(t, t0)− i~εU ε(t0, t)F (t)U ε(t, t0) (2.21)

= U ε(t0, t)[P∗(t), P∗(t)

]U ε(t, t0) +O(ε).

Mit dieser Schreibweise meinen wir ‖i~εU ε(t0, t)F (t)U ε(t, t0)‖ = O(ε). Nun wissen wiralso, dass wir (2.18) approximativ als Integral über die Zeitableitung von A(t) schreibenkönnen. Daher setzen wir (2.21) in (2.18) ein und erhalten mit partieller Integration

‖U ε(t, t0)− U εa (t, t0)‖

=∥∥∥∥U ε(t, t0)

∫ t

t0

( ddt′A(t′)

)U ε(t0, t′)U εa (t′, t0) + i~εU ε(t0, t)F (t)U εa (t′, t0) dt′

∥∥∥∥≤∥∥∥∥∫ t

t0

( ddt′A(t′)

)U ε(t0, t′)U εa (t′, t0) dt′

∥∥∥∥+ ~ε|t− t0| sups∈[t0,t]

‖F (t)‖

≤∥∥∥∥∥[A(t′)U ε(t0, t′)U εa (t′, t0)

∣∣∣∣tt0

+∫ t

t0A(t′)

( ddt′U

ε(t0, t′)U εa (t′, t0))

dt′∥∥∥∥∥+O(ε|t− t0|)

≤ ‖A(t)‖+ ‖A(t0)‖+∥∥∥∥∫ t

t0A(t′)

( ddt′U

ε(t0, t′)U εa (t′, t0))

dt′∥∥∥∥+O(ε|t− t0|)

≤ O(ε (1 + |t− t0|)).

16

Page 25: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

2.2. Raumadiabatische Störungstheorie

Dafür haben wir genutzt, dass ‖A(t′)‖ = O(ε) und∥∥∥∥ ddt′U

ε(t0, t′)U εa (t′, t0)∥∥∥∥ =

∥∥∥U ε(t0, t′)(H(t)−Ha(t))U εa (t′, t0)

∥∥∥=∥∥∥U ε(t0, t′)[P∗(t′), P∗(t′)]U εa (t′, t0)

∥∥∥uniform beschränkt ist. Damit ist die Aussage des Theorems bewiesen.

Mit der Definition des adiabatischen Hamiltonoperators Ha(t) entsteht die Mög-lichkeit eine effektive Schrödingergleichung aufzustellen, die eine simplere Form als dieSchrödingergleichung des Gesamtsystems hat und uns gestattet, eine Lösung approxi-mativ zu bestimmen. Dies geschieht, indem wir zuvor einen effektiven HamiltonoperatorHeff(t) aus Ha(t) herleiten, der die Dynamik auf dem Unterraum P∗(t)H beschreibt. Aufdie nachfolgende raumadiabatische Theorie, welche den relevanten Hintergrund dieserArbeit bildet, ist diese Taktik ebenso anwendbar. In diesem Setting werden wir späterdie effektive Dynamik noch einmal genauer betrachten.

2.2 Raumadiabatische StörungstheorieIn Abschnitt 2.1 haben wir zeitabhängige Hamiltonoperatoren betrachtet. In der Praxisgibt es jedoch auch viele relevante physikalische Systeme, die nicht direkt von der Zeitabhängen. Somit sind die jeweiligen korrespondierenden Hamiltonoperatoren zeitunab-hängig. In der raumadiabatischen Störungstheorie übernehmen Ortsvektoren die Rolleder Zeit t aus der gewöhnlichen zeitadiabatischen Theorie. Eine Skalierung des Konfigu-rationsraums erlaubt es dennoch häufig, unterschiedliche Zeitskalen zu nutzen. Ist diesder Fall, ist eine adiabatische Entkopplung möglich. Folglich gelingt es, Freiheitsgradezu finden, die sich langsamer mit der Zeit entwickeln als andere, sodass diese voneinan-der getrennt werden können. Ein bekanntes Beispiel dafür ist die Born-Oppenheimer-Approximation, auf die wir hinarbeiten und in Kapitel 3 genauer eingehen werden. Siedient zur Vereinfachung von Systemen mit mehreren Teilchen, um Aussagen über de-ren Dynamik zu erhalten. Eine essentielle Rolle spielt dabei das Massenverhältnis vonKernen und Elektronen in Molekülen und die damit verbundenen unterschiedlichenkinetischen Energien. Daher berücksichtigen wir im Folgenden nur Quantensysteme,deren Gesamtenergie beschränkt ist, um eine Trennung der Skalen zu erzwingen.

2.2.1 Rahmenbedingungen

Ein sich in Rd+q entwickelndes Quantensystem können wir mithilfe der zeitabhängigenSchrödingergleichung mit zeitunabhängigem Potential beschreiben. Wir fassen diesedazu als partielle Differentialgleichung mit den Ortsvektoren x ∈ Rd und y ∈ Rq sowieder mikroskopischen Zeit s ∈ R als Argumente einer Wellenfunktion Ψ auf:

i~ ∂∂sΨ(s, x, y) = HΨ(s, x, y). (2.22)

Hierbei ist der selbstadjungierte Hamiltonoperator von der Form

H = − ~2

2m∆x −~2

2m∆y + 1ε2V

(x,y

ε

). (2.23)

17

Page 26: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 2. Adiabatische Störungstheorie für die Schrödinger Gleichung

Die ersten beiden Terme von H beschreiben die kinetischen Energien des Systems in x-und y-Richtung. Der hintere Term mit dem Potential V (x, y/ε) schränkt die Bewegungdes Teilchens ein. Mittels einer Skalierung durch ε in y-Richtung wird es gezwungen,nahe dem Unterraum y = 0 zu bleiben und die Wellenfunktion oszilliert entlang dery-Richtung stärker als entlang der x-Richtung. Das System kann dadurch später inlangsame, durch x repräsentierte und schnelle, durch y repräsentierte Freiheitsgradegetrennt werden. Fürs Erste bleiben wir in dem allgemeinen Setting, in dem xDimensiond und y Dimension q hat. In Kapitel 4 wenden wir diese Ergebnisse dann auf einzweidimensionales Modell an.

Die langsamen Freiheitsgrade bewegen sich voraussichtlich mit Geschwindigkeitender Ordnung ε, womit wir deren Bewegung über mikroskopische Zeiten der Ordnung ε−1

folgen müssen, um in endlicher Zeit Veränderungen feststellen zu können. Wir wechselndaher zu makroskopischen Zeiteinheiten t = εs, indem wir (2.22) und (2.23) mit τ1(s) =s/ε und τ2(y) = y/ε transformieren, y = y/ε ersetzen und mit ε2 multiplizieren. InVerbindung mit Ψ ε(t, x, y) := Ψ(t/ε, x, y) erhalten wir die Gleichung:

i~ε ∂∂tΨ ε(t, x, y) = HεΨ ε(t, x, y) (2.24)

mit dem Hamiltonoperator

Hε = − ~2

2mε2∆x −~2

2m∆y + V (x, y). (2.25)

Das Trennen des Systems in langsame und schnelle Freiheitsgrade entspricht der Auftei-lung des Hilbertraums H = L2(Rd)⊗Hf = L2(Rd,Hf). Der Hilbertraum der langsamenFreiheitsgrade ist L2(Rd) und der der schnellen Freiheitsgrade Hf . Wir lassen die kon-krete Form von Hf außer Acht, da sie für uns nicht von Relevanz ist.

Bei der Dynamik von Atomteilchen, bei der sich die Kerne langsam im Vergleichzu den Elektronen bewegen, ist die Aufteilung in langsame und schnelle Freiheitsgradeoffensichtlich. Es gibt jedoch Fälle, in denen die entsprechende Form (2.25) des Hamil-tonoperators erst noch erarbeitet werden muss, vergleiche [HST01] und [TS02].

Den auf Hf wirkende Teil des Hamiltonoperators Hε bezeichnen wir mit

Hf(x) = − ~2

2m∆y + V (x, y). (2.26)

Hervorzuheben ist die Abhängigkeit von dem langsamen Parameter x ∈ Rd durch dasPotential V . Für jedes x ist der Hamiltonoperator Hf(x) selbstadjungiert mit Spektrumσ(Hf(x)). Der operatorwertige Multiplikationsoperator

Hf : Rd → L(DHf ,Hf) (2.27)x 7→ Hf(x)

definiert dann mit(Hfϕ) (x) := Hf(x)ϕ(x)

einen Hamiltonoperator, der auf L2(Rd,Hf) wirkt. Die Domäne DHf sei dabei ausge-stattet mit der Graphennorm von Hf(x0) mit dem Anfangswert x0.

18

Page 27: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

2.2. Raumadiabatische Störungstheorie

Die in Abschnitt 2.1.2 verwendete und für uns notwendige Lückenbedingung müssenwir nun auch auf die raumadiabatische Theorie übertragen. Diese vorauszusetzen, istfür uns unabdinglich und sowohl in dieser Arbeit als auch generell bei der Problem-lösung in der adiabatischen Störungstheorie ein wichtiger Bestandteil. Für Problemeohne Lückenbedingung sei auf den Ausblick in Kapitel 5 verwiesen.

Definition 2.5 (Lückenbedingung). Für x ∈ Rd sei die Lückenbedingung erfüllt, fallsσ∗(x) ⊂ σ(Hf(x)) und zwei Funktionen f± ∈ Cb(Rd,R) existieren, die ein IntervallI(x) = [f−(x), f+(x)] definieren mit

σ∗(x) ⊂ I(x) und infx∈Rd

dist(I(x), σ(Hf(x))\σ∗(x)

)=: g > 0. (2.28)

f+(x)

σ∗(x)

f−(x)

σ(Hf(x))

x

Abbildung 2.2: Das Spektrum σ∗(x) ist lokal durch eine Lücke (weiß) isoliert.

Auf ganzH definiert analog zu oben P∗ ∈ L(L2(Rd,Hf)

)mit (P∗ϕ) (x) := P∗(x)ϕ(x)

eine Projektion. Dabei sei P∗(x) die spektrale Projektion auf den zu σ∗(x) korrespondie-renden Eigenraum von Hf(x). Sowohl P∗(x) als auch P∗ sind als orthogonale Projektio-nen selbstadjungiert. Mit P⊥∗ (x) und P⊥∗ bezeichnen wir wiederum deren entsprechendeorthogonale Projektionen. Mit Blick auf die zeitadiabatische Theorie zu Beginn kom-mutieren auch Hf(x) und P∗(x),

[Hf(x), P∗(x)] = 0.

Obwohl P∗ im Allgemeinen keine spektrale Projektion von Hf sein muss, gilt dennoch

[Hf , P∗] =[P⊥∗ , Hf

]= 0.

Insbesondere ist P∗H unter der von Hf erzeugten Dynamik invariant, also ist auch[e−

i~Hf

tε , P∗

]= 0.

19

Page 28: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 2. Adiabatische Störungstheorie für die Schrödinger Gleichung

Für Hε und den unitären Propagator e−i/~Hεt/ε gilt dies allerdings nicht. Hier ist

‖[Hε, P∗]‖ =∥∥∥∥∥[− ~2

2mε2∆x, P∗

]∥∥∥∥∥ = O(ε). (2.29)

Um uns das klar zu machen, betrachten wir[~2

2mε2∆x, P∗

]= ~2

2mε2(∆xP∗ − P∗∆x)

= ~2

2mε2(∇x((∇xP∗) + P∗∇x

)− P∗∆x

)= ~2

2mε2((∆xP∗) + 2(∇xP∗) · ∇x + P∗∆x − P∗∆x

)= ~2

2mε2(∆xP∗) + ~2

mε(∇xP∗) · ε∇x.

Wir müssen beachten, dass ‖ε∇xΨ ε‖2 ε-unabhängig beschränkt ist, anstatt von derOrdnung O(ε) zu sein. Dies ist der Fall, da wir nur Situationen betrachten, in denen dieGesamtenergie ε-unabhängig beschränkt ist, also insbesondere die kinetische Energieder Wellenfunktion Ψ ε:

Ekin(Ψ ε) =⟨Ψ ε,− ~2

2mε2∆xΨ

ε⟩ = ~2

2m‖ε∇xΨε‖2 = O(1).

Hierauf werden wir im nachfolgenden Abschnitt 2.2.2 noch einmal genauer eingehen.A priori ist es nicht klar, dass der spektrale Unterraum P∗H für ausreichend lange

Zeiten approximativ invariant unter der von Hε erzeugten Zeitentwicklung ist. Umaber nicht triviale Bewegungen der Freiheitsgrade analysieren zu können, also zumBeispiel die der Kerne, müssen wir deren Dynamik über lange Zeiten verfolgen, da dieGesamtenergie beschränkt ist und sich die Kerne infolgedessen langsam bewegen. Nunist zwar der spektrale Unterraum P∗H zumindest für mikroskopische Zeiten s = t/εund kleine ε approximativ klein, jedoch nicht für endliche makroskopische Zeiten t, dieuns primär interessieren, ∥∥∥[e− i

~Hε t

ε , P∗]∥∥∥ = O(|t|). (2.30)

Dieser Umstand ist vergleichbar mit (2.18). Wir starten zunächst mit einer Gleichungder Ordnung O(1). Mittels analoger Argumentation wie in dem zeitadiabatischen Fallverbessern wir (2.30) zu

∥∥∥[e− i~H

ε tε , P∗

]∥∥∥ = O(ε (1 + |t|)).

Voraussetzung ist auch hier, dass σ∗(x) durch eine Lücke vom restlichen Teil des Spek-trums σ(Hf(x)) getrennt ist, das heißt die Lückenbedingung aus Definition 2.5 musserfüllt sein. Dies ist eines der Resultate der raumadiabatischen Theorie und steht inVerbindung mit dem nachfolgenden raumadiabatischen Theorem.

20

Page 29: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

2.2. Raumadiabatische Störungstheorie

2.2.2 Das raumadiabatische Theorem

Aufgrund der Lücke im Spektrum von Hf erwarten wir eine Entkopplung des Unter-raums P∗H von seinem orthogonalen Kompliment für kleine ε respektive für langsameFreiheitsgrade. Deshalb führen wir den diagonalen Hamiltonoperator

Hεdiag = P∗H

εP∗ + P⊥∗ HεP⊥∗

ein. Nach Konstruktion ist Hεdiag selbstadjungiert, kommutiert mit P∗ und es gilt∥∥∥[e− i

~Hεdiag

tε , P∗

]∥∥∥ = 0.

Dass es Sinn macht den Operator Hεdiag anstelle von Hε zu verwenden, werden wir im

Folgenden zeigen. Hierfür muss die von Hεdiag erzeugte Dynamik für Startwerte in P∗H

für ausreichend lange Zeiten nahe der von Hε erzeugten Dynamik liegen. Das wiederumbedeutet, folgender Term muss klein sein:(

e−i~H

ε tε − e−

i~H

εdiag

)P∗

= −e−i~H

ε tε

(e

i~H

ε tε e−

i~H

εdiag

tε − Id

)P∗

= −e−i~H

ε tε

[e

i~H

εse−i~H

εdiags

] tε

0P∗

= − i~

e−i~H

ε tε

∫ tε

0e

i~H

εsHεe−i~H

εdiags − e

i~H

εsHεdiage−

i~H

εdiags ds P∗

= i~

e−i~H

ε tε

∫ tε

0e

i~H

εs(Hεdiag −Hε)e−

i~H

εdiags ds P∗

= i~

e−i~H

ε tε

∫ tε

0e

i~H

εs(Hεdiag −Hε)P∗e−

i~H

εdiags ds

= i~

e−i~H

ε tε

∫ tε

0e

i~H

εs[P∗, Hε]P∗e−i~H

εdiags ds . (2.31)

Mit Verweis auf (2.29) ist der Kommutator in (2.31) von der Ordnung ε. Wegen derwachsenden oberen Integrationsgrenze wird der ordnungsmäßig kleine Integrand ausge-glichen und somit ist auf den ersten Blick der ganze Term von der Ordnung 1. Es istalso zunächst keineswegs klar, dass die von Hε

diag erzeugte Dynamik auch nach längerenZeiten der tatsächlichen Dynamik in etwa entspricht. Erneut kommt die Argumentationzum Tragen, welche wir bereits im Beweis des zeitadiabatischen Theorems 2.3 verwen-det haben: Der Integrand oszilliert und dadurch vergrößert sich der Fehler nicht, trotzdes sich verändernden Integrationsbereichs. Dies wurde in [ST01] detailliert gezeigt undin [Teu03] überarbeitet dargestellt. Als Resultat geht das raumadiabatische Theoremhervor.

Die adiabatische Entkopplung basiert auf einer Trennung der (Zeit-)Skalen. Bei un-serem Vorhaben, das System in schnelle und langsame Freiheitsgrade zu unterteilen,müssen wir daher noch mehr Details berücksichtigen. Momentan beruht diese Tren-nung beispielsweise bei Molekülen darauf, dass sich die Kerne langsam im Vergleichzu den Elektronen bewegen. Im Allgemeinen kann dies allerdings nicht vorausgesetzt

21

Page 30: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 2. Adiabatische Störungstheorie für die Schrödinger Gleichung

werden. Es trifft nur zu, wenn deren kinetische Energien ungefähr von der gleichenGrößenordnung sind. Ist dies nicht der Fall, können die Kerne trotz einer viel größerenMasse, also unabhängig davon wie klein ε ist, hohe Geschwindigkeiten erreichen, dieim Vergleich zu denen der Elektronen nicht mehr vernachlässigbar sind. Die adiabati-sche Theorie fällt dadurch in sich zusammen. Um dies zu verhindern, betrachten wirnur Fälle mit endlichen Gesamtenergien, indem wir E(Hε) := χ(−∞,E](Hε) als Projek-tion auf Gesamtenergien kleiner als E definieren. Hier bezeichne χ die charakteristischeFunktion.

Mit diesen Voraussetzungen sind wir nun in der Lage, das raumadiabatischeTheorem konkret zu formulieren:

Theorem 2.6 (Teufel 2003). Sei Hε der in (2.25) definierte Hamiltonoperator. Derin (2.26) und (2.27) definierte selbstadjungierte Operator Hf ∈ C∞b (Rd,L(DHf ,Hf))mit einer geeigneten dicht definierten Domäne DHf , wobei DHf ausgestattet sei mit derGraphennorm von Hf(x0), genüge der Lückenbedingung aus Definition 2.5. Dann istHε

diag selbstadjungiert auf der Domäne von Hε und es existiert eine von ε unabhängigeKonstante C <∞, sodass für alle t ∈ R gilt∥∥∥(e−

i~H

ε tε − e−

i~H

εdiag

)E(Hε)

∥∥∥L(H)

≤ εC (1 + |E|) (1 + |t|) . (2.32)

Diese Aussage lässt sich analog zu dem zeitadiabatischen Fall deuten. Das Systembleibt unter der adiabatischen Entwicklung von Hε

diag approximativ, das heißt bis aufeinen Fehler der Ordnung ε, in dem Unterraum P∗H. Es sei angemerkt, dass (2.32)ebenso für E

(− ~2

2mε2∆xΨ

ε)gilt.

Bemerkung. Essenzielle Grundlage für die gesamte Theorie und deren Anwendung in-nerhalb dieser Arbeit ist die Forderung nach der Lückenbedingung im Spektrum. Wasist jedoch, wenn diese Forderung nicht erfüllt und keine Lücke vorhanden ist? Mit dieserFrage haben sich unter anderem Bornemann [Bor98], Avron und Elgart [AE99] sowieTeufel [Teu03] beschäftigt. Unlängst ist man sogar noch davon ausgegangen, dass dieLückenbedingung in den adiabatischen Theoremen unumgänglich und unersetzlich ist.Dies stellte sich jedoch zum Teil als falsch heraus.

2.2.3 Effektive Dynamik

Sei es aus physikalischen Gründen oder aufgrund der Wahl des Startwertes, in vie-len Situationen sind wir in erster Linie nur an der von PHε

diagP erzeugten Dynamikin P∗H interessiert. Wie bereits am Ende von Abschnitt 2.1 erwähnt, ist es möglich,einen effektiven Hamiltonoperator abzuleiten, der eine einfachere Form hat als der desGesamtsystems und besagte Dynamik beschreibt.

Die Eigenfunktionen φj(x, ·) ∈ Hf zu den Eigenwerten Ej(x) von Hf(x) lösen dasEigenwertproblem

Hf(x)φj(x, ·) = Ej(x)φj(x, ·). (2.33)

Liegt unser Startwert Ψ ε0 = ψε0φj zu Beginn in dem spektralen Unterraum P∗H, danngilt nach dem raumadiabatischen Theorem 2.6 approximativ für spätere Zeitpunkte

Ψ ε(t, x, y) =(e−

i~H

ε tεΨ ε0

)(x, y) ≈

(e−

i~H

εdiag

tεΨ ε0

)(x, y) =: ψε(t, x)φj(x, y). (2.34)

22

Page 31: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

2.2. Raumadiabatische Störungstheorie

Die zeitabhängige Funktion ψε(t) ∈ L2(Rd) ist dabei die Lösung der effektiven Schrö-dingergleichung

i~ε ∂∂tψε(t, x) = Hε

effψε(t, x), ψε(t0) = ψε0 ∈ L2(Rd). (2.35)

In anderen Worten heißt das also, löst ψε die Gleichung der langsamen Freiheitsgra-de (2.35), dann erhalten wir die Lösung der Schrödingergleichung des Gesamtsystems(2.24) durch Multiplizieren von ψε mit dem zugehörigen Eigenzustand φj der schnellenFreiheitsgrade aus (2.33).

Der effektive Hamiltonoperator Hεeff wirkt auf dem Konfigurationsraum L2(Rd) der

langsamen Freiheitsgrade. Er kann mithilfe des Identifikationsoperators

U : P∗H → L2(Rd) (2.36)ψεφj 7→ ψε,

welcher den relevanten spektralen Unterraum P∗H, in dem die quantenmechanischeDynamik stattfindet, und den Hilbertraum L2(Rd) verknüpft, angegeben werden:

Hεeff = UP∗Hε

diagU∗ = UP∗HεP∗U∗. (2.37)

Die effektive Dynamik innerhalb des zum Eigenwert Ej(x) korrespondierenden Unter-raums P∗H wird also durch den Hamiltonoperator Hε

eff beschrieben und ψε(t, x) kannausgedrückt werden durch

ψε(t, x) = e−i~H

εeff

tεψε(0, x).

Wie wir später sehen werden, gehen die schnellen Freiheitsgrade nur über die EigenwerteEj(x) und ein zusätzliches Potential, welches ausschließlich Ableitungen der zugehörigenEigenfunktionen enthält, in Hε

eff ein. Dadurch bedingt ist seine Form in vielen Fällenum einiges simpler als die des Hamiltonoperators des Gesamtsystems. Wenn uns die Ei-genwerte bekannt sind, können wir die effektive Schrödingergleichung lösen und damitdann schließlich eine Aussage über die effektive Dynamik im Gesamtsystem ableiten.Diese Vorgehensweise ermöglicht uns die Analyse von größeren und eventuell kompli-zierteren Systemen. In Kapitel 3 werden wir Hε

eff explizit berechnen und in Kapitel 4diesen Lösungsansatz schließlich für ein zweidimensionales Beispielmodell numerischuntersuchen. Wir werden bemerken, dass der Unterschied bezüglich der Laufzeit dernumerischen Berechnungen immens ist.

23

Page 32: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls
Page 33: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

3 Die zeitabhängige Born-Oppenheimer-Approximation

Die Born-Oppenheimer-Approximation bezeichnet eine Näherung, mittels der man dieSchrödingergleichung von Systemen aus mehreren Teilchen vereinfachen kann. In einemSystem mit schweren und leichten Teilchen ändern diese ihre Bewegungseigenschaften,wie zum Beispiel Richtung und Geschwindigkeit auf unterschiedlichen Zeitskalen. DieGleichung, die die Dynamik der schweren, langsamen Teilchen beschreibt, kann daherohne Berücksichtigung der Bewegung der leichten, schnellen approximativ gelöst werdenund umgekehrt.

Der Grundstein zur Born-Oppenheimer-Approximation wurde 1927 von Max Bornund J. Robert Oppenheimer mit [BO27] gelegt. Diese erkannten, dass in Molekülendas Verhältnis zwischen der elektronischen und der nuklearen Masse me/mn klein istund daher als Parameter für adiabatische Berechnungen genutzt werden kann. Ver-einfacht gesagt kann die Dynamik der Kerne aus der Perspektive der Elektronen ver-nachlässigt werden, da sie aus deren Sicht praktisch still stehen. Demgegenüber wirddie Bewegung der Kerne von den leichten Elektronen nahezu nicht beeinflusst. DurchTrennung der Freiheitsgrade ist eine Zerlegung der molekularen Schrödingergleichungin eine Gleichung für die Elektronen und eine für die Kerne möglich. Diese sind getrenntvoneinander nun einfacher zu lösen.

Mit [PST07] als Grundlagen wollen wir in diesem Kapitel die effektive Dynamiknäher untersuchen, indem wir die konventionelle zeitabhängige Born-Oppenheimer-Ap-proximation betrachten und den Born-Oppenheimer-Hamiltonoperator erster Ordnungherleiten. Darüber hinaus geben wir ein systematisches Schema für die Korrektur inhöherer Ordnungen an und präsentieren eine elementare Herleitung des zeitabhängigenBorn-Oppenheimer-Hamiltonoperators zweiter Ordnung.

3.1 Die konventionelle Born-Oppenheimer-Approximation

3.1.1 Anwendung auf molekulare Dynamik

Zu Beginn dieses Kapitels wollen wir die bisherige Theorie auf ein molekulares Modellmit Kernen und Elektronen anwenden. Wir konkretisieren die bisher verwendete Massem durch die Kernmasse mn und die Elektronenmasse me und geben den molekularen

Page 34: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 3. Die zeitabhängige Born-Oppenheimer-Approximation

Hamiltonoperator wie folgt an

Hmol = − ~2

2mn∆x −

~2

2me∆y + V (x, y).

Die Positionen der K Kerne und l Elektronen sei dabei gegeben durch x = {x1, . . . , xK}und y = {y1, . . . , yl}. Das System entwickelt sich also in Rd+q = R3(K+l). Mit V (x, y)bezeichnen wir das Wechselwirkungspotential. Intuitiv ist dies das Coulomb-Potentialfür Elektronen und Kerne. Damit unsere Ergebnisse mathematisch korrekt sind, gehenwir jedoch von dem physikalischen Modell aus, dass Kerne keine Punktladungen sind,sondern ihre Ladungsverteilungen „verschmiert“ sind. Um die Notation zu vereinfachen,haben die Kerne alle die gleiche Masse mn und es sei ~ = 1 = me. Der adiabatischenParameter sei nun definiert durch

ε =√memn

.

Auf diese Weise erhalten wir den Hamiltonoperator

Hε = −12ε

2∆x +Hf(x) (3.1)

mit dem Hamiltonoperator der schnellen Freiheitsgrade beziehungsweise im weiterenVerlauf auch elektronischer Hamiltonoperator genannt

Hf(x) = −12∆y + V (x, y).

Hat das gegebene Molekül eine sehr komplexe Struktur, ist es oft nötig, nur einen vor-sichtig ausgewählten Teil der Freiheitsgrade zu betrachten. Im Folgenden nehmen wiran, dass der Hamiltonoperator schlussendlich von der Form (3.1) ist. Dies geschieht ge-gebenenfalls auf Kosten einer entsprechenden Anpassung von Hf(x). Ebenso setzen wirwie in Kapitel 2.2 voraus, dass die kinetische Energie der anfänglichen Wellenfunktionunabhängig von ε beschränkt ist und die Kerne sich folglich langsam im Vergleich zuden Elektronen bewegen.

Hinter der elektronischen Schrödingergleichung, genauer gesagt dem Eigenwertpro-blem der Elektronen (2.33), steckt eine teils aufwändige und nicht triviale Lösungstheo-rie. Für unsere Zwecke sehen wir daher das Problem bereits als gelöst und die Eigenwertesowie die zugehörigen Eigenvektoren als gegeben an. Die Eigenvektoren aus (2.33) kön-nen wir so wählen, dass sie eine Orthonormalbasis bilden, da Hf(x) selbstadjungiert ist.Sie sind also bezüglich des Skalarprodukts in Hf durch

〈φj(x, ·), φj′(x, ·)〉Hf= δjj′ (3.2)

normiert, jedoch nur bis auf eine Phase θj(x) eindeutig bestimmt. Die Eigenwerte iden-tifizieren wir einschließlich Vielfachheit durch

E1(x) ≤ E2(x) ≤ . . . .

Der Graph von Ej wird als j-te Born-Oppenheimer-Energieoberfläche bezeichnet. Nach(2.34) sind Zustände, deren Elektronen sich im j-te Eigenzustand befinden, von derForm

ψε(t, x)φj(x, y). (3.3)

26

Page 35: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

3.1. Die konventionelle Born-Oppenheimer-Approximation

Nach unseren Vorbereitungen in Kapitel 2.2.1 wissen wir, dass wir (3.3) entweder alseine Wellenfunktion in dem gesamten Hilbertraum H = L2(Rd)⊗Hf = L2(Rd,Hf) odernur als eine Wellenfunktion für die Elektronen, welche von x ∈ Rd abhängt, auffassenkönnen. Im zweiten Fall wird x jeweils fest gehalten und wir führen die abkürzende No-tation ψε(t, x)φj(x) ∈ Hf ein. Auf diese Weise kann, zusammen mit der physikalischenBra-Ket-Notation, die Projektion auf Zustände der Form (3.3) angegeben werden durch

Pj(x) = |φj(x)〉 〈φj(x)| .

Aus Darstellungsgründen verzichten wir ab sofort auf den zuvor stets verwendeten Sternim Index, behalten allerdings den Bezug zu dem isolierten Teil des Spektrums im Hin-terkopf. In Abschnitt 2.2.1 haben wir bereits erwähnt, dass es sich dabei um orthogonaleProjektionen handelt. Hier sehen wir, dass dies darauf zurückzuführen ist, dass die φj(x)orthonormal sind. Im Hinblick auf die Lückenbedingung (2.28) ist

P (x) =∑j∈I|φj(x)〉 〈φj(x)|

die Projektion auf den relevanten spektralen Unterraum von H. Benachbarte Born-Oppenheimer-Energieoberflächen werden dabei durch die Indexmenge I beschrieben.Eine vom stetigen Spektrum und anderen Energieoberflächen isolierende Lücke setzenwir voraus. Im weiteren Verlauf dieser Arbeit nehmen wir den speziellen aber wich-tigsten Fall an, dass I nur aus einem Element besteht, ergo I = {j}. Wir gehen also,wie in Abbildung 3.1 dargestellt, von einer einzelnen isolierten Energieoberfläche aus.Entsprechend ergibt sich demnach für uns

P (x) = Pj(x) = |φj(x)〉 〈φj(x)| .

Ei(x)f+(x)Ej(x)f−(x)

σ(Hf(x))

x

Abbildung 3.1: Einzelne isolierte Energieoberfläche Ej(x), wobei j ∈ I, i /∈ I.

27

Page 36: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 3. Die zeitabhängige Born-Oppenheimer-Approximation

3.1.2 Der effektive Hamiltonoperator

Wir wollen nun den in Abschnitt 2.2.3 angesprochenen effektiven HamiltonoperatorHεeff

konkret berechnen. Zu diesem Zweck ermitteln wir

(Hεψεφj) (t, x, y) =(−1

2ε2∆x +Hf(x)

)ψε(t, x)φj(x, y)

= −12ε

2∇x((∇xψε(t, x)

)φj(x, y) + ψε(t, x)

(∇xφj(x, y)

))+Hf(x)ψε(t, x)φj(x, y)

= −12ε

2((

∆xψε(t, x)

)φj(x, y) + 2

(∇xψε(t, x)

)(∇xφj(x, y)

)+ ψε(t, x)

(∆xφj(x, y)

))+ Ej(x)ψε(t, x)φj(x, y)

=(−1

2ε2∆x + Ej(x)

)ψε(t, x)φj(x, y)

− ε2(∇xψε(t, x))(∇xφj(x, y)

)− 1

2ε2ψε(t, x)

(∆xφj(x, y)

).

Um die Dynamik in dem spektralen Unterraum zu erhalten, projizieren wir darauf unterVerwendung des Skalarprodukts auf Hf :

〈φj(x), Hεψε(t, x)φj(x)〉Hf

=(−1

2ε2∆x + Ej(x)

)ψε(t, x)

− ε2(∇xψε(t, x))〈φj(x),∇xφj(x)〉Hf

− 12ε

2ψε(t, x)〈φj(x),∆xφj(x)〉Hf

=(1

2(−iε∇x − εA(x)

)2+ Ej(x) + 1

2ε2VBH(x)

)ψε(t, x). (3.4)

In dem letzten Term wird A(x) Koeffizient des Berry-Zusammenhangs genannt undangegeben durch

A(x) = i〈φj(x),∇xφj(x)〉Hf. (3.5)

Das zusätzliche Potential VBH wird in diesem Kontext als Born-Huang-Potential be-zeichnet und lautet

VBH(x) = 〈∇xφj(x),∇xφj(x)〉Hf−A(x)2

= 〈∇xφj(x),∇xφj(x)〉Hf− 〈∇xφj(x), φj(x)〉Hf

〈φj(x),∇xφj(x)〉Hf

= 〈∇xφj(x),∇xφj(x)〉Hf− 〈∇xφj(x), P (x)∇xφj(x)〉Hf

= 〈∇xφj(x), P⊥(x)∇xφj(x)〉Hf. (3.6)

Wie bereits oben angemerkt, sind die φj(x) nur bis auf eine Phase θj(x) eindeutigbestimmt. Durch eine Transformation φj(x) 7→ φj(x) = eiθj(x)φj(x) sind wir in derLage, reelle Eigenvektoren zu wählen und wir erhalten Im〈φj(x),∇xφj(x)〉 = 0. Mit

28

Page 37: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

3.2. Korrekturen für die Approximation höherer Ordnungen

(3.2) gilt zusätzlich

Re〈φj(x),∇xφj(x)〉Hf= 1

2(〈φj(x),∇xφj(x)〉Hf

+ 〈∇xφj(x), φj(x)〉Hf

)= 1

2∇x〈φj(x), φj(x)〉Hf

= 0.

Bei reellen Eigenvektoren gilt für das Born-Huang-Potential also

VBH(x) = 〈∇xφj(x),∇xφj(x)〉Hf. (3.7)

Den Berry-Zusammenhangskoeffizienten (3.5) können wir dadurch ebenso wegeichen,das heißt, wir setzen A(x) = 0. Mit dem Impulsoperator p = −iε∇x haben wir deneffektiven Hamiltonoperator für diesen Fall schlussendlich definiert durch

Hεeff(x) := −1

2ε2∆x + Ej(x) + 1

2ε2VBH(x)

= 12p

2 + Ej(x) + 12ε

2VBH(x). (3.8)

Er wird teilweise auch als der konventionelle Born-Oppenheimer-Hamiltonoperator ers-ter Ordnung bezeichnet. Die elektronischen Freiheitsgrade gehen in ihn nur über dieeffektive potentielle Energie, genauer gesagt den Eigenwert Ej(x) und das Born-Huang-Potential VBH(x) ein.

Unter Zuhilfenahme des Identifikationsoperators (2.36) kann nun, bis auf einen Feh-ler der Ordnung ε, auf die Dynamik des Gesamtsystems geschlossen werden. Zielt manauf eine höhere Genauigkeit der Kerndynamik ab, stößt obiger konventioneller Ansatzallerdings an seine Grenzen. Mithilfe von Korrekturtermen kann jedoch der konventio-nelle Born-Oppenheimer-Hamiltonoperator weiter verbessert werden, um die Approxi-mation höherer Ordnungen, das heißt Fehler kleiner der Ordnung ε, zu ermöglichen. Imnachfolgenden Abschnitt werden wir die allgemeine Idee dafür grob erläutern und andem Beispiel der zweiten Ordnung demonstrieren.

3.2 Korrekturen für die Approximationhöherer Ordnungen

Die Theorie zur korrekten Approximation höherer Ordnungen wurde in zahlreichenPublikationen entwickelt. Dazu gehören unter anderem [MS02,NS04, PST03]. Allesamtbedienen sie sich der Pseudodifferentialrechnung. Im Gegensatz dazu orientieren wir unsan der Vorgehensweise in [PST07], indem wir im Folgenden keine Pseudodifferentialope-ratoren verwenden, sondern auf elementarere Mittel zurückgreifen. Dabei beschränkenwir uns auf das Wesentliche und deuten das allgemeine Schema ohne die nötigen tech-nischen Details nur an. Ein größeres Augenmerk legen wir auf die explizite Herleitungdes Born-Oppenheimer-Hamiltonoperators zweiter Ordnung.

Für eine Approximation höherer Ordnung ist die Projektion auf Eigenzustände derForm (3.3) nicht mehr geeignet. Ihr Bild ist nur bis auf Fehler der Ordnung ε invariantunter der von Hε erzeugten Dynamik. Deshalb müssen wir eine modifizierte Projektion

29

Page 38: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 3. Die zeitabhängige Born-Oppenheimer-Approximation

P ε finden, sodass der sogenannte superadiabatische Unterraum P εH auch für kleinereFehler invariant bleibt. Als ein Ergebnis in [PST03] hat es sich herausgestellt, dass diesmöglich ist und es solch eine orthogonale Projektion P ε gibt, die sich nicht zu sehrvon der ursprünglichen unterscheidet und für die es für jedes n ∈ N eine von ε und tunabhängige Konstante Cn <∞ gibt mit∥∥∥(e−iHε t

ε − e−iP εHεP ε tε

)P εE(Hε)

∥∥∥L(H)

≤ Cn εn |t|. (3.9)

An dieser Stelle haben wir die Zusammensetzung der effektive Dynamik aus (2.37)adaptiert. Unser Ziel ist es P εHεP ε nach Potenzen von ε zu entwickeln:

P εHεP ε =n∑i=0

εiHi +O(εn+1

)=: Hε

(n) +O(εn+1

).

Aus (3.9) wird folglich∥∥∥∥(e−iHε tε − e−iHε

(n)tε

)P εE(Hε)

∥∥∥∥L(H)

≤ Cn εn |t|. (3.10)

Ebenso kann P ε als Taylorreihe entwickelt und abhängig von der Ordnung an der ent-sprechenden Stelle abgeschnitten werden. Wir könnten nunHε

(n) als Born-Oppenheimer-Hamiltonoperator n-ter Ordnung bezeichnen, allerdings wirkt Hε

(n) auf einem ε-unab-hängigen Unterraum von H. Ein Nutzen der Standard-Born-Oppenheimer-Approxima-tion besteht jedoch auch darin, dass der effektive Hamiltonoperator auf Wellenfunktio-nen wirkt, die nur von nuklearen und, wenn überhaupt, von einigen wenigen diskretenelektronischen Freiheitsgraden abhängen.

In der Zwischenzeit wollen wir uns aber die Konstruktion der Projektion P ε überle-gen. Dabei folgen wir weiterhin dem Schema von Panati, Spohn und Teufel [PST07], diesich wiederum unter anderem an den in [NS04,MS02] entwickelten Methoden orientierthaben. Wir versuchen, die Koeffizienten Pj in der Reihendarstellung P ε(n) =

∑ni=0 ε

iPjzu bestimmen. Approximativ soll P ε(n) dabei eine Projektion sein und mit Hε bis zuTerme der Ordnung εn+1 kommutieren. Offensichtlich ist dies für P0 = P (x) mit Bezugauf (2.29) erfüllt. Hier und im Folgenden verzichten wir größtenteils auf die Darstellungder Abhängigkeiten von x, um die Übersichtlichkeit und Leserlichkeit der Gleichungenzu fördern. Induktiv können wir nun Pj ermitteln, indem wir∥∥∥∥(P ε(n)

)2− P ε(n)

∥∥∥∥ = O(εn+1

)(3.11)

sowie ∥∥∥[P ε(n)), Hε]∥∥∥ = O

(εn+1

)(3.12)

fordern und annehmen, dass selbiges schon für P ε(n−1) gilt. Wie zu Beginn dieses Ab-schnitts erwähnt, nutzt die allgemeine Konstruktion an dieser Stelle in der Regel ε-Dif-ferentialoperatoren. Demgegenüber leiten wir P1 elementar her. Wir berechnen also(

P ε(1)

)2− P ε(1) = (P0 + εP1)2 − (P0 + εP1)

= P0 + εP0P1 + εP1P0 + ε2P 21 − P0 − εP1

= ε(P0P1 + P1P0 − P1) +O(ε2). (3.13)

30

Page 39: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

3.2. Korrekturen für die Approximation höherer Ordnungen

Damit der Term der Ordnung ε in (3.13) verschwindet, setzen wir

P1 = P0P1 + P1P0 +O(ε). (3.14)

Die gleiche Prozedur führen wir durch für[P ε(1), H

ε]

= [P0 + εP1, Hε]

= [P0, Hε] + [εP1, Hf ] +O(ε2)

= iε(∇xP0) · (−iε∇x) + ε[P1, Hf ] +O(ε2).

In diesem Fall muss P1 der Gleichung

[P1, Hf ] = i(∇xP0) · (iε∇x) +O(ε)

genügen. Wir multiplizieren diese von beiden Seiten mit P0 beziehungsweise P⊥0 . Fürdie linke Seite erhalten wir unter anderem

P0[P1, Hf ]P⊥0 = P0(P1Hf − EjP1)P⊥0 = P0P1(Hf − Ej)P⊥0

undP⊥0 [P1, Hf ]P0 = P⊥0 (P1Ej −HfP1)P0 = P⊥0 (Ej −Hf)P1P0.

Auf dem Bild von P⊥0 können wir (Hf − Ej) invertieren und folgern

P0P1P⊥0 = iP0(∇xP0)P⊥0 (Hf − Ej)−1 · (iε∇x) +O(ε)

sowieP⊥0 P1P0 = −i(Hf − Ej)−1P⊥0 (∇xP0)P0 · (iε∇x) +O(ε).

Die Blockdiagonalterme P0P1P0 und P⊥0 P1P⊥0 sind von der Ordnung ε, da

P0(∇xP0)P0 = P⊥0 (∇xP0)P⊥0 = 0.

Für P1 ergibt sich infolgedessen

P1 = P0P1P0 + P0P1P⊥0 + P⊥0 P1P0 + P⊥0 P1P

⊥0

= P0P1P⊥0 + P⊥0 P1P0 +O(ε)

= iP0(∇xP0)(Hf − Ej)−1P⊥0 · (iε∇x) + adj. +O(ε)= i(iε∇x) · P0(∇xP0)(Hf − Ej)−1P⊥0 + adj. +O(ε)

= (−iε∇x) ·(−iP0(∇xP0)(Hf − Ej)−1P⊥0

)+ adj. +O(ε),

wobei wir, um eine mögliche Lösung zu erhalten, den hintersten Term der OrdnungO(ε) weglassen können. Diese Lösung erfüllt Bedingung (3.14) und wir bezeichnen sieebenso mit P1:

P1 = (−iε∇x) ·(−iP0(∇xP0)(Hf − Ej)−1P⊥0

)+ adj.

=: p ·B +B∗ · p. (3.15)

31

Page 40: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 3. Die zeitabhängige Born-Oppenheimer-Approximation

Die Glattheit von Hf(x) und die Lückenbedingung machen es uns möglich, in obigerRechnung auszunutzen, dass die spektrale Projektion P0(x) und die reduzierte Resol-vente

(Hf(x)−Ej(x)

)−1P0(x)⊥ glatte von x abhängige Funktionen mit Werten in dem

Raum der beschränkten Operatoren über dem elektronischen Hilbertraum Hf sind. In(2.20) haben wir dies bereits für den zeitadiabatischen Fall angedeutet gehabt. Auch imFolgenden werden wir uns besagte Eigenschaft häufiger zu Nutze machen, ebenso dieFolgerung, dass das Kommutieren des Impulsoperators p = −iε∇x mit solchen Opera-toren nur zu einem zusätzlichen Term der Ordnung ε führt, da hierdurch ε frei wirdund sich nicht wegen ‖∇xψε‖ = O(1/ε) vorher heraus kürzt. Des Weiteren sei an dieserStelle angemerkt, dass P1 aufgrund des Impulsoperators nicht von x abhängt.

Nachdem wir jetzt P ε(1) so konstruiert haben, dass die Bedingungen (3.11) und (3.12)erfüllt sind, lässt sich in einem weiteren Schritt ein Term der Ordnung ε2 so hinzufügen,dass die Projektion orthogonal wird. Diese bezeichnen wir ebenfalls mit P ε(1). Für mehrInformationen sei an dieser Stelle auf [NS04] verwiesen. Eine bisher offen gebliebeneFrage ist, ob das Bild von P ε(1) bis zu Fehlern der Ordnung O(ε2) invariant unter dervon Hε erzeugten Dynamik ist. Wieder einmal würde ein naiver Ansatz wie in (2.31)nur zu einem Fehler der Ordnung ε führen, welcher sich aber mit der raumadiabatischenVorgehensweise aus [ST01] beziehungsweise [Teu03] verbessern ließe. Stattdessen ziehenwir direkt [MS02, PST03] heran und mit zwei Konstanten C1 und C2 erhalten wiranstelle der Abschätzung des raumadiabatischen Theorems 2.6∥∥∥∥(e−iHε t

ε − e−iP ε(1)H

εP ε(1)

)P ε(1)E(Hε)

∥∥∥∥L(H)

≤ ε2 (C1 + |t|C2) . (3.16)

Mit dieser Abschätzung ermitteln wir nun einen effektiven Hamiltonoperator auf demBild von P ε(1), indem wir P ε(1)H

εP ε(1) nach ε-Potenzen bis zu Terme der Ordnung ε2 entwi-ckeln. Die konventionelle Born-Oppenheimer-Approximation wird auf diese Weise durcheine weitere Ordnung verbessert. Das Bild von P ε(1) wird indessen nicht von Wellenfunk-tionen der Form ψε(t, x)φεj(x, y) wie in (3.3) aufgespannt, weil der Impulsoperator inder Darstellung von P1 in (3.15) auftaucht. Folglich kann P ε(1)H

εP ε(1) nicht als effektiverHamiltonoperator, der ausschließlich auf Wellenfunktionen der Kerne wirkt, angesehenwerden. Genau das ist aber die entscheidende Eigenschaft, die die Born-Oppenheimer-Approximation überhaupt so nützlich macht. Um diese Charakteristik auch für höhe-re Ordnungen aufrecht zu erhalten, verknüpfen wir das Bild von P ε(1) mit dem Raumder nuklearen Wellenfunktionen L2(Rd). Auf diese Art und Weise sind wir bereits inAbschnitt 2.2.3 vorgegangen. Der unitäre Operator U ε : P εH → L2(Rd) definiert deneffektiven Born-Oppenheimer-Hamiltonoperator mit P ε = P0 + εP1 +O(ε2) durch

HεBO = U εP εHεP εU ε ∗.

Die Konstruktion von U ε wird in [PST03] spezifiziert. Aus (3.9) wird∥∥∥(e−iHε tε − U ε ∗e−iHε

BOtεU ε

)P εE(Hε)

∥∥∥L(H)

≤ Cn εn |t|. (3.17)

Der Ausdruck e−iHεBOt/ε bezeichnet den effektiven Born-Oppenheimer-Propagator der

Kerne innerhalb der isolierten Born-Oppenheimer-Energieoberfläche bis zu einer belie-bigen Ordnung. Den gewöhnlichen nullten und ersten Born-Oppenheimer-Hamiltonope-rator liefert die Reihenentwicklung von Hε

BO nach ε als führende Ordnungen. Zusätzlich

32

Page 41: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

3.2. Korrekturen für die Approximation höherer Ordnungen

sind aber auch die korrekten Terme höherer Ordnungen berücksichtigt. Später werdenwir auch den Born-Oppenheimer-Hamiltonoperator zweiter Ordnung mit dem SymbolHε

BO bezeichnen.An dieser Stelle sei noch einmal darauf hingewiesen, dass wir uns auf den Fall ei-

ner einzelnen isolierten Energieoberfläche Ej beschränkten. Auf die in [PST03] dar-gestellte allgemeine Konstruktion mit allen technischen Details verzichten wir. Umden Born-Oppenheimer-Hamiltonoperator zweiter Ordnung durch elementare Rech-nungen für unseren Spezialfall zu berechnen, nutzen wir nun unser ZwischenergebnisP ε = P0 + εP1 +O(ε2) mit P0(x) = |φj(x)〉 〈φj(x)| und zudem P1 = p ·B +B∗ · p aus(3.15). Analog zu P ε ist für U ε unser Ansatz U ε = U0 + εU1 +O(ε2), wobei

U0(x) = 〈φj(x)|

der unitäre Operator von P0H nach L2(Rd) aus (2.36) ist. Da die Wahl der Basis φj(x)eine Rolle spielt, ist U ε in dieser Form nicht eindeutig.

Im nächsten Schritt wollen wir U1 berechnen. Dazu nehmen wir ohne Beschränkungder Allgemeinheit an, dass U1 = U0A für einen Operator A auf H ist. Bis zur Ordnungε2 soll U ε ∗ unitär sein. Daher fordern wir

U εU ε ∗ =(U0 + εU1 +O(ε2)

) (U∗0 + εU∗1 +O(ε2)

)= U0U

∗0 + ε(U0U

∗1 + U1U

∗0 ) +O(ε2)

= Id + ε(U0A∗U∗0 + U0AU

∗0 ) +O(ε2) (3.18)

!= Id +O(ε2).

Des Weiteren soll das Bild von U ε ∗ dem Bild von P ε bis zu O(ε2) entsprechen:

P ε⊥U ε ∗ =(Id− P0 − εP1 +O(ε2)

) (U∗0 + εU∗1 +O(ε2)

)= U∗0 + εU∗1 − P0U

∗0 − εP0U

∗1 − εP1U

∗0 +O(ε2)

= U∗0 + εA∗U∗0 − P0U∗0 − εP0A

∗U∗0 − εP1U∗0 +O(ε2)

= (Id− P0 + εA∗ − εP0A∗ − εP1)U∗0 +O(ε2)

= P⊥0 U∗0 + ε (A∗ − P0A

∗ − P1)U∗0 +O(ε2)= ε (A∗ − P0A

∗ − P1)U∗0 +O(ε2) (3.19)!= O(ε2).

Der O(ε)-Term in (3.19) verschwindet für

(Id− P0)A∗ = P1P0 +O(ε), (3.20)

da dann

(A∗ − P0A∗ − P1)U∗0 = (P1P0 − P1 +O(ε))U∗0 = −P1P

⊥0 U

∗0 +O(ε) = O(ε).

Aus Gleichung (3.20) zusammen mit P1P0 = B∗ · p folgern wir A∗ = B∗ · p als einemögliche Lösung. Diese eliminiert zudem den Term der Ordnung O(ε) in (3.18). Wirerhalten

U1 = U0A = U0 p ·B. (3.21)

33

Page 42: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 3. Die zeitabhängige Born-Oppenheimer-Approximation

So definiert, ist der Operator U ε(1) = U0 + εU1 aber nur fast unitär und verknüpft P ε(1)Hmit L2(Rd) ebenso nur fast. Die Unitarität von U ε(1) ist für das Nachfolgende allerdingsessenziell. Mit einem zusätzlichen Term der Ordnung ε2, genauer ε2U ε2 mit U ε2 = O(1),können wir U ε(1) zu einem richtigen unitären Operator verbessern. Dieser verknüpftP ε(1)H mit L2(Rd) exakt. Den verbesserten Operator U ε(1) = U0 +εU1 +ε2U ε2 bezeichnenwir, wie wir es bereits bei P ε(1) gehandhabt haben, mit dem gleichen Symbol. Abgesehenvon der Unitarität von U ε(1) ist die exakte Form von U ε2 für uns nicht weiter relevant.Die genaue Konstruktion findet sich aber auch in [PST03]. Die Unitarität kombiniertmit

U1U∗0 + U0U

∗1 = U0 p ·BU∗0 + U0B

∗ · pU∗0 = U0P1U∗0 = 0

liefert uns

U ε(1)Uε ∗(1) =

(U0 + εU1 + ε2U ε2

) (U∗0 + εU∗1 + ε2U ε ∗2

)= Id + ε2 (U0U

ε ∗2 + U1U

∗1 + U ε2U

∗0 ) +O(ε3)

und bedeutet, dass sich die letzten beiden Summanden gegenseitig aufheben und derTerm zwischen den Klammern von der Ordnung O(ε) ist:

U0Uε ∗2 + U1U

∗1 + U ε2U

∗0 = O(ε). (3.22)

Gleichermaßen wissen wir

e−iP ε(1)H

εP ε(1)

tεP ε(1) = U ε ∗(1)U

ε(1)e−iP ε

(1)HεP ε

(1)tεU ε ∗(1)U

ε(1)P

ε(1)

= U ε ∗(1)e−iUε

(1)Pε(1)H

εP ε(1)U

ε ∗(1)

tεU ε(1)P

ε(1).

Bis zu Fehlern der Ordnung ε2 ist dies eine Approximation der tatsächlichen Zeitent-wicklung, denn zusammen mit (3.16) gilt∥∥∥∥(e−iHε t

ε − U ε ∗(1)e−iUε

(1)Pε(1)H

εP ε(1)U

ε ∗(1)

tεU ε(1)

)P ε(1)E(Hε)

∥∥∥∥L(H)

≤ ε2 (C1 + |t|C2) .

Um den Born-Oppenheimer-Hamiltonoperator zweiter Ordnung zu ermitteln, entwi-ckeln wir nun U ε(1)P

ε(1)H

εP ε(1)Uε ∗(1) nach Potenzen von ε. Im Gegensatz zu P ε(1)H

εP ε(1)wirkt dieser nun ausschließlich auf dem von ε unabhängigen nuklearen Wellenfunktio-nenraum L2(Rd) und hat die Form

U ε(1)Pε(1)H

εP ε(1)Uε ∗(1) = U ε(1)H

εU ε ∗(1) +O(ε3)

=(U0 + εU1 + ε2U ε2

)Hε

(U∗0 + εU∗1 + ε2U ε ∗2

)+O(ε3)

= U0HεU∗0 + εU0H

εU∗1 + ε2U0HεU ε ∗2

+ εU1HεU∗0 + ε2U1H

εU∗1 + ε2U ε2HεU∗0 +O(ε3)

= U0HεU∗0 (3.23)

+ ε (U0HεU∗1 + U1H

εU∗0 ) (3.24)+ ε2U1H

εU∗1 (3.25)+ ε2 (U0H

εU ε ∗2 + U ε2HεU∗0 ) +O(ε3). (3.26)

34

Page 43: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

3.2. Korrekturen für die Approximation höherer Ordnungen

Die vier Terme (3.23)-(3.26) werten wir einzeln aus, beginnend mit (3.23):

U0HεU∗0 = U0

(−1

2ε2∆x +Hf

)U∗0

= −12ε

2∆x + U0

[−1

2ε2∆x, U

∗0

]+ Ej

= p2

2 −12ε

2U0(∆xU

∗0 − U∗0 ∆x

)+ Ej

= p2

2 −12ε

2U0(∇x((∇xU∗0 ) + U∗0∇x

)− U∗0 ∆x

)+ Ej

= p2

2 −12ε

2U0((∆xU

∗0 ) + 2(∇xU∗0 ) · ∇x + U∗0 ∆x − U∗0 ∆x

)+ Ej

= p2

2 −12ε

2U0(∆xU∗0 )− ε2U0(∇xU∗0 ) · ∇x + Ej

= p2

2 −12ε

2〈φj ,∆xφj〉 − iε〈φj ,∇xφj〉 · p+ Ej

= p2

2 −12ε

2(iA · ∇x − i∇xA− 〈∇xφj ,∇xφj〉

)− εA · p+ Ej

= p2

2 −12εp · A −

12εA · p+ 1

2ε2〈∇xφj ,∇xφj〉+ Ej

= 12(p− εA)2 + Ej + 1

2ε2VBH

mit dem Berry-Zusammenhangskoeffizienten A(x) aus (3.5) und dem Born-Huang-Po-tential VBH(x) aus (3.6). Wie zu erwarten war, ist das gerade der konventionelle Born-Oppenheimer-Hamiltonoperator (3.8) respektive (3.4). Für einen korrekten Born-Op-penheimer-Hamiltonoperator zweiter Ordnung müssen wir beachten, dass nicht PH,sondern P ε(1) der adiabatisch invariante Unterraum bezüglich dieser Ordnung ist. Ergoist es also notwendig, zusätzlich die übrigen Terme (3.24)-(3.26) zu betrachten. Mit

(∇xB)P0 =(∇x

(−iP0(∇xP0)(Hf − Ej)−1P⊥0

))P0

= −iP0(∇xP0)(Hf − Ej)−1∇x(Id− P0)P0

= iP0(∇xP0)(Hf − Ej)−1(∇xP0)P0

und analog

P0∇xB∗ = −iP0(∇xP0)(Hf − Ej)−1(∇xP0)P0

35

Page 44: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 3. Die zeitabhängige Born-Oppenheimer-Approximation

erhalten wir für (3.24)

U0HεU∗1 + U1H

εU∗0 = U0HεB∗ · pU∗0 + U0 p ·BHεU∗0

= U0

(p2

2 B∗ · p+ p ·Bp

2

2

)U∗0 + U0

(HfB

∗ · p+ p ·BHf)U∗0

= U0

(p2

2 B∗ · p+ p ·Bp

2

2

)U∗0

= U0

(p2

2 B∗ · p+B∗

p2

2 · p+ p ·Bp2

2 − p3 ·B

)U∗0

= U0(p2 ·B∗ · p− p2 ·B · p+O(ε)

)U∗0

= −iεU0(p · ∇xB∗ · p− p · ∇xB · p

)U∗0 +O(ε2)

= −2εU0 p · P0(∇xP0)(Hf − Ej)−1(∇xP0)P0 · pU∗0 +O(ε2)=: −2εM(x, p) +O(ε2).

Auch hier haben wir ein weiteres Mal genutzt, dass das Kommutieren von p = −iε∇xmit x-abhängigen Operatoren nur zu einem zusätzlichen Term der Ordnung O(ε) führt.Deshalb wäre auch eine andere Reihenfolge der Operatoren möglich. Bei dieser Wahlerhalten wir

M(x, p) = U0 p · P0(∇xP0)(Hf − Ej)−1(∇xP0)P0 · pU∗0= p · U0P0(∇xP0)(Hf − Ej)−1(∇xP0)P0U0 · p+O(ε)= p · 〈φj | (∇xP0)(Hf − Ej)−1(∇xP0) |φj〉 · p+O(ε)= p · 〈φj | (∇xP0)(Id− P0)(Hf − Ej)−1(Id− P0)(∇xP0) |φj〉 · p+O(ε)= p · 〈∇xφj | (Hf − Ej)−1 |∇xφj〉 · p+O(ε)=⟨∇xφj · p, (Hf − Ej)−1P⊥0 p · ∇xφj

⟩Hf

+O(ε), (3.27)

indem wir∇xP0(x) = |∇xφj(x)〉 〈φj(x)|+ |φj(x)〉 〈∇xφj(x)|

verwenden. Ferner ist⟨∇xφj · p, (Hf − Ej)−1P⊥0 p · ∇xφj

⟩Hf

=d∑

i,k=1

⟨(∂xiφj)(−iε∂xi), (Hf − Ej)−1P⊥0 (−iε∂xk

)(∂xkφj)⟩Hf.

Eine naheliegende symmetrische Wahl fürM wäre also

(Mψε) (x) =d∑

i,k=1

((−iε∂xi

)vik(x)

(−iε∂xk

))ψε(x) (3.28)

mit der x-abhängigen d×d-Matrix v mit den Einträgen

vik(x) =⟨∂xiφj(x),

(Hf(x)− Ej(x)

)−1P⊥0 (x)∂xk

φj(x)⟩Hf. (3.29)

36

Page 45: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

3.2. Korrekturen für die Approximation höherer Ordnungen

Die Schreibweise ∂x meint dabei ∂∂x .

Mit (3.25) fahren wir weiter fort und ermitteln

U1HεU∗1 = U0 p ·B

(p2

2 +Hf

)B∗ · pU∗0

= U0 p ·B (Hf − Ej)B∗ · pU∗0 + U0 p ·B(p2

2 + Ej

)B∗ · pU∗0

= U0 p · P0(∇xP0)(Hf − Ej)−1(∇xP0)P0 · pU∗0

+ U0 p ·BB∗ · pU∗0

(p2

2 + Ej

)+O(ε)

=M(x, p) + U1U∗1

(p2

2 + Ej

)+O(ε). (3.30)

Der zweite Summand in (3.30) wird unter Verwendung von (3.22) durch den viertenund letzten zu berechnenden Term (3.26) eliminiert:

U0HεU ε ∗2 + U ε2H

εU∗0 = U0

(p2

2 + Ej

)U ε ∗2 + U ε2

(p2

2 + Ej

)U∗0

= (U0Uε ∗2 + U ε2U

∗0 )(p2

2 + Ej

)+O(ε)

= −U1U∗1

(p2

2 + Ej

)+O(ε)

Zu guter Letzt setzen wir nun alles zusammen und erhalten für eine einzelne isolierteEnergieoberfläche Ej als Born-Oppenheimer-Hamiltonoperator zweiter Ordnung

HεBO(x) = 1

2(p− εA(x)

)2 + Ej(x) + 12ε

2VBH(x)− ε2M(x, p). (3.31)

mit A(x), VBH(x) undM(x, p) definiert in (3.5), (3.6) und (3.27). Im Vergleich zu demkonventionellen Ausdruck enthält dieser den zusätzlichen Term ε2M(x, p). Mittels der x-abhängigen Matrix v wird durchM(x, p) eine geschwindigkeitsabhängige Komponentehinzugefügt.

Der Fehler der Approximation der tatsächlichen Zeitentwicklung durch die von demBorn-Oppenheimer-Hamiltonoperator zweiter Ordnung Hε

BO generierte Dynamik istjetzt von der Ordnung ε2:

∥∥∥(e−iHε tε − U ε ∗(1)e

−iHεBO

tεU ε(1)

)P ε(1)E(Hε)

∥∥∥L(H)

≤ ε2 (C1 + |t|C2) . (3.32)

37

Page 46: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 3. Die zeitabhängige Born-Oppenheimer-Approximation

Dies ergibt sich mittels (3.10), (3.17) und einer Dreiecksungleichung aus∥∥∥∥(e−iHε tε − e−iHε

(n)tε

)P εE(Hε)

∥∥∥∥L(H)

+∥∥∥∥(e−iHε

(n)tε − U ε ∗e−iHε

BOtεU ε

)P εE(Hε)

∥∥∥∥L(H)

+∥∥∥(U ε ∗e−iHε

BOtεU εP ε − U ε ∗(1)e

−iHεBO

tεU ε(1)P

ε(1)

)E(Hε)

∥∥∥L(H)

≤ Cn εn |t|+ C2 ε2 |t|+ C1 ε

2

≤ ε2 (C1 + |t|C2)≤ ε2C (1 + |t|) .

Der konstante Faktor besteht also aus einem zeitabhängigen und einem zeitunabhängi-gen Teil. Die Qualität des Hamiltonoperators Hε

BO beeinflusst den zeitabhängigen Teil,während die Wahl der Startwerte aus dem adiabatischen beziehungsweise superadiaba-tischen Unterraum sich auf den zeitunabhängigen Teil auswirkt.

Aus der Lösung der effektiven Born-Oppenheimer-Schrödingergleichung der Kerne

iε ∂∂tψε(t, x) = Hε

BOψε(t, x), ψε(t) ∈ L2(Rd)

können wir eine approximative Lösung Ψ ε(t) der molekularen Schrödingergleichung desGesamtsystems

iε ∂∂tΨ ε(t, x, y) = HεΨ ε(t, x, y), Ψ ε(t) ∈ H = L2(Rd,Hf)

erhalten, indem wir Ψ εBO(t) := U ε ∗(1)ψε(t) definieren. Nach (3.32) ist die Abweichung von

Ψ εBO(t) von der tatsächlichen Lösung Ψ ε(t) bei dem selben Startwert Ψ ε(0) = Ψ εBO(0) ∈P ε(1)E(Hε)H in der Norm von H von der Ordnung ε2. In vielen Fällen ist es jedoch nichtzwingend nötig, auf die Wellenfunktion bezüglich des Gesamtsystems zurück abzubil-den, da es ausreicht ψε(t) zu kennen. Für die Wahrscheinlichkeitsverteilung bezüglichder Positionen der Kerne gilt zum Beispiel⟨

Ψ εBO(t, x, ·), Ψ εBO(t, x, ·)⟩Hf

=⟨U ε ∗(1)ψ

ε(t, x), U ε ∗(1)ψε(t, x)

⟩Hf

= 〈ψε(t, x), ψε(t, x)〉Hf

=∫|ψε(t, x)|2 dy

= |ψε(t, x)|2 +O(ε2).

Abschließend sei angemerkt, dassHεBO invariant unter der Wahl der Phase von φj(x)

ist. Der Berry-Zusammenhangskoeffizient A(x) hängt von dieser Wahl ab und kann beiuns, wie gezeigt, durch eine geeignete Wahl der Phase weggeeicht werden. Die Termeder Ordnung ε2 hingegen sind zwar durch φj(x) definiert, aber dennoch invariant unterdieser Wahl. Klar wird das, indem wir VBH(x) undM(x, p) mittels

TrHf (A) =∑

i∈I∪Ic

〈φi(x)|A|φi(x)〉

38

Page 47: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

3.2. Korrekturen für die Approximation höherer Ordnungen

in Abhängigkeit von P0(x) ausdrücken:

VBH(x) = TrHf

(∇xP0(x) · ∇xP0(x)P⊥0 (x)

),

M(x, p) = TrHf

((p · ∇xP0(x)

)2(Hf(x)− Ej(x)

)−1P⊥0 (x)

).

39

Page 48: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls
Page 49: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

4 Numerische Simulation

Die mathematische Theorie der vorangegangenen Kapitel wollen wir nun auf ein konkre-tes Modell anwenden. Dazu betrachten wir jedoch keinen molekularen Hamiltonopera-tor, sondern ersetzen den elektronischen Hamiltonoperator durch einen harmonischenOszillator. Die Wellenfunktion breitet sich in zwei Dimensionen aus und wird durchein Potential eingeschränkt. Der effektive Hamiltonoperator lebt auf einem eindimen-sionalen Raum. Dieses Modell wurde bereits in [Böl09] thematisiert und an der dortverwendeten Vorgehensweise werden wir uns in diesem Kapitel orientieren. Ziel ist es,die jeweils bei der Born-Oppenheimer-Approximation mit und ohne Korrekturterm ent-standenen Fehler miteinander zu vergleichen. Außerdem werfen wir einen kurzen Blickauf die benötigten Rechenzeiten. Wir werden sehen, dass eine beachtliche qualitati-ve Verbesserung durch die Korrektur gelingt, deren zeitlicher Mehraufwand sich aberdennoch im Rahmen hält. Die gesamten numerischen Berechnungen werden wir mit~ = m = 1 durchführen, sie wurden aber variabel implementiert.

4.1 Umsetzung der Anwendung

Zunächst wollen wir in den folgenden Abschnitten die Umsetzung der numerischenBerechnungen beschreiben. Angefangen mit dem Setting erläutern wir anschließend dieverwendeten numerischen Methoden, bevor wir auf die explizite Implementation unddie Konvergenz der Simulation bei unterschiedlicher Diskretisierung eingehen.

4.1.1 Rahmenbedingungen

Das System, das wir studieren, wird durch die zweidimensionale Schrödingergleichung

iε∂tΨ ε(t, x, y) = −ε2

2 ∂2xΨ

ε(t, x, y)− 12∂

2yΨ

ε(t, x, y) + V (x, y)Ψ ε(t, x, y) (4.1)

mit x, y ∈ R beschrieben. Das Teilchen wird durch das harmonische Potential

V (x, y) = −12x

2 + 12ω(x)2y2

eingesperrt. Die x-abhängige Frequenz ist dabei gegeben durch ω(x) =√

2(1 + x4). DieSchrödingergleichung (4.1) hat ihren Ursprung in

i∂tΨ ε(t, x, y) =(−1

2∂2x −

12∂

2y −

12(εx)2 + 1

2ω(εx)2y2)Ψ ε(t, x, y)

Page 50: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 4. Numerische Simulation

mit x = εx. In x- und y-Richtung wird das Potential also asymmetrisch skaliert. Al-ternativ kann man das Potential auch in y-Richtung statt in x-Richtung skalieren. Indiesem Fall kann man jedoch die Gesamtenergie trotz beliebig kleinem ε nicht beschrän-ken und es ist nötig weitere zusätzliche Anpassungen vorzunehmen. In Abbildung 4.1ist die Gestalt des Potentials mittels eines Kontur- und eines HSL-Plots illustriert. Zurbesseren graphischen Darstellung haben wir dabei im Konturplot Werte größer Zweiabgeschnitten. Beim HSL-Plot werden komplexe Werte mit Farben dargestellt. Dem Ar-gument, also dem Winkel beziehungsweise der Phase, wird ein Farbton zugeordnet undder Betrag wird durch die Helligkeit veranschaulicht. Das Potential ist so konstruiert,dass es für wachsendes |x| sowohl schmaler als auch steiler wird und dadurch eine Art„Tal“ bildet. Hier und im Folgenden sei bei numerischen Betrachtungen x ∈ Ix = [−5, 5]und y ∈ Iy = [−4, 4]. Zudem bezeichnen wir, angelehnt an die geometrischen Eigen-schaften vieler physikalischer Systeme, die Richtung der langsamen Freiheitsgrade alslongitudinale und die der schnellen Freiheitsgrade als transversale Richtung. In Abbil-dung 4.2 ist das Potential vereinfacht dreidimensional dargestellt, wobei wir dieses Malbei Zehn abgeschnitten haben.

Der Hamiltonoperator der schnellen Freiheitsgrade Hf = −12∂

2y + V (x, y) hat die

Form eines harmonischen Oszillators. Dieser ist ein Standardproblem aus der Literatur,weshalb uns die Lösung des Eigenwertproblems der schnellen Freiheitsgrade

Hf(x)φj(x, ·) = Ej(x)φj(x, ·)

mit j ∈ N bekannt ist. Das Spektrum ist diskret und die Eigenfunktionen sind gegebendurch

φj(x, y) =

√√√√ ω(x)12

2j j!π12Hj

(ω(x)

12 y)

e−ω(x)

2 y2, (4.2)

wobei Hj die Hermiteschen Polynome bezeichnet, deren Definition und allgemeine Ei-genschaften im Appendix A.1 zu finden sind. Die zu (4.2) korrespondierenden Eigen-werte sind gegeben durch

Ej(x) = ω(x)(j + 1

2

)− 1

2x2 (4.3)

und jeweils durch eine Lücke von den anderen Eigenwerten getrennt. Insbesondere istalso E0(x) ein einzelner isolierter Eigenwert. Im Hinblick auf Kapitel 3 wählen wir alsoj = 0 und betrachten den Fall I = {0}. Die Wellenfunktion Ψ ε zum Zeitpunkt t = 0 seigegeben durch

Ψ ε0 (x, y) := Ψ ε(0, x, y) = ψε(0, x)φ0(x, y).

Präziser formuliert sei sie also zusammengesetzt aus einer gaußförmigen Wellenfunktionψε0(x) := ψε(0, x) mit

ψε(0, x) = 1(2π)

14

e−12x

2 eix·px/ε

und dem Grundzustand φ0(x, y) des harmonischen Oszillators

φ0(x, y) =(ω(x)π

) 14

e−ω(x)

2 y2.

42

Page 51: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

4.1. Umsetzung der Anwendung

-4 -2 0 2 4

Longitudinale x

-4

-2

0

2

4

Tra

nsv

ersa

ley

-12 -10 -8 -7 -5 -3 -2 0 2 3

(a) Konturplot

-4 -2 0 2 4Longitudinale x

-4

-2

0

2

4

Tra

nsv

ersa

ley

(b) HSL-Plot

Abbildung 4.1: Darstellung des Potentials V .

Der Anfangsimpuls in die longitudinale x-Richtung wird dabei mit px bezeichnet. Fürunsere Zwecke werden wir, falls nicht anders angegeben, stets px = 1 verwenden.

Sind ε und px hinreichend klein, oszilliert das Teilchen innerhalb des Potentials undwird von diesem im Bereich des Sattelpunkts festgehalten. Der Konturplot in Abbil-dung 4.3 zeigt dieses Verhalten für ε = 1/8. Das Betragsquadrat der Wellenfunktion

43

Page 52: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 4. Numerische Simulation

Longitudinale x

-4-2

02

4Tra

nsver

saley

-4

-2

0

2

4

V(x, y

)

-10

-5

0

5

Abbildung 4.2: Dreidimensionaler Plot des Potentials V .

|Ψ ε|2 ist zu unterschiedlichen Zeiten zusammen mit dem Potential V abgebildet. DasTeilchen bewegt sich zunächst nach rechts und wird dann an der Engstelle des Potentialsreflektiert. Das Oszillationsverhalten lässt sich auch in dem HSL-Plot in Abbildung 4.4beobachten.

Wir gehen davon aus, dass die Wellenfunktion Ψ ε in der transversalen Richtungadiabatisch im Grundzustand φ0 verbleibt, was sich mit unserem raumadiabatischenVerständnis deckt, bei dem die schnellen Freiheitsgrade von dieser Eigenfunktion be-schrieben werden. Infolgedessen können wir die raumadiabatische Theorie anwenden.Die Ergebnisse der vorangegangenen Kapitel 2 und 3 geben uns als Schrödingerglei-chung der langsamen Freiheitsgrade

iε∂tψε(t, x) = HεBOψ

ε(t, x)

mit dem effektiven Born-Oppenheimer-Hamiltonoperator (3.31) in Form von

HεBO(x) = −1

2ε2∂2x + E0(x) + 1

2ε2VBH(x)− ε2M(x, p), (4.4)

wobei der Berry-Zusammenhangskoeffizient aufgrund der Reellwertigkeit der Eigenfunk-tion φ0(x, y) weggeeicht wurde undM(x, p) der Korrekturterm ist, den wir testen wol-len. Der isolierte Eigenwert im Grundzustand ist E0(x) =

(ω(x)− x2)/2 und mit

∇xφ0(x, y) = ω′(x)2π

14

(1

2ω(x)34− ω(x)

14 y2

)e−

ω(x)2 y2 (4.5)

44

Page 53: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

4.1. Umsetzung der Anwendung

-4 -2 0 2 4-2

-1

0

1

2T

ran

sver

saley

t = 0

-4 -2 0 2 4-2

-1

0

1

2t = 1

-4 -2 0 2 4-2

-1

0

1

2

Tra

nsv

ersa

ley

t = 2

-4 -2 0 2 4-2

-1

0

1

2t = 3

-4 -2 0 2 4

Longitudinale x

-2

-1

0

1

2

Tra

nsv

ersa

ley

t = 4

-4 -2 0 2 4

Longitudinale x

-2

-1

0

1

2t = 5

Abbildung 4.3: Konturplot des Betragsquadrats der Wellenfunktion |Ψ ε|2 zu-sammen mit dem Potential V für ε = 1/8 zu unterschiedlichenZeitpunkten.

sowie (3.7) ermitteln wir das Born-Huang-Potential

VBH(x) = 〈∇xφ0(x, y),∇xφ0(x, y)〉Hf=∫R

(∇xφ0(x, y)

)2 dy

=∫R

ω′(x)2

4√πω(x)

( 14ω(x) − y

2 + ω(x)y4)

e−ω(x)y2 dy

= ω′(x)2

4√πω(x)

(1

4ω(x) ·√π√

ω(x)−√π

2ω(x)32

+ ω(x) 3√π

4ω(x)52

)

= ω′(x)2

8ω(x)

( 12ω(x) −

1ω(x) + 3

2ω(x)

)= ω′(x)2

8ω(x)2 . (4.6)

45

Page 54: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 4. Numerische Simulation

-4 -2 0 2 4

-2

0

2

Tra

nsv

ersa

ley

t = 0

-4 -2 0 2 4

-2

0

2

t = 1

-4 -2 0 2 4

-2

0

2

Tra

nsv

ersa

ley

t = 2

-4 -2 0 2 4

-2

0

2

t = 3

-4 -2 0 2 4

Longitudinale x

-2

0

2

Tra

nsv

ersa

ley

t = 4

-4 -2 0 2 4

Longitudinale x

-2

0

2

t = 5

Abbildung 4.4: HSL-Plot der Wellenfunktion Ψ ε für ε = 1/8 zu unterschiedli-chen Zeitpunkten.

Für die Berechnung des Korrekturterms M aus (4.4) benötigen wir die expliziteDarstellung von ⟨

∇xφ0(x, y),(Hf(x)− E0(x)

)−1P⊥0 (x)∇xφ0(x, y)

⟩Hf

(4.7)

aus (3.29) mit der reduzierten Resolvente

R(x) :=(Hf(x)− E0(x)

)−1P⊥0 (x) =

∞∑i=1

|φi(x, y)〉 〈φi(x, y)|Ei(x)− E0(x) . (4.8)

Dazu machen wir uns zunächst klar, dass wir ∇xφ0(x, y) aus (4.5) auch in Abhängigkeitvon φ0(x, y) schreiben können:

∇xφ0(x, y) = ω′(x)4ω(x)φ0(x, y)− ω′(x)

2 y2φ0(x, y).

46

Page 55: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

4.1. Umsetzung der Anwendung

In Kombination mit φi(x, y) = φi(x, y) ist dann

〈φi(x, y),∇xφ0(x, y)〉Hf=∫Rφi(x, y)∇xφ0(x, y) dy

= ω′(x)4ω(x)

∫Rφi(x, y)φ0(x, y) dy − ω′(x)

2

∫Rφi(x, y)y2φ0(x, y) dy .

Mithilfe von (A.1) und (A.2) können wir diesen Term für unterschiedliche i ∈ N aus-werten und stellen fest, dass er nur für i = 2 ungleich Null ist. Dementsprechend ist

〈φ2(x, y),∇xφ0(x, y)〉Hf= −ω

′(x)2

∫Rφ2(x, y)y2φ0(x, y) dy = − ω′(x)

232ω(x)

. (4.9)

Indem wir beachten, dass φi(x, y) und ∇xφi(x, y) selbstadjungiert sind und alle Termemit i 6= 2 wegfallen, sind wir infolgedessen in der Lage, (4.7) zu ermitteln:

⟨∇xφ0(x, y), R(x)∇xφ0(x, y)

⟩Hf

=∞∑i=1

⟨∇xφ0(x, y), φi(x)

⟩⟨φi(x),∇xφ0(x, y)

⟩Ei(x)− E0(x)

=⟨φ2(x),∇xφ0(x, y)

⟩2E2(x)− E0(x)

= 12ω(x)

(− ω′(x)

232ω(x)

)2

= ω′(x)2

16ω(x)3 .

Der Korrekturterm aus (3.28) hat schließlich in unserem konkreten Fall die Form

M(x) = −ε2∂xω′(x)2

16ω(x)3∂x. (4.10)

Alles zusammen ist der effektive Hamiltonoperator mit (4.6) und (4.10) sowie derFrequenz ω(x) =

√2(1 + x4) und ihrer Ableitung

ω′(x) = 2√

2x3√

1 + x4,

explizit gegeben durch

HεBO(x) = −1

2ε2∂2x + E0(x) + 1

2ε2VBH(x)− ε2M(x) (4.11)

= −12ε

2∂2x + ω(x)− x2

2 + 12ε

2 ω′(x)2

8ω(x)2 + ε4∂xω′(x)2

16ω(x)3∂x

= −12ε

2∂2x +

√2(1 + x4)− x2

2 + 12ε

2 x6

2(1 + x4)2 + ε4∂xx6

4√

2(1 + x4)52∂x.

47

Page 56: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 4. Numerische Simulation

Um die Born-Oppenheimer-Approximation zweiter Ordnung vollständig und korrektdurchführen zu können, benötigen wir darüber hinaus die Adjungierte des in Kapitel3.2 eingeführten Operators U ε(1). Dieser hat mit (3.21) die Form

U ε ∗(1) = U∗0 + εU∗1 = U∗0 + εB∗ · pU∗0 .

Mittels der Definition aus (3.15) gilt für B∗ zusammen mit der Resolvente (4.8)

B∗ = iP⊥0 (x)(Hf(x)− Ej(x)

)−1(∇xP0(x))P0(x)

= i∞∑i=1

|φi(x, y)〉 〈φi(x, y)|Ei(x)− E0(x)

(|∇xφ0(x, y)〉 〈φ0(x, y)|

+ |φ0(x, y)〉 〈∇xφ0(x, y)|)|φ0(x, y)〉 〈φ0(x, y)|

= i∞∑i=1

|φi(x, y)〉 〈φi(x, y)|∇xφ0(x, y)〉 〈φ0(x, y)|Ei(x)− E0(x)

= i |φ2(x, y)〉 〈φ2(x, y)|∇xφ0(x, y)〉 〈φ0(x, y)|Ei(x)− E0(x)

= −i ω′(x)2

52ω(x)2

|φ2(x, y)〉 〈φ0(x, y)| .

Dabei haben wir (3.2), das Wegfallen von Summanden mit i 6= 2 sowie (4.9) verwendet.Nach (4.3) ist zudem Ei(x)− E0(x) = 2ω(x). Mit B∗ können wir nun U ε ∗(1) berechnen.

U ε ∗(1) = U∗0 + εB∗ · pU∗0= |φ0(x, y)〉 − iε2B∗

(|∇xφ0(x, y)〉+ |φ0(x, y)〉 · ∇x

)= |φ0(x, y)〉 − ε2 ω′(x)

252ω(x)2

|φ2(x, y)〉 · ∇x.

Dementsprechend gilt also angewandt auf ψε:

U ε ∗(1)ψε(t, x) = ψε(t, x) |φ0(x, y)〉 − ε2 ω′(x)

252ω(x)2

|φ2(x, y)〉 · ∇xψε(t, x)

= ψε(t, x) |φ0(x, y)〉 − ε2 x3

4(1 + x4)32|φ2(x, y)〉 · ∇xψε(t, x).

Unser Ziel ist es im Folgenden, die Qualität der Approximation in Abhängigkeit von demM-Term, also dem hintersten Term in (4.11), zu untersuchen und die theoretischen Re-sultate aus Kapitel 3 an dem obigen Beispiel zu veranschaulichen. Dies bewerkstelligenwir, indem wir den Fehler der konventionellen Born-Oppenheimer-Approximation, dasheißt die normierte Differenz von der numerischen Approximation ohne den M-Termund der numerischen Lösung der zweidimensionalen Schrödingergleichung des Gesamt-systems, mit dem Fehler der Born-Oppenheimer-Approximation zweiter Ordnung, alsounter Berücksichtigung des M-Terms, vergleichen. Bezüglich ersterem verwenden wiralso Hε

eff aus (3.8) und berechnen numerisch∥∥∥e−iHε tεΨ ε0 − U ε ∗(1)e

−iHεeff

tεψε0

∥∥∥L2. (4.12)

48

Page 57: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

4.1. Umsetzung der Anwendung

Demgegenüber reicht für die vollständige Born-Oppenheimer-Approximation zweiterOrdnung mittels (4.11) die Verwendung des adiabatischen Unterraums für die Anfangs-daten nicht aus. Zwar werden wir sehen, dass im vorliegenden Fall die Approximationdurch den Korrekturterm alleine schon etwas verbessert wird, jedoch nicht um eine gan-ze Ordnung. Dafür ist es indessen nötig, dass die Startwerte aus dem superadiabatischenUnterraum stammen. Dies erzwingen wir, indem wir U ε ∗(1)ψ

ε0 als Startwert verwenden.

Mit U ε(1)Uε ∗(1)ψ

ε0 = ψε0 berechnen wir dann numerisch

∥∥∥e−iHε tεU ε ∗(1)ψ

ε0 − U ε ∗(1)e

−iHεBO

tεψε0

∥∥∥L2. (4.13)

4.1.2 Numerische Methoden

An dieser Stelle wollen wir kurz die numerischen Methoden ansprechen, auf die wir fürunsere Berechnungen zurückgreifen. Eine dieser Methoden ist das sogenannte Strang-Splitting [Str68], mit dessen Hilfe wir die Zeitentwicklung des Systems approximieren.Allgemein wird dabei der Zeitentwicklungsoperator eines gegebenen HamiltonoperatorsH = H0 + V auf L2(Rd) für ein hinreichend großes n ∈ N mittels

e−iHt ≈(e−iV t

2n e−iH0te−iV t2n

)nin kleinere Zeitabschnitte unterteilt. Der dabei entstehende Fehler ist von der Ord-nung O

(1/n2). Für das Potential V ist e−iV t ein Multiplikationsoperator im Ortsraum.

Enthält H0 den ableitenden Teil des Operators, zum Beispiel in Form eines Laplace-Operators, können wir diesen durch eine Fouriertransformation in einen Multiplikati-onsoperator im Fourierraum umwandeln. Der Gradient ∇ wird auf diese Weise zu einerMultiplikation mit ik und ∆ dementsprechend zu einer Multiplikation mit −k2, wennk die Variable im Fourierraum bezeichnet:

e−iH0t = F−1e−i k22 tF .

Die Fouriertransformation lässt sich numerisch sehr effizient berechnen, indem der so-genannte FFT-Algorithmus angewendet wird. Die Bezeichnung stammt aus dem Engli-schen und steht für „fast Fourier transform“. Der Algorithmus steht in vielen Program-miersprachen bereits durch diverse Programmbibliotheken zur Verfügung. Besonderseffektiv und schnell ist dieser, wenn die Anzahl der verwendeten Gitterpunkte eineZweierpotenz ist. Dies werden wir bei der Wahl des Gitters berücksichtigen.

Unter Berücksichtigung der angesprochenen Punkte ist unsere Vorgehensweise zurApproximation der Zeitentwicklung also die Folgende: Der Anfangswert wird zunächstmit dem Potentialteil des Propagators eines halben Zeitschritts e−iV t/2n multipliziert.Dies entspricht der Multiplikation mit einer komplexen Exponentialfunktion. Anschlie-ßend führen wir eine Fouriertransformation durch, um den ableitenden Teil des Pro-pagators in eine Multiplikation im Fourierraum zu überführen. Dadurch müssen wirnur eine weitere Multiplikation mit einer komplexen Exponentialfunktion durchführen.Darauf folgt eine Fourierrücktransformation, um die Wellenfunktion im ursprünglichenOrtsraum zu erhalten. Als letzter Schritt wird das Ganze erneut mit dem Potentialteil

49

Page 58: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 4. Numerische Simulation

des Propagators eines halben Zeitschritts multipliziert. Dieses Prozedere wird hinrei-chend oft wiederholt:

e−iHt =(

e−iV t2nF−1e−i k2

2tnFe−iV t

2n

)n+O

(1n2

).

Die Zeitentwicklung gegeben durch den Operator (4.11) können wir also analog durch

e−iHεBO

tεψε0(x) =

(e−if(x) t

ε1

2nF−1e−iε2 k22

12n

(Id + ikF

(iεv(x) tn

)F−1ik

)e−iε2 k2

2tε

12nFe−if(x) t

ε1

2n

)nψε0(x) +O

(1n2

)approximieren, indem wir zusätzlich den Propagator des Korrekturteils mit den Gra-dienten im Exponenten entwickeln. Außerdem verwenden wir die beiden Terme f(x) =E0(x) + 1

2ε2VBH(x) und

v(x) = −ε2 x6

4√

2(1 + x4)52. (4.14)

Es ist M(x) = ∂xv(x)∂x. Die Abbildung 4.5 soll einen Eindruck vermitteln wie derEinfluss von (4.6) und (4.14) in (4.11) verläuft.

−6 −4 −2 0 2 4 60.000.020.040.060.080.100.120.140.160.18

x

VB

H(x

)

(a) Born-Huang-Potential VBH

−6 −4 −2 0 2 4 6−0.0006

−0.0005

−0.0004

−0.0003

−0.0002

−0.0001

0.0000

x

v(x

)

(b) v-Term

Abbildung 4.5: Das Born-Huang-Potential aus (4.6) und der v-Term (4.14) fürε = 1/8 im Vergleich.

4.1.3 Implementierung

Den groben Aufbau der Implementierung wollen wir im Folgenden andeuten. Für denvollständigen Quellcode sei auf Appendix A.2 verwiesen. Das Programm ist in der uni-versellen Programmiersprache Python [Pyt, Ros95] implementiert. Zusätzlich nutzenwir unter anderem die Programmbibliotheken NumPy [Oli06] und Matplotlib [Hun07].Der Programmcode selbst besteht aus drei Teilen: Dem Hauptprogramm run.py, ei-ner Sammlung diverser Funktionen function_data_base.py und der Parameterda-tei parameter.py. Ersteres ist ausführbar und steuert die Anzeige und Initiierung al-ler verfügbaren Berechnungen. Dabei werden zu Beginn function_data_base.py und

50

Page 59: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

4.1. Umsetzung der Anwendung

parameter.py importiert. Je nach Auswahl werden die entsprechenden Funktionen ausder Sammlung aufgerufen, um die gewünschten Berechnungen durchzuführen. Sind kei-ne expliziten Angaben gemacht, werden die Werte verwendet, die in der Parameterdateihinterlegt sind. Dort findet sich unter anderem das reduzierte Plancksche Wirkungs-quantum ~, die Teilchenmasse m, der Anfangsimpuls px und der adiabatische Para-meter ε wieder. Der Parameter h = t/n gibt die Schrittgröße im Strang-Splitting an.Des Weiteren wird hier die Größe der Intervalle Ix und Iy festgelegt. Unsere Wahl mitIx = [−5, 5] und Iy = [−4, 4] beruht darauf, dass die Funktionswerte der WellenfunktionΨ ε zum Rand hin hinreichend klein und damit vernachlässigbar werden. Das IntervallIx× Iy wird diskretisiert, indem wir es in ein Gitter mit 2nx×2ny Gitterpunkten unter-teilen. Die Gittergröße wird dabei über die Parameter nx und ny geregelt. Die Abständeder Gitterpunkte zueinander sind äquidistant.

4.1.4 Konvergenztests

Es ist wichtig für den Vergleich der analytischen Fehler der konventionellen Born-Oppenheimer-Approximation erster Ordnung (4.12) und der Born-Oppenheimer-Ap-proximation zweiter Ordnung (4.13), dass diese nicht durch numerisch bedingte Fehlergestört oder verfälscht werden. Wir müssen also darauf achten, dass die numerischenFehler kleiner sind als die analytischen Fehler, an denen wir eigentlich interessiert sind,ohne dabei einen enormen Rechenaufwand zu betreiben, der zu langen Rechenzeitenführt. Dazu betrachten wir zum Einen die Abhängigkeit von der Gittergröße der Diskre-tisierung und versuchen geeignete Werte für die Parameter nx und ny unseres (2nx , 2ny )-Gitters zu finden. Zum Anderen ist die Zeitschrittweite h beim Strang-Splitting einwichtiger Faktor.

Wir beginnen mit der Untersuchung der Gittergröße. Dazu lassen wir den Korrek-turterm erst einmal außen vor und berechnen die durch (4.11) ohne Korrekturtermgegebene Zeitentwicklung von Ψ ε0 bis zu der makroskopischen Gesamtzeit T = 12 fürunterschiedliche Gittergrößen. Dies entspricht in etwa der Zeit einer Periode mit jeweilseiner Reflexion rechts und links, siehe Abbildung 4.6. Da bei unterschiedlichen Gitter-größen die Anzahl ausgewerteter Punkte variiert, berücksichtigen wir nur Punkte, diemit den jeweiligen Punkten des anderen Gitters vergleichbar sind. Zu jedem Zeitpunkt tberechnen wir also für verschiedene nx = ny punktweise die Differenz der in den diversenGittern entwickelten Wellenfunktionen und betrachten anschließend die L2-Norm, dasheißt konkret ‖Ψ εeff,mx

− Ψ εeff,nx‖L2 . Analog gehen wir bezüglich der tatsächlichen Zeit-

entwicklung im zweidimensionalen Gesamtsystem vor und werten ‖Ψ εges,mx− Ψ εges,nx

‖L2

aus. Je kleiner ε ist, desto präziser muss die Numerik sein. Deshalb verwenden wir füralle nachfolgenden Tests, insbesondere der Gitteranalyse, ε = 1/16, entsprechend denkleinsten von uns beabsichtigten Wert für ε. Als Schrittgröße h verwenden wir vorläufig10−5, wobei sich herausstellen wird, dass eine größere Schrittweite ausreichend ist. DasErgebnis des Konvergenztests sieht man in Abbildung 4.7. Sie zeigt, dass die Werte fürΨ εeff schneller stabil werden als für Ψ εges. Erhöhen wir die Gittergröße von nx = 6 suk-zessiv auf nx = 8, verringert sich der numerische Fehler der effektiven Approximationjeweils um mehrere Zehnerpotenzen. Ab nx = 8 stagniert der Fehler im Bereich von10−5 und eine weitere Verfeinerung des Gitters führt nur zu einer Verbesserung der Ge-nauigkeit von nicht mehr als 10−5. Problematischer ist hingegen die zweidimensionale

51

Page 60: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 4. Numerische Simulation

-4 -2 0 2 4-2

-1

0

1

2

Tra

nsv

ersa

ley

t = 0

-4 -2 0 2 4-2

-1

0

1

2t = 2

-4 -2 0 2 4-2

-1

0

1

2

Tra

nsv

ersa

ley

t = 4

-4 -2 0 2 4-2

-1

0

1

2t = 6

-4 -2 0 2 4-2

-1

0

1

2

Tra

nsv

ersa

ley

t = 8

-4 -2 0 2 4-2

-1

0

1

2t = 10

-4 -2 0 2 4

Longitudinale x

-2

-1

0

1

2

Tra

nsv

ersa

ley

t = 11

-4 -2 0 2 4

Longitudinale x

-2

-1

0

1

2t = 12

Abbildung 4.6: Konturplot des Betragsquadrats der Wellenfunktion |Ψ ε|2 zu-sammen mit dem Potential V für ε = 1/16 zu verschiedenenZeitpunkten.

Art der Berechnung. Hier wird der Fehler langsamer klein und ist erst ab einem Sprungvon 29 auf 210 Gitterpunkten in besagtem Bereich. Für die Kombination der Zeitinter-valle zusammen mit den ε-Werten, die wir betrachten, erwarten wir einen ungefährenanalytischen Fehler von minimal 10−4. Dieser Wert hat sich in vorangegangenen Ver-suchen, auf die wir hier aber nicht weiter eingehen, abgezeichnet. Die durch nx = 9gegebene Genauigkeit sollte dementsprechend angemessen sein.

Es sei angemerkt, dass die Programmlaufzeiten bei einer größeren Anzahl an Gitter-punkten exponentiell ansteigen. Dies ist vor allem für die zweidimensionale Berechnungvon großer Bedeutung. Während sie für eine Gittergröße von nx = 6 keine 5 Minutendauert, ist für nx = 10 bereits eine Rechenzeit von mehr als 40 Stunden nötig. Um-so wichtiger ist es also, sinnvolle Werte für die Durchführung der Approximation zu

52

Page 61: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

4.1. Umsetzung der Anwendung

0 2 4 6 8 10 12

Zeit t

0

10−9

10−8

10−7

10−6

10−5

10−4

10−3

10−2

10−1

100

101

‖Ψε mx−

Ψε nx‖ L

2

eff : mx = 6, nx = 7

ges : mx = 6, nx = 7

eff : mx = 7, nx = 8

ges : mx = 7, nx = 8

eff : mx = 8, nx = 9

ges : mx = 8, nx = 9

eff : mx = 9, nx = 10

ges : mx = 9, nx = 10

Abbildung 4.7: Vergleich verschiedener Gittergrößen mit nx = ny für ε = 1/16und h = 10−5.

verwenden. Auf die Unterschiede in der Laufzeit und den messbaren Geschwindigkeits-vorteil der Born-Oppenheimer-Approximation gehen wir in Abschnitt 4.2.2 nochmalgenauer ein.

Mit festem nx = 9 testen wir nun diverse Werte für ny. Abbildung 4.8 zeigt dasErgebnis. Die Abweichung zwischen ny = 6 und ny = 7 bei der zweidimensionalen Be-rechnung vergrößert sich im Verlauf der dargestellten Zeit um eine Zehnerpotenz. DerSprung von ny = 9 auf ny = 10 verursacht in derselben Zeit nur eine Verschlechterungum etwa eine halbe Zehnerpotenz. Insgesamt verbessert sich der numerische Fehler mitjeder Verfeinerung des Gitters. Da die getesteten Werte der effektiven Berechnung soklein sind, dass sie an dieser Stelle nicht weiter berücksichtigt werden müssen, entschei-den wir uns für ny = 8. Hier bleibt der Fehler für Zeiten, die für uns relevant sind, stetsunter 10−5. Verwenden wir in die transversale Richtung also 28 Gitterpunkte, sollte dernumerische Fehler unter dem erwarteten analytischen Fehler bleiben. Zusammengefasstkönnen wir folglich für alle ε, die nicht größer sind als 1/16, die Werte nx = 9 undny = 8, das heißt ein 29 × 28-Gitter, verwenden.

Kommen wir nun zu der beim Strang-Splitting verwendeten Zeitschrittweite h. Demvorangegangenen Konvergenztest bezüglich der Gittergröße entsprechend berechnen wirjeweils die ein- und die zweidimensionale Lösung bis zur makroskopischen GesamtzeitT = 12 für variierende Zeitschrittweiten. In Abbildung 4.9 ist zu sehen, dass im ef-fektiven Fall schon für kleine h eine hohe Genauigkeit gegeben ist. Das Augenmerkliegt also erneut auf der zweidimensionalen Berechnung im Gesamtsystem. Die Diffe-

53

Page 62: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 4. Numerische Simulation

0 2 4 6 8 10 12

Zeit t

0

10−9

10−8

10−7

10−6

10−5

10−4

‖Ψε my−

Ψε ny‖ L

2

eff : my = 6, ny = 7

ges : my = 6, ny = 7

eff : my = 7, ny = 8

ges : my = 7, ny = 8

eff : my = 8, ny = 9

ges : my = 8, ny = 9

eff : my = 9, ny = 10

ges : my = 9, ny = 10

Abbildung 4.8: Vergleich verschiedener Gittergrößen mit festem nx = 9 undvariablen ny für ε = 1/16 und h = 10−5.

0 2 4 6 8 10 12

Zeit t

0

10−9

10−8

10−7

10−6

10−5

10−4

‖Ψε h1−

Ψε h2‖ L

2 eff : h1 = 140000, h2 = 1

50000

ges : h1 = 140000, h2 = 1

50000

eff : h1 = 150000, h2 = 1

60000

ges : h1 = 150000, h2 = 1

60000

eff : h1 = 160000, h2 = 1

70000

ges : h1 = 160000, h2 = 1

70000

eff : h1 = 170000, h2 = 1

80000

ges : h1 = 170000, h2 = 1

80000

Abbildung 4.9: Vergleich verschiedener Zeitschrittweiten h mit nx = 9, ny = 8und ε = 1/16.

54

Page 63: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

4.2. Durchführung der Approximation

renz ‖Ψ εges,h1− Ψ εges,h2

‖L2 wird für die gewählten Schrittweiten in gleichmäßigen Sprün-

gen kleiner. Im Vergleich zueinander haben die Kurven einen ähnlichen Verlauf mitsteigender Tendenz. Nach einem starken Anstieg am Anfang wachsen sie ab ungefährt = 4 nur noch langsam. Abermals müssen wir uns für eine Kurve entscheiden, die imBereich von 10−5 liegt, weshalb wir h = 6 · 10−5 für unsere Zwecke wählen.

Zusammenfassend halten wir als geeignete Werte für die Parameter in dem nachfol-genden Abschnitt nx = 9, ny = 8 und h = 6 · 10−5 fest.

4.2 Durchführung der ApproximationDie Ergebnisse der numerischen Berechnungen wollen wir nun präsentieren. Als abschlie-ßendes Resultat soll die Theorie an dem zuvor eingeführten Praxisbeispiel bestätigt undveranschaulicht werden. Auf unterschiedliche Arten werden wir dabei den analytischenFehler der Approximation darstellen und deuten. Dazu approximieren wir sowohl beifestem ε über längere Zeiten hinweg als auch mit verschiedenen Werten für ε, wo wirdie Ergebnisse zu einem festen Zeitpunkt betrachten.

4.2.1 Vergleich der Fehler

Die Theorie über die Born-Oppenheimer-Approximation aus den vorangegangenen Ka-piteln hält, was sie verspricht. In Abbildung 4.10 haben wir den Logarithmus des ana-lytischen Fehlers in Abhängigkeit von log(ε) dargestellt. Dadurch können wir an der

-4.0 -3.5 -3.0 -2.5 -2.0 -1.5 -1.0log2(ε)

-14

-12

-10

-8

-6

-4

log

2(‖

Ψε−

Ψε eff‖ L

2)

ohne VBH und Mohne Mmit M, ohne U ε ∗

(1)

mit M und U ε ∗(1)

Abbildung 4.10: Vergleich der verschiedenen Arten der Approximation in Ab-hängigkeit von ε zum Zeitpunkt t = 3.

55

Page 64: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 4. Numerische Simulation

Steigung der Kurve die Ordnung der Approximation ablesen. Schon auf den erstenBlick werden großen Unterschiede ersichtlich. Die Approximation in führender Ordnung(rot), also mit nur den ersten beiden Summanden aus (4.11), liefert das schlechteste Er-gebnis. Die konventionelle Born-Oppenheimer-Approximation (blau) verbessert diesesResultat. Der Fehler wird aber auch hier mit kleineren ε nur linear kleiner. Die graugestrichelte Linie zeigt eine Vergleichsgerade mit der Steigung Eins. Entsprechend demraumadiabatischen Theorem (2.32) bildet sie eine obere Schranke, unter der alle abge-bildeten Fehlerwerte liegen. Nehmen wir den Korrekturterm hinzu, jedoch nicht U ε ∗(1),erzielen wir erneut eine kleine Verbesserung (gelb). Zum späteren Zeitpunkt t = 12 istdiese Verbesserung sogar noch deutlicher, vergleiche dazu Abbildung 4.11. Dies ist be-

-4.0 -3.5 -3.0 -2.5 -2.0 -1.5 -1.0log2(ε)

-14

-12

-10

-8

-6

-4

-2

log

2(‖

Ψε−

Ψε eff‖ L

2)

ohne VBH und Mohne Mmit M, ohne U ε ∗

(1)

mit M und U ε ∗(1)

Abbildung 4.11: Vergleich der verschiedenen Arten der Approximation in Ab-hängigkeit von ε zum Zeitpunkt t = 12.

merkenswert und nicht selbstverständlich, da uns die Theorie aus Kapitel 3.2 nur eineVerbesserung durch denM-Term in Verbindung mit U ε(1) respektive U

ε ∗(1) verspricht. Die

Qualität und der damit verbundene Mehrwert der Born-Oppenheimer-Approximationzweiter Ordnung (grün) ist offensichtlich und springt sofort ins Auge. Die Steigung dergrünen Kurve ist signifikant steiler. Sie verhält sich ähnlich zu der grau gestricheltenVergleichsgerade mit der Steigung 2. Der Fehler verringert sich also quadratisch, wasmit (3.32) konform geht. Den kleinsten Fehler erhalten wir hier für ε = 1/16 und dieserbeträgt wie ursprünglich erwartet in etwa 10−4.

Beeindruckend ist auch der Blick auf den zeitabhängigen Fehlervergleich. Mit demadiabatischen Parameter ε = 1/8 ist der Verlauf in Abbildung 4.12 dargestellt. In dem

56

Page 65: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

4.2. Durchführung der Approximation

0 2 4 6 8 10 12

Zeit t

0.000

0.005

0.010

0.015

0.020

0.025

0.030

0.035

0.040

0.045

‖Ψε−

Ψε eff‖ L

2

ohne VBH und Mohne Mmit M, ohne U ε ∗

(1)

mit M und U ε ∗(1)

Abbildung 4.12: Vergleich der verschiedenen Arten der Approximation in Ab-hängigkeit von der Zeit t mit ε = 1/8.

abgebildeten Zeitintervall von 0 bis 12 bleibt der Fehler der Born-Oppenheimer-Appro-ximation zweiter Ordnung augenscheinlich nahezu konstant klein. Er beträgt durch-schnittlich etwa 4,4 · 10−4 für ε = 1/8 und 1,1 · 10−4 für ε = 1/16. Der Fehler derkonventionellen Approximation steigt hingegen sichtlich fortlaufend an. Zum Zeitpunktt = 5 beträgt er bereits ein Hundertstel und bei t = 12 mehr als 0,015. Das mit Ab-stand schlechteste Resultat liefert die einfache Approximation in führender Ordnung.Über längere Zeiten hinweg ist diese zunehmend unbrauchbar. Der Fehler bei t = 12beträgt mehr als das 2,5-fache der konventionellen und das 77-fache der Approximationzweiter Ordnung. Letzterer ist für ε = 1/16 sogar fast um das 170-fache größer. Be-gutachtet man die Born-Oppenheimer-Approximation zweiter Ordnung, verbessert sichder Fehler zum Zeitpunkt t = 12 gegenüber der konventionellen Approximation um denFaktor 30 für ε = 1/8 und 70 für ε = 1/16.

Ein weiteres interessantes Kriterium ist die Fehlerentwicklung in Abhängigkeit vomAnfangsimpuls. In Abbildung 4.13 testen wir die verschiedenen Arten der Approxi-mation bei jeweils unterschiedlichen Werten für px. Auffallend ist, dass für große An-fangsimpulse der Fehler der Approximation in führender Ordnung kleiner wird. Dem-gegenüber vergrößert das hinzukommende Born-Huang-Potential der konventionellenBorn-Oppenheimer-Approximation den Fehler zusätzlich. Der Korrekturterm und dasBorn-Huang-Potential haben entgegengesetztes Vorzeichen und eliminieren sich hierteilweise gegenseitig, sodass die Kombination der beiden besser ist. Für kleine Impulsesieht es anders aus. Das Born-Huang-Potential verbessert das Ergebnis deutlich. Die

57

Page 66: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 4. Numerische Simulation

0 2 4 6 8 10 12Zeit t

0.00

0.01

0.02

0.03

0.04

0.05

0.06

‖Ψε−

Ψε eff‖ L

2

ohne VBH undMohneMmitM, ohne U ε ∗

(1)

mitM und U ε ∗(1)

(a) Anfangsimpuls px = 1/4

0 2 4 6 8 10 12Zeit t

0.00

0.01

0.02

0.03

0.04

0.05

‖Ψε−

Ψε eff‖ L

2

ohne VBH undMohneMmitM, ohne U ε ∗

(1)

mitM und U ε ∗(1)

(b) Anfangsimpuls px = 12

0 2 4 6 8 10 12Zeit t

0.000

0.005

0.010

0.015

0.020

0.025

0.030

‖Ψε−

Ψε eff‖ L

2

ohne VBH undMohneMmitM, ohne U ε ∗

(1)

mitM und U ε ∗(1)

(c) Anfangsimpuls px = 3/2

0 2 4 6 8 10 12Zeit t

0.000

0.005

0.010

0.015

0.020

0.025

0.030

0.035

0.040

‖Ψε−

Ψε eff‖ L

2

ohne VBH undMohneMmitM, ohne U ε ∗

(1)

mitM und U ε ∗(1)

(d) Anfangsimpuls px = 2

Abbildung 4.13: Vergleich der verschiedenen Arten der Approximation bei je-weils unterschiedlichen Anfangsimpulsen für ε = 1/8.

Kurve verläuft flacher und damit steigt auch der Fehler im Verlauf der Zeit langsameran. Dennoch kommt keine dieser Varianten an die Born-Oppenheimer-Approximationzweiter Ordnung heran. Bei allen untersuchten Anfangsimpulsen erzielt stets die Ver-wendung des Korrekturterms innerhalb des superadiabatischen Unterraums das besteResultat. Dies rührt daher, dass M(x, p) proportional zu p2 und U ε ∗(1) proportional zup ist. Der Impuls wird also stärker berücksichtigt und fließt bei der Approximation mitein. Dies ist gerade bei großen Impulsen und insbesondere an den Stellen wichtig, andenen VBH undM(x, p) stark einwirken, vergleiche dazu Abbildung 4.5.

4.2.2 Vergleich der Laufzeiten

Die Born-Oppenheimer-Approximation bietet nicht nur einen qualitativen Vorteil, son-dern ist gegenüber der mehrdimensionalen Berechnung im Gesamtsystem um ein Viel-faches schneller. Ein Eindruck, um wie viel genau und wie die unterschiedlichen Ord-nungen sich auf den Rechenaufwand auswirken, soll hier vermittelt werden.

58

Page 67: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

4.2. Durchführung der Approximation

Es sei angemerkt, dass der Programmcode nicht bezüglich der Laufzeit optimiertwurde. Hier ist durchaus Verbesserungspotential vorhanden und es kann zum Beispielbei dem Strang-Splitting mit mehreren Operatoren, wie es bei der Approximation zwei-ter Ordnung auftritt, angesetzt werden. Für unsere Zwecke war eine ausführliche Op-timierung jedoch nicht nötig, da die effektiven Berechnungen für die Zeiten, die wirbetrachtet haben, ausreichend schnell waren. Alle genannten Berechnungen wurden aufeinem Intel R© CoreTM i5-7500 Prozessor mit vier Kernen und jeweils 3,40 GHz durchge-führt. Darüber hinaus verfügte der Rechner über 8 GiB DDR4 2400MHz Arbeitsspei-cher. Die folgende Tabelle 4.1 soll eine Übersicht über die Laufzeiten geben. Angefangenbeim Zeitpunkt t0 = 0 wurde die Zeit bis zum Zeitpunkt T entwickelt.

Typ t0 T Laufzeit

Gesamtsystem 0 12 2:14 Std.Approx. 1. Ord. 0 12 0:18 Min.Approx. 2. Ord. 0 12 1:06 Min.

Tabelle 4.1: Laufzeiten.

Die Laufzeiten wachsen mit der Zeit und der Anzahl der Zeitschritte jeweils linearan und exponentiell mit der Gittergröße. Im vorliegenden Fall wurden als Parameterh = 6 · 10−5, nx = 9 und ny = 8 verwendet. Passt man diese in Abhängigkeit vomadiabatischen Parameter an, ändert sich entsprechend dem Wert von ε auch die Lauf-zeit. Die zweidimensionale Berechnung im Gesamtsystem benötigt mit Abstand ammeisten Rechenzeit. Für komplexere Modelle ist sie daher kaum geeignet. Viel schnel-ler sind im Vergleich dazu die Born-Oppenheimer-Approximationen erster und zweiterOrdnung. Mit 18 Sekunden ist die konventionelle Approximation unschlagbar schnellund unterstreicht noch einmal die Vorzüge der adiabatischen Theorie. Etwas länger alsdie konventionelle Approximation benötigt die Approximation zweiter Ordnung. Behältman allerdings die Resultate des vorangegangenen Abschnitts 4.2.1 im Hinterkopf, istdie Laufzeit noch beeindruckender. Obwohl sie zusätzlich den Fehler verringert, beträgtdie benötigte Rechenzeit bei dieser Größenordnung nur unwesentlich mehr als die derkonventionellen und ein Bruchteil der zweidimensionalen Variante.

4.2.3 Verhalten für längere Zeiten

Ohne im Folgenden erneut Konvergenztests durchzuführen, wollen wir die Entwick-lung auch noch für längere Zeiten betrachten. Dabei muss uns bewusst sein, dass dernumerische Fehler nicht mehr zwingend kleiner als der analytische Fehler sein mussund etwaige Beobachtungen dadurch eventuell etwas getrübt sein könnten. Um diesesProblem in Grenzen zu halten, verwenden wir größere Werte für den adiabatischen Pa-rameter ε als bisher. Der numerische Fehler verliert dadurch etwas an Bedeutung, dader Gesamtfehler größer wird und zusätzlich die Relevanz der Nachkommastellen einwenig abnimmt.

In Abbildung 4.14 ist die Entwicklung der Fehler bis zu der Gesamtzeit T = 50 dar-gestellt. Der unterschiedlich starke Anstieg der Approximation in führender Ordnung,der sich bereits in Abbildung 4.12 angedeutet hat, setzt sich auch für längere Zeiten

59

Page 68: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 4. Numerische Simulation

0 10 20 30 40 50

Zeit t

0.00

0.05

0.10

0.15

0.20

0.25

0.30

0.35

‖Ψε−

Ψε eff‖ L

2

ohne VBH und Mohne Mmit M, ohne U ε ∗

(1)

mit M und U ε ∗(1)

Abbildung 4.14: Vergleich der verschiedenen Arten der Approximation in Ab-hängigkeit von der Zeit t bis T = 50 mit ε = 1/4.

fort. In Bereichen, in denen der Betrag des Born-Huang-Potentials sowie des Korrek-turterms gering sind, wächst der Fehler im Hinblick auf Abbildung 4.5 langsamer, alsan den Stellen, an denen sie einen größeren Anteil an der Fehlerkorrektur haben.

Der anfängliche Sprung der Fehler in Abbildung 4.12 und 4.14 und die geringe Stei-gung der gelben Kurve in Abbildung 4.14 legen die Unterteilung der Konstanten aus(3.32) in einen zeitabhängigen und einen zeitunabhängigen Teil noch einmal offen. Dieentsprechende Wahl des Hamiltonoperators wirkt sich auf die zeitabhängige Kompo-nente und somit auf die Steigung der Kurve aus. Die Wahl der Startwerte aus demadiabatischen beziehungsweise superadiabatischen Unterraum beeinflusst den Achsen-abschnitt, also den statischen zeitunabhängigen Teil.

Die grüne Kurve in Abbildung 4.14, genauer gesagt die Fehlerkurve der Born-Oppenheimer-Approximation zweiter Ordnung, sieht wie bereits in Abbildung 4.12 auchfür längere Zeiten augenscheinlich nahezu konstant aus. Zumindest scheint der Fehlerlangsamer zu wachsen als erwartet. Wir betrachten diesen daher in Abbildung 4.15isoliert für verschiedene Werte des adiabatischen Parameters ε. Trotz Oszillation derKurven erkennen wir nun einen leichten Gesamtzuwachs des Fehlers. Abhängig von εsteigt dieser entsprechend schnell oder langsam an und ist somit nicht konstant. DieKurve mit ε = 1/4 oszilliert am stärksten. Es sei jedoch angemerkt, dass alle drei Kur-ven relativ zu ihren Größenordnungen oszillieren, sodass die Oszillationen bei der Skalain Abbildung 4.15 unterschiedlich ausgeprägt erscheinen. Ein Vergrößern der unteren

60

Page 69: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

4.2. Durchführung der Approximation

0 10 20 30 40 50

Zeit t

0.0000

0.0005

0.0010

0.0015

0.0020

0.0025

0.0030

0.0035

0.0040

‖Ψε−

Ψε BO‖ L

2

ε = 14

ε = 18

ε = 116

Abbildung 4.15: Fehlerkurven der Born-Oppenheimer-Approximation zweiterOrdnung bei verschiedenen Werten für den Parameter ε.

beiden Kurven würde ein ähnlich stark oszillierendes Bild ergeben. Halbieren wir ε, ver-kleinert sich der Fehler um den Faktor vier, oder anders ausgedrückt, er verändert sichum den Faktor (1/2)2. Dies entspricht dem Verhalten der in Abschnitt 3.2 beschriebe-nen Theorie bezüglich einer oberen Schranke. Der Fehler scheint sich also, zumindest inunserem Beispielmodell, ungefähr so wie die Schranke zu verhalten und sich nicht nurasymptotisch ihr anzunähern. Aus diesem Grund interpolieren wir in Abbildung 4.16die Fehlerkurven entsprechend der Gestalt der Schranke aus (3.32) durch eine Geradeder Form ε2 (C1 + C2 · t). Den Nullpunkt lassen wir dabei außer Acht. Die Interpolationliefert uns folgende Koeffizienten:

ε Interpolationsgerade C1 C2

14 1,918 · 10−3 + 2,902 · 10−5 · t 0,03069 4,642 · 10−4

18 4,459 · 10−4 + 1,083 · 10−6 · t 0,02854 6,933 · 10−5

116 1,029 · 10−4 + 1,464 · 10−6 · t 0,02634 3,747 · 10−4

Tabelle 4.2: Koeffizienten der Interpolation.

Die statische zeitunabhängige Komponente ist hier also gegeben durch C1 und derzeitabhängige Teil wird durch C2 beschrieben. Das Ergebnis der Interpolation decktsich mit den Beobachtungen aus den vorangegangenen Grafiken. Die Konstante C2

61

Page 70: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 4. Numerische Simulation

0 10 20 30 40 50

Zeit t

0.0000

0.0005

0.0010

0.0015

0.0020

0.0025

0.0030

0.0035

0.0040

‖Ψε−

Ψε BO‖ L

2

ε = 14

ε = 18

ε = 116

Abbildung 4.16: Interpolationsgeraden zu den Fehlerkurven der Born-Oppen-heimer-Approximation zweiter Ordnung mit den Koeffizientenaus Tabelle 4.2.

ist in der Tat um ein Vielfaches kleiner als C1. Der Fehler steigt im Verlauf der Zeitdadurch vergleichsweise langsam an und vermittelt daher den konstanten Eindruck inden Abbildungen.

Zum Abschluss dieses Kapitels können wir also festhalten, dass auch in unserem Bei-spielmodell nachweislich ein signifikanter Mehrwert durch die Verbesserung der Born-Oppenheimer-Approximation um eine weitere Ordnung entsteht. Je länger die zu appro-ximierende Gesamtzeit ist, desto mehr fallen der Korrekturterm und die Verwendungdes superadiabatischen Unterraums ins Gewicht und vergrößern den Qualitätsunter-schied der beiden Ordnungen zunehmend.

62

Page 71: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

5 Fazit und kurzer Ausblick

Eine numerische Berechnung der Zeitentwicklung einer Wellenfunktion direkt in ihremmehrdimensionalen Gesamtsystem ist aufgrund des hohen Rechenaufwands und der da-mit verbundenen langen Laufzeit gerade bei komplexeren Modellen in der Praxis meistunbrauchbar. Die Born-Oppenheimer-Approximation liefert eine vielversprechende Al-ternative dazu. Ziel dieser Arbeit war, die adiabatische Störungstheorie und die daraufaufbauende Theorie der Born-Oppenheimer-Approximation einzuführen und anschlie-ßend Korrekturen zu präsentieren, die die Approximation um eine weitere Ordnung ver-bessern. Eine Gegenüberstellung der Approximationen verschiedener Ordnungen solltedie Qualität sowie den gewonnenen Mehrwert ausleuchten.

In Kapitel 2 stand das Schaffen der nötigen mathematischen Grundlagen im Vorder-grund. Wir haben den Rahmen der zeitadiabatischen Störungstheorie abgesteckt sowiedie Bedeutsamkeit der Lückenbedingung erklärt. Dadurch dass ein Teil des Spektrumsdes Hamiltonoperators durch eine Lücke vom restlichen Spektrum getrennt ist, konn-ten wir den Begriff der adiabatischen Entkopplung einführen sowie das zeitadiabatischeTheorem formulieren und beweisen. Als nächstes haben wir uns der raumadiabatischenTheorie zugewandt und Quantensysteme in diesem Setting betrachtet. Wir haben ge-lernt, dass eine Skalierung des Konfigurationsraums zu einer Entwicklung des Systemsin schnellen und langsamen Freiheitsgraden führt. Dies ermöglicht eine adiabatischeEntkopplung. Wichtige Voraussetzung dafür ist auch hier die Existenz eines durch ei-ne Lücke vom Rest isolierten Teil des Spektrums. Mittels der adiabatischen Entkopp-lung gelingt es, Freiheitsgrade zu finden, die sich langsamer mit der Zeit entwickeln alsandere und diese dadurch voneinander zu trennen. Darauf aufbauend haben wir dasraumadiabatische Theorem vorgestellt, bevor wir abschließend eine effektive Schrödin-gergleichung abgeleitet haben, die einfacher zu lösen ist als die Schrödingergleichungdes Gesamtsystems.

Die zeitabhängige Born-Oppenheimer-Approximation stand in Kapitel 3 im Mit-telpunkt. Sie nutzt die adiabatische Theorie, um die Schrödingergleichung von Syste-men aus mehreren Teilchen zu vereinfachen. Indem Freiheitsgrade, die sich auf unter-schiedlichen Zeitskalen entwickeln, voneinander getrennt werden, kann die Dynamikder Teilchen durch unterschiedliche Schrödingergleichungen beschrieben werden. DieseGleichungen können dann separat voneinander gelöst werden. Nachdem wir die adia-batische Theorie auf ein molekulares Modell angewandt haben, konnten wir durch dieBetrachtung der konventionellen zeitabhängigen Born-Oppenheimer-Approximation dieeffektive Dynamik näher untersuchen und den Born-Oppenheimer-Hamiltonoperatorerster Ordnung herleiten. Darüber hinaus haben wir ein systematisches Schema fürdie Korrektur in höheren Ordnungen angegeben und eine elementare Herleitung des

Page 72: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Kapitel 5. Fazit und kurzer Ausblick

zeitabhängigen Born-Oppenheimer-Hamiltonoperators zweiter Ordnung präsentiert.Die Durchführung der Born-Oppenheimer-Approximation an einem konkreten Bei-

spielmodell in Kapitel 4 hat den Mehrwert der Approximation im Allgemeinen sowieinsbesondere den der Approximation zweiter Ordnung noch einmal verdeutlicht. Selbstin dem zweidimensionalen Anwendungsbeispiel wurde deutlich, dass eine direkte Be-rechnung der Wellenfunktion in ihrem Konfigurationsraum impraktikabel und unhand-lich ist. Die Laufzeiten für die Berechnungen waren um ein Vielfaches höher als dieder Born-Oppenheimer-Approximation. Das macht sie zu einem mächtigen Werkzeug,das durch die vorgestellten Korrekturen noch zusätzlich enorm an Genauigkeit gewinnt.Trotz dieser Genauigkeit, die den Fehler auch über längere Zeiten hinweg klein hält, hältsich der zusätzliche Rechenaufwand im Rahmen und die Laufzeiten steigen gegenüberder konventionellen Born-Oppenheimer-Approximation vergleichsweise wenig an.

Abschließend können wir festhalten, dass unsere Ergebnisse in dieser Arbeit die Be-deutung der Born-Oppenheimer-Approximation unterstreichen und eine Verbesserungder Ordnung der konventionellen Approximation, wenn möglich, vorgezogen werden soll-te. In [Mát18], einem aktuellen Forschungsergebnis aus dem Bereich der Chemie, wurdeder Born-Oppenheimer-Hamiltonoperator zweiter Ordnung bereits unter Berücksichti-gung der ortsabhängigen Masse, das heißt inklusive des vorgestellten Korrekturterms,erfolgreich auf die stationäre Schrödingergleichung eines 4He+

2 -Ions angewandt.

Ausblick

Auf einige weiterführende Punkte sind wir bisher gar nicht oder nur knapp eingegangen,welche wir abschließend nicht unbeachtet lassen und kurz ansprechen wollen.

Ein sich zur Validierung an Kapitel 4 anschließendes Modell wäre folgendes Szenario:Das Teilchen wird vom Potential nicht festgehalten, sondern durchläuft den kritischenBereich nur einmal. Das Energieniveau ist für sehr kleine und sehr große Zeiten kon-stant. Hier stimmen der adiabatische und der superadiabatische Unterraum dann fastüberein. Eine numerische Simulation sollte daher ergeben, dass für Anfangs- und End-daten außerhalb des kritischen Bereichs die Qualität des Ergebnisses der Approximationkaum von der Verwendung von U ε ∗(1) abhängt.

In dem von uns verwendeten Beispielmodell haben wir als Hamiltonoperator derschnellen Freiheitsgrade einen harmonischen Oszillator mit unendlich vielen Eigenwer-ten verwendet. Eine alternative und numerisch einfachere Fragestellung wäre nun: Waspassiert, wenn wir diesen durch eine Matrix mit nur endlich vielen Eigenwerten ersetzenund je nach Wahl der Matrix die Breite der Spektrallücke variiert? Inwiefern ändert sichdas Ergebnis der Approximation in Abhängigkeit von der Größe der Lücke? Ein Ansatzfür weitere Untersuchungen wäre also, die naheliegende Vermutung zu überprüfen, obdie Korrekturterme umso mehr an Bedeutung gewinnen, je kleiner die Lücke ist.

64

Page 73: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

A Appendix

A.1 Hermitesche PolynomeDie Hermiteschen Polynome sind sowohl in der Mathematik als auch in der Physik weitverbreitet und in vielen Lehrbüchern zu finden. Im Folgenden soll ein kurzer Überblicküber die Definition sowie auch einige Eigenschaften gegeben werden. Zu diesem Zweckhaben wir [Wei] herangezogen.

Die Hermiteschen Polynome sind für n ∈ N Polynome vom Grad n und definiertdurch

Hn(x) := (−1)nex2 dn

dxn e−x2.

Sie lösen die Hermitesche Differentialgleichung, eine linearen Differentialgleichung zwei-ter Ordnung:

H ′′n(x)− 2x ·H ′n(x) + 2n ·Hn(x) = 0.

Des Weiteren genügen sie der rekursiven Darstellung

Hn+1(x) = 2xHn(x)− 2nHn−1(x),H ′n(x) = 2nHn−1(x).

Die ersten sechs Hermiteschen Polynome lauten also

H0(x) = 1,H1(x) = 2x,H2(x) = 4x2 − 2,H3(x) = 8x3 − 12x,H4(x) = 16x4 − 48x2 + 12,H5(x) = 32x5 − 160x3 + 120x.

Bezüglich der Gewichtsfunktion e−x2 ist die Orthogonalitätseigenschaft∫ +∞

−∞Hm(x)Hn(x)e−x2 dx = 2n n!

√π δmn

erfüllt. Für Funktionen der Form

fn(x) =√ √

ω

2n n!π12Hn(√ωx)e−

ω2 x

2

Page 74: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Anhang A. Appendix

gilt außerdem ∫ +∞

−∞fm(x)fn(x) dx = δmn, (A.1)

∫ +∞

−∞fm(x)xfn(x) dx =

n2ω , falls m = n− 1,√n+12ω , falls m = n+ 1,

0, sonst,

sowie

∫ +∞

−∞fm(x)x2fn(x) dx =

√n(n−1)2ω , falls m = n− 2,

2n+12ω , falls m = n,√

(n+1)(n+2)2ω , falls m = n+ 2,

0, sonst.

(A.2)

66

Page 75: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

A.2. Quelltext

A.2 Quelltext

Der nachfolgende Code wurde in Python 3.5.2 geschrieben. Zusätzlich wurden die Mo-dule Matplotlib und NumPy in den Versionen 1.5.1 beziehungsweise 1.11.0 verwendet.Der vollständige Code ist digital verfügbar unter https://github.com/dnlwbr/Born-Oppenheimer-Approximation.git.

A.2.1 run.py

1 " " "2 Hauptprogramm .34 Copyright ( c ) 2018 Danie l Weber .5 D i s t r i bu ted under the MIT License .6 ( See accompanying LICENSE f i l e or copy at https : // github . com/dnlwbr/Born−

Oppenheimer−Approximation/ blob / master /LICENSE)7 " " "89 import numpy as np

10 from parameter import eps , xrange , yrange , x , y11 from function_data_base import p o t e n t i a l , p lo t_potent ia l ,

plot_potentialHSL , plot_potent ia l3D , plot_potential_wave , \12 plot_HSL_wave , potential_VBH , correct ion_term , plot_potential_VBH ,

plot_correct ion_term , \13 mult iple_propagation_times , compare_grid , plot_grid ,

compare_spl i t t ing_steps , p l o t _ s p l i t t i n g , compare_error_eps , \14 plot_error_eps , compare_error_time , plot_error_time ,

compare_error_time_VMU , plot_error_time_VMU , ask_data_handling151617 r e t r i e s = 318 whi le True :1920 todo = s t r ( input ( " " "21 Was s o l l berechnet werden?22 [ 1 ] Konturplot des P o t e n t i a l s23 [ 2 ] HSL−Plot des P o t e n t i a l s24 [ 3 ] 3D−Plot des P o t e n t i a l s25 [ 4 ] Konturplot der Wel lenfunkt ion im P o t e n t i a l26 [ 5 ] HSL−Plot der Wel lenfunkt ion27 [ 6 ] Plot Born−Huang−P o t e n t i a l28 [ 7 ] Plot Korrekturterm29 [ 8 ] Konvergenztest − G i t t e r g r öße30 [ 9 ] Konvergenztest − Z e i t s c h r i t t w e i t e beim Strang−S p l i t t i n g31 [ 1 0 ] Qual i t ä t s t e s t nach e p s i l o n32 [ 1 1 ] Qual i t ä t s t e s t nach Ze i t33 [ 1 2 ] Qual i t ä t s t e s t nach Ze i t ( nur 2 . Ordnung )34 [Q] Beenden35 >> " " " ) )3637 i f todo == " 1 " :38 # Konturplot des P o t e n t i a l s39 p r i n t ( " P o t e n t i a l wird berechnet . . . " )40 V = p o t e n t i a l (x , y )

67

Page 76: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Anhang A. Appendix

4142 p r i n t ( " Graph wird gez e i chne t . . . " )43 p l o t_poten t i a l (x , y , V)44 break4546 e l i f todo == " 2 " :47 # HSL−Plot des P o t e n t i a l s48 p r i n t ( " P o t e n t i a l wird berechnet . . . " )49 V = p o t e n t i a l (x , y )5051 p r i n t ( " Graph wird gez e i chne t . . . " )52 plot_potentialHSL (V, xrange , yrange )53 break5455 e l i f todo == " 3 " :56 # 3D−Plot des P o t e n t i a l s57 p r i n t ( " P o t e n t i a l wird berechnet . . . " )58 V = p o t e n t i a l (x , y )5960 p r i n t ( " Graph wird gez e i chne t . . . " )61 plot_potent ia l3D (x , y , V)62 break6364 e l i f todo == " 4 " :65 # Konturplot der Wel lenfunkt ion im P o t e n t i a l66 p r i n t ( " P o t e n t i a l wird berechnet . . . " )67 V = p o t e n t i a l (x , y )6869 p r i n t ( " Wel lenfunkt ion wird berechnet . . . " )70 t imes = [ t f o r t in range (0 , 6 , 1) ]71 a l l P s i = mult iple_propagation_times (x , y , t imes )7273 p r i n t ( " Graph wird gez e i chne t . . . " )74 plot_potential_wave (x , y , V, a l l P s i )75 break7677 e l i f todo == " 5 " :78 # HSL−Plot der Wel lenfunkt ion79 p r i n t ( " Wel lenfunkt ion wird berechnet . . . " )80 t imes = [ t f o r t in range (0 , 6 , 1) ]81 a l l P s i = mult iple_propagation_times (x , y , t imes )8283 p r i n t ( " Graph wird gez e i chne t . . . " )84 plot_HSL_wave ( a l l P s i , xrange , yrange )85 break8687 e l i f todo == " 6 " :88 # Plot Born−Huang−P o t e n t i a l89 p r i n t ( " Born−Huang−P o t e n t i a l wird berechnet . . . " )90 VBH = potential_VBH ( x )9192 p r i n t ( " Graph wird gez e i chne t . . . " )93 plot_potential_VBH (x , VBH)94 break9596 e l i f todo == " 7 " :

68

Page 77: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

A.2. Quelltext

97 # Plot Korrekturterm98 p r i n t ( " Korrekturterm wird berechnet . . . " )99 M = correct ion_term (x , eps )

100101 p r i n t ( " Graph wird gez e i chne t . . . " )102 plot_correct ion_term (x , M)103 break104105 e l i f todo == " 8 " :106 # Konvergenztest − G it t e rg r öße107 t imes = np . l i n s p a c e (0 , 12 , 25)108 nx1 , nx2 , nx3 , nx4 , nx5 , ny = 6 , 7 , 8 , 9 , 10 , None109 # nx1 , nx2 , nx3 , nx4 , nx5 , ny = 9 , 9 , 9 , 9 , 9 , [ 6 , 7 , 8 , 9 , 10 ]110111 handl ing = ask_data_handling ( )112 i f handl ing i s None :113 cont inue114115 p r i n t ( " \ nWel lenfunktionen werden berechnet . . . " )116 normDiff = compare_grid ( xrange , yrange , times , nx1 , nx2 , nx3 , nx4 ,

nx5 , ny=ny ,117 potV=p o t e n t i a l , e f f e k t i v=False , handl ing=

handl ing )118119 p r i n t ( " e f f e k t i v e Wel lenfunkt ionen werden berechnet . . . " )120 normDif f_ef f = compare_grid ( xrange , yrange , times , nx1 , nx2 , nx3 ,

nx4 , nx5 , ny=ny )121122 p r i n t ( " Graph wird gez e i chne t . . . " )123 p lot_gr id ( normDiff , normDiff_eff , times , nx1 , nx2 , nx3 , nx4 , nx5 ,

ny=ny )124 break125126 e l i f todo == " 9 " :127 # Konvergenztest − Z e i t s c h r i t t w e i t e beim Strang−S p l i t t i n g128 t imes = np . l i n s p a c e (0 , 12 , 25)129 h1 , h2 , h3 , h4 , h5 = 1/40000 , 1/50000 , 1/60000 , 1/70000 , 1/80000130131 handl ing = ask_data_handling ( )132 i f handl ing i s None :133 cont inue134135 p r i n t ( " \ nPotent i a l wird berechnet . . . " )136 V = p o t e n t i a l (x , y )137138 p r i n t ( " Wel lenfunkt ionen werden berechnet . . . " )139 normDiff = compare_spl i t t ing_steps (x , y , times , h1 , h2 , h3 , h4 , h5

, potV=V, e f f e k t i v=False , handl ing=handl ing )140141 p r i n t ( " e f f e k t i v e Wel lenfunkt ionen werden berechnet . . . " )142 normDif f_ef f = compare_spl i t t ing_steps (x , y , times , h1 , h2 , h3 , h4

, h5 )143144 p r i n t ( " Graph wird gez e i chne t . . . " )145 p l o t _ s p l i t t i n g ( normDiff , normDiff_eff , times , h1 , h2 , h3 , h4 , h5 )146 break

69

Page 78: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Anhang A. Appendix

147148 e l i f todo == " 10 " :149 # Qual i t ä t s t e s t nach e p s i l o n150 t = 3151 e p s i l o n s = [2∗∗( −1) , 2∗∗( −1.5) , 2∗∗( −2) , 2∗∗( −2.5) , 2∗∗( −3) ,

2∗∗( −3.5) , 2∗∗( −4) ]152153 handl ing = ask_data_handling ( )154 i f handl ing i s None :155 cont inue156157 p r i n t ( " \ nPotent i a l wird berechnet . . . " )158 V = p o t e n t i a l (x , y )159160 p r i n t ( " Wel lenfunkt ionen werden berechnet . . . " )161 normDiff_eff , normDiff_effV , normDiff_effVM , normDiff_effVMU =

compare_error_eps (x , y , t , e p s i l o n s , V, handl ing )162163 p r i n t ( " Graph wird gez e i chne t . . . " )164 plot_error_eps ( normDiff_eff , normDiff_effV , normDiff_effVM ,

normDiff_effVMU , t , e p s i l o n s )165 break166167 e l i f todo == " 11 " :168 # Qual i t ä t s t e s t nach Ze i t169 t imes = np . l i n s p a c e (0 , 12 , 25)170171 handl ing = ask_data_handling ( )172 i f handl ing i s None :173 cont inue174175 p r i n t ( " \ nPotent i a l wird berechnet . . . " )176 V = p o t e n t i a l (x , y )177178 p r i n t ( " Wel lenfunkt ionen werden berechnet . . . " )179 normDiff_eff , normDiff_effV , normDiff_effVM , normDiff_effVMU =

compare_error_time (x , y , times , V, handl ing )180181 p r i n t ( " Graph wird gez e i chne t . . . " )182 plot_error_time ( times , normDiff_eff , normDiff_effV , normDiff_effVM

, normDiff_effVMU )183 break184185 e l i f todo == " 12 " :186 # Qual i t ä t s t e s t nach Ze i t ( nur 2 . Ordnung )187 t imes = np . l i n s p a c e (0 , 50 , 101)188 e p s i l o n s = [1/4 , 1/8 , 1/16 ]189190 handl ing = ask_data_handling ( )191 i f handl ing i s None :192 cont inue193194 p r i n t ( " \ nPotent i a l wird berechnet . . . " )195 V = p o t e n t i a l (x , y )196197 p r i n t ( " Wel lenfunkt ionen werden berechnet . . . " )

70

Page 79: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

A.2. Quelltext

198 normDiff_effVMU_eps = compare_error_time_VMU (x , y , times , e p s i l o n s, V, handl ing )

199200 p r i n t ( " Graph wird gez e i chne t . . . " )201 plot_error_time_VMU ( times , normDiff_effVMU_eps , i n t e r p o l a t e=True )202 break203204 e l i f todo in ( "Q" , " q " ) :205 p r i n t ( " Beendet . " )206 break207208 e l s e :209 p r i n t ( "Ungü l t i g e Eingabe ! " )210211 r e t r i e s −= 1212 i f r e t r i e s == 0 :213 p r i n t ( " Abbruch aufgrund von d r e i f e h l e r h a f t e n Versuchen . " )214 break

A.2.2 parameter.py

1 " " "2 Zu verwendende Parameter .34 Die h i e r angegebenen Parameter werden von den Funktionen verwendet , f a l l s

ke ine anderen Fest legungen in dem5 Hauptskr ipt vorgenommen wurden .67 Copyright ( c ) 2018 Danie l Weber .8 D i s t r i bu ted under the MIT License .9 ( See accompanying LICENSE f i l e or copy at https : // github . com/dnlwbr/Born−

Oppenheimer−Approximation/ blob / master /LICENSE)10 " " "1112 import numpy as np131415 # a d i a b a t i s c h e r Parameter16 eps = 1/81718 # r e d u z i e r t e s Planck ’ s ches Wirkungsquantum19 h_red = 12021 # Auswertungsbereich22 Nx = 2 ∗∗ 9 # Anzahl der Gitterpunkte

in x−Richtung23 Ny = 2 ∗∗ 8 # Anzahl der Gitterpunkte

in y−Richtung24 xrange = 5 # I n t e r v a l l g r e n z e x−

Richtung25 yrange = 4 # I n t e r v a l l g r e n z e y−

Richtung26 x = np . arange(−xrange , xrange , 2 ∗ xrange / Nx) # I n t e r v a l l x−Richtung27 y = np . arange(−yrange , yrange , 2 ∗ yrange / Ny) # I n t e r v a l l y−Richtung2829 # Anfangswerte

71

Page 80: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Anhang A. Appendix

30 px = 1 # Anfangsimpuls in x−Richtung

3132 # Tei lchenmasse33 m = 13435 # S c h r i t t w e i t e im Strang−S p l i t t i n g36 h = 1/6000037 # I s t e in Zeitpunkt n i cht e r r e i chbar , wird automatisch d i e nächstmö g l i c h e

k l e i n e r e S c h r i t t w e i t e verwendet

72

Page 81: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

A.2. Quelltext

A.2.3 function_data_base.py

1 " " "2 Sammlung von Funktionen .34 Copyright ( c ) 2018 Danie l Weber .5 D i s t r i bu ted under the MIT License .6 ( See accompanying LICENSE f i l e or copy at https : // github . com/dnlwbr/Born−

Oppenheimer−Approximation/ blob / master /LICENSE)7 " " "89 import os

10 import p i c k l e11 import sys12 import time13 from c o l o r s y s import hls_to_rgb14 from f r a c t i o n s import Fract ion1516 import matp lo t l i b as mpl17 import matp lo t l i b . pyplot as p l t18 import mpl_too lk i t s . mplot3d19 import numpy as np20 from numpy import exp , s q r t21 from numpy . f f t import f f t , i f f t , f f t 2 , i f f t 2 , f f t f r e q2223 from parameter import eps , h_red , px , m, h242526 # Frequenz27 de f f requency ( x ) :28 " " " Berechnet d i e Funktionswerte der Frequenz . " " "29 r e s u l t = s q r t (2 ∗ (1 + x ∗∗4) )30 re turn r e s u l t313233 # Ableitung der Frequenz34 de f f requency_der ived ( x ) :35 " " " Berechnet d i e Funktionswerte der Able itung der Frequenz . " " "36 r e s u l t = 2 ∗ s q r t (2 ) ∗ x∗∗3 / s q r t (1 + x ∗∗ 4)37 return r e s u l t383940 # Anfangswert der Wel lenfunkt ion der Kerne (x−Richtung )41 de f calc_psi_0 (x , e p s i l o n=eps ) :42 " " " Berechnet d i e Funktionswerte von psi_0 f ür e in gegebenen Vektor x .

" " "43 psi_0 = (1 / (2∗np . p i ) ∗∗(1/4) ) ∗ exp(−x∗∗2 / 2) ∗ exp (1 j ∗ x ∗ px /

e p s i l o n )44 re turn psi_0454647 # Grundzustand der Eigenfunkt ionen der Elektronen48 de f calc_phi_0 (x , y ) :49 " " " Berechnet d i e Funktionswerte von phi_0 f ür zwei gegebene Vektoren x

und y . " " "50 Nx, Ny = len ( x ) , l en ( y )

73

Page 82: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Anhang A. Appendix

51 phi_0 = np . empty ( (Ny, Nx) )52 w = frequency ( x )53 f o r j in range (Ny) :54 f o r i in range (Nx) :55 phi_0 [ j , i ] = (w[ i ] / np . p i ) ∗∗(1/4) ∗ exp(−w[ i ] / 2 ∗ y [ j ] ∗ ∗ 2 )5657 return phi_0585960 # Eigenfunkt ion phi_261 de f calc_phi_2 (x , y ) :62 " " " Berechnet d i e Funktionswerte von phi_2 f ür zwei gegebene Vektoren x

und y . " " "63 Nx, Ny = len ( x ) , l en ( y )64 phi_2 = np . empty ( (Ny, Nx) )65 w = frequency ( x )66 f o r j in range (Ny) :67 f o r i in range (Nx) :68 phi_2 [ j , i ] = ( s q r t (w[ i ] ) / (8 ∗ s q r t (np . p i ) ) ) ∗∗(1/2) ∗ (4 ∗ w

[ i ] ∗ y [ j ]∗∗2 − 2) ∗ exp(−w[ i ] / 2 ∗ y [ j ] ∗ ∗ 2 )6970 return phi_2717273 # Wel lenfunkt ion im 2D−Gesamtsystem74 de f ca lc_Psi (x , y , ps i , phi_0 , useU , e p s i l o n=eps ) :75 " " " Berechnet d i e Funktionswerte der zwe id imens iona len Wel lenfunkt ion

aus den gegebenen76 Funktionswerten von p s i und phi . " " "77 Ny, Nx = np . shape ( phi_0 )78 Psi = np . empty ( (Ny, Nx) , dtype=complex )7980 i f useU i s True :81 phi_2 = calc_phi_2 (x , y )82 w = frequency ( x )83 w_derived = frequency_der ived ( x )84 k = 2 ∗ np . p i ∗ f f t f r e q (Nx, 2 ∗ np . abs ( x [ 0 ] ) / Nx)85 ps i_der ived = i f f t (1 j ∗ k ∗ f f t ( p s i ) )8687 f o r j in range (Ny) :88 f o r i in range (Nx) :89 Psi [ j , i ] = p s i [ i ] ∗ phi_0 [ j , i ] \90 − e p s i l o n ∗∗2 ∗ w_derived [ i ] / ( 2 ∗ ∗ ( 5 / 2 ) ∗ w[ i

] ∗ ∗ 2 ) ∗ phi_2 [ j , i ] ∗ ps i_der ived [ i ]9192 e l s e :93 f o r j in range (Ny) :94 f o r i in range (Nx) :95 Psi [ j , i ] = p s i [ i ] ∗ phi_0 [ j , i ]9697 return Psi9899

100 # P o t e n t i a l101 de f p o t e n t i a l (x , y ) :102 " " " Berechnet das P o t e n t i a l aus den Vektoren x und y . " " "

74

Page 83: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

A.2. Quelltext

103 Nx, Ny = len ( x ) , l en ( y )104 V = np . empty ( (Ny, Nx) )105 w = frequency ( x )106 f o r j in range (Ny) :107 f o r i in range (Nx) :108 V[ j , i ] = −0.5 ∗ x [ i ]∗∗2 + 0 .5 ∗ w[ i ]∗∗2 ∗ y [ j ]∗∗2109110 return V111112113 # P o t e n t i a l oben abschneiden114 de f potent ia l_cut (V, c ) :115 " " " Schne idet a l l e Werte des gegebenen P o t e n t i a l s oberhalb von zwei ab .

" " "116 Vy, Vx = V. shape117 Vcut = np . copy (V)118 f o r j in range (Vy) :119 f o r i in range (Vx) :120 i f c < Vcut [ j , i ] :121 Vcut [ j , i ] = None122123 return Vcut124125126 # Konturplot des P o t e n t i a l s127 de f p l o t_pot en t i a l (x , y , V) :128 " " " E r s t e l l t e in Konturplot e i n e s gegebenen P o t e n t i a l s . " " "129 Vcut = potent ia l_cut (V, 4)130131 p l t . f i g u r e ( )132 plot_norm = mpl . c o l o r s . SymLogNorm( l i n t h r e s h =1.5 , vmin=−20, vmax=4)133 p l t . contour f (x , y , Vcut , 500 , cmap=p l t . cm . get_cmap ( ’ j e t ’ ) , norm=

plot_norm , extend=’max ’ )134 p l t . c o l o r b a r ( o r i e n t a t i o n=’ h o r i z o n t a l ’ , format=’% . 0 f ’ , extend=’max ’ ,

ex t end f rac =0.02 , shr ink =0.85)135 p l t . y t i c k s (np . l i n s p a c e (−4 , 4 , 5) )136 p l t . t i t l e ( ’ Konturplot von $V$ ’ )137 p l t . x l a b e l ( ’ Long i tud ina l e $x$ ’ )138 p l t . y l a b e l ( ’ Transver sa l e $y$ ’ )139 p l t . show ( )140141142 # HLS−Transformation143 de f complex_to_HSL ( z ) :144 " " " Trans formier t komplexe Funktionswerte in den HSL−Farbraum . " " "145 H = np . ang le ( z )146 L = (1 − 2 ∗∗ (−np . abs ( z ) ) ) ∗ 1147 S = 1148149 c = np . v e c t o r i z e ( hls_to_rgb ) (H, L , S)150 c = np . array ( c ) # (3 , n ,m) zu (n ,m, 3 ) zu (m, n , 3 )151 c = c . swapaxes (0 , 2)152 c = c . swapaxes (0 , 1)153154 return c155

75

Page 84: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Anhang A. Appendix

156157 # HSL−Plot des P o t e n t i a l s158 de f plot_potentialHSL (V, xrange , yrange ) :159 " " " E r s t e l l t e in HSL−Plot e i n e s gegebenen P o t e n t i a l s . " " "160 p l t . f i g u r e ( )161 img = complex_to_HSL (V)162 p l t . imshow ( img , aspect =0.68 , extent=[−xrange , xrange , −yrange , yrange

] )163 p l t . y t i c k s (np . l i n s p a c e (−4 , 4 , 5) )164 p l t . t i t l e ( ’HSL−Plot von $V$ ’ )165 p l t . x l a b e l ( ’ Long i tud ina l e $x$ ’ )166 p l t . y l a b e l ( ’ Transver sa l e $y$ ’ )167 p l t . show ( )168169170 # Plot des P o t e n t i a l s171 de f p lot_potent ia l3D (x , y , V) :172 " " " E r s t e l l t e in 3D−Plot e i n e s gegebenen P o t e n t i a l s . " " "173 Vcut = potent ia l_cut (V, 10)174175 X, Y = np . meshgrid (x , y )176 f i g = p l t . f i g u r e ( )177 ax = f i g . gca ( p r o j e c t i o n=’ 3d ’ )178 ax . plot_wireframe (X, Y, Vcut , r s t r i d e =32, c s t r i d e =30, alpha =0.5)179 zmin = Vcut [ l en ( Vcut ) //2 , −1]180 ax . contour (X, Y, Vcut , z d i r=’ z ’ , o f f s e t=zmin , alpha =0.6)181 p l t . y t i c k s (np . l i n s p a c e (−4 , 4 , 5) )182 p l t . t i t l e ( ’ 3D−Plot von $V$ ’ )183 ax . s e t_x labe l ( ’ Long i tud ina l e $x$ ’ )184 ax . s e t_y labe l ( ’ Transver sa l e $y$ ’ )185 ax . s e t _ z l a b e l ( ’$V(x , y ) $ ’ )186 p l t . show ( )187188189 # Konturplot der Wel lenfunkt ion im P o t e n t i a l190 de f plot_potential_wave (x , y , V, a l l P s i ) :191 " " " E r s t e l l t e i n z e l n e Konturplots e i n e r Wel lenfunkt ion zu ver sch i edenen

Zeitpunkten in einem gegebenen P o t e n t i a l . " " "192 p l t . f i g u r e ( f i g s i z e =(8 , 8) )193 l v l s = len ( a l l P s i )194195 i f l v l s == 1 :196 Psi = np . abs ( a l l P s i [ 0 ] [ 1 ] ) ∗∗ 2197 p l t . contour (x , y , V, l e v e l s=np . arange ( 0 . 2 , 2 , 0 . 3 ) , cmap=p l t . cm .

get_cmap ( ’YlOrRd ’ ) )198 p l t . contour (x , y , Psi , cmap=p l t . cm . get_cmap ( ’ b inary ’ ) )199 p l t . yl im (−3 , 3)200 p l t . t i t l e ( ’ $t = {}$ ’ . format ( a l l P s i [ 0 ] [ 0 ] ) )201 p l t . x l a b e l ( ’ Long i tud ina l e $x$ ’ )202 p l t . y l a b e l ( ’ Transver sa l e $y$ ’ )203 e l s e :204 z = l v l s // 2 + ( l v l s % 2 > 0)205 f o r i in range (1 , l v l s +1) :206 p l t . subp lot ( z , 2 , i )207 p l t . contour (x , y , V, l e v e l s=np . arange ( 0 . 2 , 2 , 0 . 3 ) , cmap=p l t .

cm . get_cmap ( ’YlOrRd ’ ) )

76

Page 85: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

A.2. Quelltext

208 plot_data = np . abs ( a l l P s i [ i − 1 ] [ 1 ] ) ∗∗ 2209 p l t . contour (x , y , plot_data , cmap=p l t . cm . get_cmap ( ’ b inary ’ ) )210 p l t . yl im (−2 , 2)211 p l t . y t i c k s (np . l i n s p a c e (−2 , 2 , 5) )212 p l t . t i t l e ( ’ $t = {}$ ’ . format ( a l l P s i [ i − 1 ] [ 0 ] ) )213 i f i == l v l s :214 p l t . x l a b e l ( ’ Long i tud ina l e $x$ ’ )215 i f ( l v l s % 2 == 0) and ( i == l v l s −1) :216 p l t . x l a b e l ( ’ Long i tud ina l e $x$ ’ )217 i f i % 2 != 0 :218 p l t . y l a b e l ( ’ Transver sa l e $y$ ’ )219220 p l t . t ight_layout ( )221 p l t . show ( )222223224 # HSL−Plot der Wel lenfunkt ion225 de f plot_HSL_wave ( a l l P s i , xrange , yrange ) :226 " " " E r s t e l l t e i n z e l n e HSL−Plot s e i n e r Wel lenfunkt ionen zu gegebenen

Zeitpunkten . " " "227 p l t . f i g u r e ( f i g s i z e =(8 , 8) )228 l v l s = len ( a l l P s i )229230 i f l v l s == 1 :231 img = complex_to_HSL ( a l l P s i [ 0 ] [ 1 ] )232 p l t . imshow ( img , aspect=’ auto ’ , extent=[−xrange , xrange , −yrange ,

yrange ] )233 p l t . yl im (−3 , 3)234 p l t . t i t l e ( ’ $t = {}$ ’ . format ( a l l P s i [ 0 ] [ 0 ] ) )235 p l t . x l a b e l ( ’ Long i tud ina l e $x$ ’ )236 p l t . y l a b e l ( ’ Transver sa l e $y$ ’ )237 e l s e :238 z = l v l s // 2 + ( l v l s % 2 > 0)239 f o r i in range (1 , l v l s +1) :240 p l t . subplot ( z , 2 , i )241 img = complex_to_HSL ( a l l P s i [ i − 1 ] [ 1 ] )242 p l t . imshow ( img , aspect=’ auto ’ , extent=[−xrange , xrange , −

yrange , yrange ] )243 p l t . yl im (−3 , 3)244 p l t . y t i c k s (np . l i n s p a c e (−2 , 2 , 3) )245 p l t . t i t l e ( ’ $t = {}$ ’ . format ( a l l P s i [ i − 1 ] [ 0 ] ) )246 i f i == l v l s :247 p l t . x l a b e l ( ’ Long i tud ina l e $x$ ’ )248 i f ( l v l s % 2 == 0) and ( i == l v l s −1) :249 p l t . x l a b e l ( ’ Long i tud ina l e $x$ ’ )250 i f i % 2 != 0 :251 p l t . y l a b e l ( ’ Transver sa l e $y$ ’ )252253 p l t . t ight_layout ( )254 p l t . show ( )255256257 # Born−Huang−P o t e n t i a l258 de f potential_VBH ( x ) :259 " " " Berechnet d i e Funktionswerte des Born−Huang−P o t e n t i a l s . " " "260 r e s u l t = x∗∗6 / (2 ∗ (1 + x ∗∗4) ∗∗2)

77

Page 86: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Anhang A. Appendix

261 return r e s u l t262263264 # M−Term265 de f correct ion_term (x , e p s i l o n ) :266 " " " Berechnet d i e Funktionswerte des Korrekturterms . " " "267 r e s u l t = − e p s i l o n ∗∗2 ∗ x∗∗6 / (4 ∗ s q r t (2 ) ∗ (1 + x ∗∗4) ∗∗(5/2) )268 return r e s u l t269270271 # Plot Born−Huang−P o t e n t i a l272 de f plot_potential_VBH (x , VBH) :273 " " " P l o t t e t des Born−Huang−P o t e n t i a l s . " " "274 p l t . f i g u r e ( )275 p l t . p l o t (x , VBH)276 p l t . t i t l e ( r " Born−Huang−P o t e n t i a l $V_\mathrm{BH}$ " )277 p l t . x l a b e l ( " $x$ " )278 p l t . y l a b e l ( r "$V_\mathrm{BH}( x ) $ " )279 p l t . show ( )280281282 # Plot Born−Huang−P o t e n t i a l283 de f p lot_correct ion_term (x , M) :284 " " " P l o t t e t den Korrekturterm . " " "285 p l t . f i g u r e ( )286 p l t . p l o t (x , M)287 p l t . t i t l e ( r " $\ mathfrak {v}$−Term" )288 p l t . x l a b e l ( " $x$ " )289 p l t . y l a b e l ( r " $\ mathfrak {v }( x ) $ " )290 p l t . show ( )291292293 # a d i a b a t i s c h e Ze i t entwick lung nach Born−Oppenheimer (1d)294 de f time_propagation_BO (x , ps i , t , useVBH=True , useM=True , e p s i l o n=eps ,

t_delta=h) :295 " " " Berechnet d i e a d i a b a t i s c h e Ze i t entwick lung m i t t e l s dem Born−

Oppenheimer−Hamiltonoperator . " " "296 i f t == 0 :297 return p s i298299 E_0 = ( s q r t (2 + 2 ∗ x ∗∗ 4) − x ∗∗ 2) / 2300 i f useVBH i s True :301 VBH = potential_VBH ( x )302 e l s e :303 VBH = 0304305 i f abs ( round ( t / t_delta ) − ( t / t_delta ) ) <= 10 ∗∗ (−10) :306 s p l i t t i n g _ s t e p s = i n t ( round ( t / t_delta ) )307 e l s e :308 s p l i t t i n g _ s t e p s = i n t (np . c e i l ( t / t_delta ) )309310 k = 2 ∗ np . p i ∗ f f t f r e q ( l en ( x ) , 2∗np . abs ( x [ 0 ] ) / l en ( x ) )311 eA = exp(−1 j /h_red ∗ (E_0 + h_red ∗∗2/(2∗m) ∗ e p s i l o n ∗∗2 ∗ VBH) ∗ t /(

e p s i l o n ∗2∗ s p l i t t i n g _ s t e p s ) )312 i f useM i s True :313 eB = exp(−1 j /h_red ∗ e p s i l o n ∗∗2 ∗ h_red ∗∗2/(2∗m) ∗ k∗∗2 ∗ t /(

78

Page 87: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

A.2. Quelltext

e p s i l o n ∗2∗ s p l i t t i n g _ s t e p s ) )314 M = correct ion_term (x , e p s i l o n )315316 f o r i in range ( s p l i t t i n g _ s t e p s ) :317 p s i = eB ∗ f f t (eA ∗ p s i )318 pmp = i f f t (1 j ∗ k ∗ p s i )319 i f s t r (np . i snan (np . sum(pmp) ) ) == ’ True ’ :320 r a i s e ValueError ( ’ B i t t e k l e i n e r e S c h r i t t w e i t e wä hlen ! (

Aktue l l : h =’ + s t r (h) + ’ ) ’ )321 pmp = 1 j ∗ k ∗ f f t (1 j /h_red ∗ e p s i l o n ∗ M ∗ t /

s p l i t t i n g _ s t e p s ∗ pmp)322 p s i += pmp323 p s i = eA ∗ ( i f f t (eB ∗ p s i ) )324 e l s e :325 eB = exp(−1 j /h_red ∗ e p s i l o n ∗∗2 ∗ h_red ∗∗2/(2∗m) ∗ k∗∗2 ∗ t /(

e p s i l o n ∗ s p l i t t i n g _ s t e p s ) )326 f o r i in range ( s p l i t t i n g _ s t e p s ) :327 p s i = eA ∗ ( i f f t (eB ∗ f f t (eA ∗ p s i ) ) )328329 return p s i330331332 # Ze i tentwick lung (2d)333 de f time_propagation_2d (x , y , Psi , V, t , e p s i l o n=eps , t_delta=h) :334 " " " Berechnet d i e Ze i t entwick lung m i t t e l s dem zweid imens iona len

Hamiltonoperator . " " "335 i f t == 0 :336 return Psi337338 Nx, Ny = len ( x ) , l en ( y )339 kx_tmp = 2 ∗ np . p i ∗ f f t f r e q ( l en ( x ) , 2∗np . abs ( x [ 0 ] ) / l en ( x ) )340 ky_tmp = 2 ∗ np . p i ∗ f f t f r e q ( l en ( y ) , 2∗np . abs ( y [ 0 ] ) / l en ( y ) )341 kx = np . empty ( (Ny, Nx) )342 ky = np . empty ( (Ny, Nx) )343344 f o r j in range (Ny) :345 kx [ j , : ] = kx_tmp346 f o r i in range (Nx) :347 ky [ : , i ] = ky_tmp348349 i f abs ( round ( t / t_delta ) − ( t / t_delta ) ) <= 10 ∗∗ (−10) :350 s p l i t t i n g _ s t e p s = i n t ( round ( t / t_delta ) )351 e l s e :352 s p l i t t i n g _ s t e p s = i n t (np . c e i l ( t / t_delta ) )353354 eA = exp(−1 j /h_red ∗ ( e p s i l o n ∗∗2 ∗ h_red ∗∗2/(2∗m) ∗ kx∗∗2 + h_red

∗∗2/(2∗m) ∗ ky ∗∗2) ∗ t /( e p s i l o n ∗ s p l i t t i n g _ s t e p s ) )355 eB = exp(−1 j /h_red ∗ V ∗ t /( e p s i l o n ∗2∗ s p l i t t i n g _ s t e p s ) )356 f o r i in range ( s p l i t t i n g _ s t e p s ) :357 Psi = eB ∗ ( i f f t 2 (eA ∗ f f t 2 (eB ∗ Psi ) ) )358359 return Psi360361362 # Vektor aus Ze i tentwick lungen363 de f mult iple_propagat ion_times (x , y , times , ∗potV , e f f e k t i v=True ,

79

Page 88: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Anhang A. Appendix

364 useVBH=True , useM=True , useU=None , e p s i l o n=eps , t_delta=h) :

365 " " " Berechnet d i e Ze i tentwick lungen zu mehreren Zeitpunkten . " " "366 psi_0 = calc_psi_0 (x , e p s i l o n=e p s i l o n )367 phi_0 = calc_phi_0 (x , y )368 i f e f f e k t i v i s Fa l se :369 i f l en ( potV ) > 1 :370 p r i n t ( ’Warnung : Zu v i e l e Argumente ! Üb e r f l ü s s i g e werden

i g n o r i e r t . ’ )371 e l i f l en ( potV ) == 0 or potV [ 0 ] i s None :372 r a i s e IOError ( ’ Kein P o t e n t i a l angegeben ! ’ )373 V = potV [ 0 ]374 i f useU i s None :375 useU = False376 Psi = calc_Psi (x , y , psi_0 , phi_0 , useU , e p s i l o n=e p s i l o n )377 typ = " 2D"378 e l s e :379 i f l en ( potV ) > 0 and potV [ 0 ] i s not None :380 p r i n t ( ’Warnung : Zu v i e l e Argumente ! Üb e r f l ü s s i g e s P o t e n t i a l

wird i g n o r i e r t . ’ )381 i f useU i s None :382 i f useM i s True :383 useU = True384 e l s e :385 useU = False386 p s i = psi_0387 typ = " 1D"388389 a l l P s i = [ ]390 p r i n t ( ’ \ t [ ’ + typ + ’ ] Todo : ’ + ’ [ ’ + 100 ∗ ’− ’ + ’ ] ’ )391 sys . s tdout . wr i t e ( ’ \ t F o r t s c h r i t t : [ ’ )392 sys . s tdout . f l u s h ( )393 printedMark = 0394 s t a r t = time . time ( )395 f o r i in range ( l en ( t imes ) ) :396 i f e f f e k t i v i s Fa l se :397 i f i >= 1 :398 Psi = time_propagation_2d (x , y , Psi , V, t imes [ i ]− t imes [ i

−1] , e p s i l o n=eps i l on , t_delta=t_delta )399 e l s e :400 Psi = time_propagation_2d (x , y , Psi , V, t imes [ i ] , e p s i l o n=

eps i l on , t_delta=t_delta )401 e l s e :402 i f i >= 1 :403 p s i = time_propagation_BO (x , ps i , t imes [ i ]− t imes [ i −1] ,

useVBH , useM , e p s i l o n=eps i l on , t_delta=t_delta )404 e l s e :405 p s i = time_propagation_BO (x , ps i , t imes [ i ] , useVBH , useM ,

e p s i l o n=eps i l on , t_delta=t_delta )406 Psi = calc_Psi (x , y , ps i , phi_0 , useU , e p s i l o n=e p s i l o n )407 a l l P s i . append ( [ t imes [ i ] , Ps i ] )408 i f l en ( t imes )−1 == 0 :409 currentMark = 100410 e l s e :411 currentMark = 100 ∗ i // ( l en ( t imes ) −1)412 i f printedMark < currentMark :

80

Page 89: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

A.2. Quelltext

413 sys . s tdout . wr i t e ( ’ ∗ ’ ∗ ( currentMark − printedMark ) )414 sys . s tdout . f l u s h ( )415 printedMark = currentMark416 sys . s tdout . wr i t e ( ’ ] \ n ’ )417 end = time . time ( )418 t o t a l = [ i n t ( ( end − s t a r t ) % 60) , ( end − s t a r t ) // 60 ]

# [ s , min , (h) ]419 i f t o t a l [ 1 ] >= 60 :420 t o t a l . append ( t o t a l [ 1 ] // 60)421 t o t a l [ 1 ] = i n t ( t o t a l [ 1 ] % 60)422 p r i n t ( ’ \tBenö t i g t e Ze i t : { : 0 2 . 0 f } : { : 0 2 d } : { : 0 2 d} Stunden\n ’ . format (

t o t a l [ 2 ] , t o t a l [ 1 ] , t o t a l [ 0 ] ) )423 e l s e :424 p r i n t ( ’ \tBenö t i g t e Ze i t : { : 0 2 . 0 f } : { : 0 2 d} Minuten\n ’ . format ( t o t a l

[ 1 ] , t o t a l [ 0 ] ) )425426 return a l l P s i427428429 # Verg l e i che Ze i tentwick lungen be i u n t e r s c h i e d l i c h e n G i t t e r g r ößen430 de f compare_grid ( xrange , yrange , times , ∗nx , ny=None , potV=None ,431 e f f e k t i v=True , useVBH=True , useM=True , useU=False ,

handl ing="new" ) :432 " " " Berechnet d i e normierte D i f f e r e n z der Wel lenfunkt ionen zu gegebenen

G i t t e r g r ößen . " " "433 i f e f f e k t i v i s Fa l se and handl ing [ 0 ] == " load " :434 with open ( handl ing [ 1 ] , ’ rb ’ ) as i n p u t _ f i l e :435 p r i n t ( " \t2D−Daten werden ge laden . . . " )436 nx , ny , a l lPs i_nxy = p i c k l e . load ( i n p u t _ f i l e )437 e l s e :438 x_nx , y_ny , a l lPs i_nxy = [ ] , [ ] , [ ]439 V = None440441 i f ny i s None :442 empty = True443 ny = [ ]444 e l s e :445 empty = False446447 f o r i in range ( l en ( nx ) ) :448 x_nx . append (np . arange(−xrange , xrange , 2 ∗ xrange / (2 ∗∗ nx [ i

] ) ) ) # I n t e r v a l l x−Richtung449 i f empty :450 ny . append ( nx [ i ] )451 y_ny . append (np . arange(−yrange , yrange , 2 ∗ yrange / (2 ∗∗

nx [ i ] ) ) ) # I n t e r v a l l y−Richtung452 e l s e :453 i f i < l en ( ny ) :454 y_ny . append (np . arange(−yrange , yrange , 2 ∗ yrange / (2

∗∗ ny [ i ] ) ) ) # I n t e r v a l l y−Richtung455 e l s e :456 p r i n t ( "Warnung : Weniger y−G i t t e r g r ößen a l s x−G i t t e r g r ö

ßen angegeben ! Verwende ny = nx . " )457 ny . append ( nx [ i ] )458 y_ny . append (np . arange(−yrange , yrange , 2 ∗ yrange / (2

∗∗ ny [ i ] ) ) ) # I n t e r v a l l y−Richtung

81

Page 90: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Anhang A. Appendix

459 i f e f f e k t i v i s Fa l se :460 i f potV i s None :461 r a i s e IOError ( ’ Keine P o t e n t i a l f u n k t i o n angegeben ! ’ )462 e l s e :463 V = potV (x_nx [ i ] , y_ny [ i ] )464 al lPs i_nxy . append ( mult iple_propagation_times (x_nx [ i ] , y_ny [ i ] ,

t imes , V, e f f e k t i v=e f f e k t i v ,465 useVBH=useVBH ,

useM=useM ,useU=useU ) )

466 i f e f f e k t i v i s Fa l se and handl ing [ 0 ] == " save " :467 with open ( handl ing [ 1 ] , ’wb ’ ) as i n p u t _ f i l e :468 p r i n t ( " \t2D−Daten werden g e s p e i c h e r t in " + handl ing [ 1 ] + " \n"

)469 p i c k l e . dump ( [ nx , ny , a l lPs i_nxy ] , i n p u t _ f i l e )470471 normDiff = np . empty ( ( l en ( nx ) −1, l en ( t imes ) ) )472 f o r i in range ( l en ( nx ) −1) :473 f o r t in range ( l en ( t imes ) ) :474 d i f f_x = np . abs ( nx [ i +1] − nx [ i ] )475 d i f f_y = np . abs ( ny [ i + 1 ] − ny [ i ] )476 i f nx [ i ] <= nx [ i +1] and ny [ i ] <= ny [ i +1] :477 normFact = s q r t ( (4 ∗ xrange ∗ yrange ) / ( nx [ i ] ∗ ny [ i ] ) )478 normDiff [ i ] [ t ] = normFact ∗ np . l i n a l g . norm(479 al lPs i_nxy [ i ] [ t ] [ 1 ] − ( a l lPs i_nxy [ i + 1 ] [ t ] [ 1 ] ) [ : : 2 ∗ ∗

di f f_y , : : 2 ∗ ∗ d i f f_x ] )480 e l i f nx [ i ] <= nx [ i +1] and ny [ i ] > ny [ i +1] :481 normFact = s q r t ( (4 ∗ xrange ∗ yrange ) / ( nx [ i ] ∗ ny [ i +

1 ] ) )482 normDiff [ i ] [ t ] = normFact ∗ np . l i n a l g . norm(483 ( al lPs i_nxy [ i ] [ t ] [ 1 ] ) [ : : 2 ∗∗ di f f_y , : ] − ( a l lPs i_nxy [

i + 1 ] [ t ] [ 1 ] ) [ : , : : 2 ∗∗ d i f f_x ] )484 e l i f nx [ i ] > nx [ i + 1 ] and ny [ i ] <= ny [ i + 1 ] :485 normFact = s q r t ( (4 ∗ xrange ∗ yrange ) / ( nx [ i + 1 ] ∗ ny [ i

] ) )486 normDiff [ i ] [ t ] = normFact ∗ np . l i n a l g . norm(487 ( al lPs i_nxy [ i ] [ t ] [ 1 ] ) [ : , : : 2 ∗∗ d i f f_x ] − ( a l lPs i_nxy [

i + 1 ] [ t ] [ 1 ] ) [ : : 2 ∗∗ di f f_y , : ] )488 e l s e :489 normFact = s q r t ( (4 ∗ xrange ∗ yrange ) / ( nx [ i + 1 ] ∗ ny [ i

+ 1 ] ) )490 normDiff [ i ] [ t ] = normFact ∗ np . l i n a l g . norm(491 ( al lPs i_nxy [ i ] [ t ] [ 1 ] ) [ : : 2 ∗ ∗ di f f_y , : : 2 ∗ ∗ d i f f_x ] ) −

al lPs i_nxy [ i + 1 ] [ t ] [ 1 ]492493 return normDiff494495496 # Plot des G i t t e r g r ößen−V e r g l e i c h s497 de f p lot_gr id ( normDiff , normDiff_eff , times , ∗nx , ny=None ) :498 " " " P l o t t e t b i s zu v i e r Ergebn i s se des V e r g l e i c h s von den ( b i s zu f ünf )

G i t t e r g r ößen . " " "499 i f l en ( nx ) > 5 or l en ( nx ) <= 1 :500 p r i n t ( "Zu wenig/ v i e l e G i t t e r g r ößen angegeben ! (Nö t i g : 2 b i s 5) " )501 p r i n t ( " Abbruch . " )502 return

82

Page 91: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

A.2. Quelltext

503 e l s e :504 i f ny i s not None and l en ( nx ) != len ( ny ) :505 p r i n t ( "Zu wenig/ v i e l e G i t t e r g r ößen angegeben ! (Nö t i g : 2 b i s 5)

" )506 p r i n t ( " Abbruch . " )507 return508509 p l t . f i g u r e ( )510 i f ny i s None :511 l b l = r ’ $\mathrm{{ ges : }} $ ’ + r ’ $ m_x = {0} , n_x = {1}$ ’512 l b l _ e f f = r ’ $\mathrm{{ e f f : }} $ ’ + r ’ $ m_x = {0} , n_x = {1}$ ’513 p l t . y l a b e l ( r ’ $\Vert \ mathit {\ Psi }^\ vareps i lon_ {m_x} − \ mathit {\ Psi

}^\ vareps i lon_ {n_x} \Vert_{L^2}$ ’ )514 nxy = nx515 e l s e :516 l b l = r ’ $\mathrm{{ ges : }} $ ’ + r ’ $ m_y = {0} , n_y = {1}$ ’517 l b l _ e f f = r ’ $\mathrm{{ e f f : }} $ ’ + r ’ $ m_y = {0} , n_y = {1}$ ’518 p l t . y l a b e l ( r ’ $\Vert \ mathit {\ Psi }^\ vareps i lon_ {m_y} − \ mathit {\ Psi

}^\ vareps i lon_ {n_y} \Vert_{L^2}$ ’ )519 nxy = ny520521 p l t . p l o t ( times , normDif f_ef f [ 0 ] , c=’ orange ’ , marker=’ v ’ , mec=’ None ’ ,

l s= ’ ’ , l a b e l=l b l _ e f f . format ( nxy [ 0 ] , nxy [ 1 ] ) )522 p l t . p l o t ( times , normDiff [ 0 ] , c o l o r=’ orange ’ , marker=’+’ , l s=’ ’ , l a b e l=

l b l . format ( nxy [ 0 ] , nxy [ 1 ] ) )523 i f l en ( nxy ) >= 3 :524 p l t . p l o t ( times , normDif f_ef f [ 1 ] , c=’ f i r e b r i c k ’ , marker=’ v ’ , mec=’

None ’ , l s= ’ ’ ,525 l a b e l=l b l _ e f f . format ( nxy [ 1 ] , nxy [ 2 ] ) )526 p l t . p l o t ( times , normDiff [ 1 ] , c o l o r=’ f i r e b r i c k ’ , marker=’+’ , l s=’ ’ ,

l a b e l=l b l . format ( nxy [ 1 ] , nxy [ 2 ] ) )527 i f l en ( nxy ) >= 4 :528 p l t . p l o t ( times , normDif f_ef f [ 2 ] , ’ bv ’ , markeredgeco lor=’ None ’ , l s=

’ ’ , l a b e l=l b l _ e f f . format ( nxy [ 2 ] , nxy [ 3 ] ) )529 p l t . p l o t ( times , normDiff [ 2 ] , ’ b+’ , l a b e l=l b l . format ( nxy [ 2 ] , nxy

[ 3 ] ) )530 i f l en ( nxy ) >= 5 :531 p l t . p l o t ( times , normDif f_ef f [ 3 ] , ’ gv ’ , markeredgeco lor=’ None ’ , l s=

’ ’ , l a b e l=l b l _ e f f . format ( nxy [ 3 ] , nxy [ 4 ] ) )532 p l t . p l o t ( times , normDiff [ 3 ] , ’ g+’ , l a b e l=l b l . format ( nxy [ 3 ] , nxy

[ 4 ] ) )533534 p l t . y s c a l e ( ’ symlog ’ , l i n t h r e s h y =10 ∗∗ (−9) )535 p l t . l egend ( l o c=’ bes t ’ , numpoints=1, nco l =2, columnspacing =1,

handletextpad =0)536 p l t . x l a b e l ( ’ Ze i t $t$ ’ )537 p l t . show ( )538539540 # Verg l e i che Ze i tentwick lungen be i u n t e r s c h i e d l i c h e r S c h r i t t g r öße beim

Strang−S p l i t t i n g541 de f compare_spl i t t ing_steps (x , y , times , ∗hs , potV=None , e f f e k t i v=True ,542 useVBH=True , useM=True , useU=False , handl ing="

new" ) :543 " " " Berechnet d i e normierte D i f f e r e n z der Wel lenfunkt ionen zu gegebenen

S c h r i t t w e i t e n des Strang−S p l i t t i n g s " " "

83

Page 92: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Anhang A. Appendix

544 i f e f f e k t i v i s Fa l se and handl ing [ 0 ] == " load " :545 with open ( handl ing [ 1 ] , ’ rb ’ ) as i n p u t _ f i l e :546 p r i n t ( " \t2D−Daten werden ge laden . . . " )547 a l lPs i_n = p i c k l e . load ( i n p u t _ f i l e )548 e l s e :549 a l lPs i_n = [ ]550 f o r i in range ( l en ( hs ) ) :551 a l lPs i_n . append ( mult iple_propagat ion_times (552 x , y , times , potV , e f f e k t i v=e f f e k t i v , useVBH=useVBH , useM=

useM , useU=useU , t_delta=hs [ i ] ) )553 i f e f f e k t i v i s Fa l se and handl ing [ 0 ] == " save " :554 with open ( handl ing [ 1 ] , ’wb ’ ) as i n p u t _ f i l e :555 p r i n t ( " \t2D−Daten werden g e s p e i c h e r t in " + handl ing [ 1 ] + " \n"

)556 p i c k l e . dump( al lPsi_n , i n p u t _ f i l e )557558 normDiff = np . empty ( ( l en ( hs ) −1, l en ( t imes ) ) )559 normFact = s q r t ( (4 ∗ np . abs ( x [ 0 ] ) ∗ np . abs ( y [ 0 ] ) ) / ( l en ( x ) ∗ l en ( y ) ) )560 f o r i in range ( l en ( hs ) −1) :561 f o r t in range ( l en ( t imes ) ) :562 normDiff [ i ] [ t ] = normFact ∗ np . l i n a l g . norm( a l lPs i_n [ i + 1 ] [ t

] [ 1 ] − a l lPs i_n [ i ] [ t ] [ 1 ] )563564 return normDiff565566567 # Plot des Strang−S p l i t t i n g −V e r g l e i c h s568 de f p l o t _ s p l i t t i n g ( normDiff , normDiff_eff , times , ∗hs ) :569 " " " P l o t t e t b i s zu v i e r Ergebn i s se des V e r g l e i c h s von den ( b i s zu f ünf )

S c h r i t t w e i t e n im Strang−S p l i t t i n g . " " "570 i f l en ( hs ) > 5 or l en ( hs ) <= 1 :571 p r i n t ( "Zu wenig/ v i e l e S c h r i t t g r ößen angegeben ! (Nö t i g : 2 b i s 5) " )572 p r i n t ( " Abbruch . " )573 return574575 hs_frac = [ Fract ion ( hs [ i ] ) . l imit_denominator ( ) f o r i in range ( l en ( hs ) )

]576 hs_lbl = [ r ’ \ f r a c {{{0}}}{{{1}}} ’ . format ( hs_frac [ i ] . numerator , hs_frac [

i ] . denominator ) f o r i in range ( l en ( hs_frac ) ) ]577 l b l = r ’ $\mathrm{{ ges : }} $ ’ + r ’ $ h_1 = {0} , h_2 = {1}$ ’578 l b l _ e f f = r ’ $\mathrm{{ e f f : }} $ ’ + r ’ $ h_1 = {0} , h_2 = {1}$ ’579580 p l t . f i g u r e ( )581 p l t . p l o t ( times , normDif f_ef f [ 0 ] , c=’ orange ’ , marker=’ v ’ , mec=’ None ’ ,

l s= ’ ’ ,582 l a b e l=l b l _ e f f . format ( hs_lbl [ 0 ] , hs_lbl [ 1 ] ) )583 p l t . p l o t ( times , normDiff [ 0 ] , c o l o r=’ orange ’ , marker=’+’ , l s=’ ’ , l a b e l=

l b l . format ( hs_lbl [ 0 ] , hs_lbl [ 1 ] ) )584 i f l en ( hs ) >= 3 :585 p l t . p l o t ( times , normDif f_ef f [ 1 ] , c=’ f i r e b r i c k ’ , marker=’ v ’ , mec=’

None ’ , l s= ’ ’ ,586 l a b e l=l b l _ e f f . format ( hs_lbl [ 1 ] , hs_lbl [ 2 ] ) )587 p l t . p l o t ( times , normDiff [ 1 ] , c o l o r=’ f i r e b r i c k ’ , marker=’+’ , l s=’ ’ ,

l a b e l=l b l . format ( hs_lbl [ 1 ] , hs_lbl [ 2 ] ) )588 i f l en ( hs ) >= 4 :589 p l t . p l o t ( times , normDif f_ef f [ 2 ] , ’ bv ’ , markeredgeco lor=’ None ’ , l s=

84

Page 93: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

A.2. Quelltext

’ ’ ,590 l a b e l=l b l _ e f f . format ( hs_lbl [ 2 ] , hs_lbl [ 3 ] ) )591 p l t . p l o t ( times , normDiff [ 2 ] , ’ b+’ , l a b e l=l b l . format ( hs_lbl [ 2 ] ,

hs_lbl [ 3 ] ) )592 i f l en ( hs ) == 5 :593 p l t . p l o t ( times , normDif f_ef f [ 3 ] , ’ gv ’ , markeredgeco lor=’ None ’ , l s=

’ ’ ,594 l a b e l=l b l _ e f f . format ( hs_lbl [ 3 ] , hs_lbl [ 4 ] ) )595 p l t . p l o t ( times , normDiff [ 3 ] , ’ g+’ , l a b e l=l b l . format ( hs_lbl [ 3 ] ,

hs_lbl [ 4 ] ) )596597 p l t . y s c a l e ( ’ symlog ’ , l i n t h r e s h y =10 ∗∗ (−9) )598 p l t . l egend ( l o c=’ bes t ’ , numpoints=1, nco l =2, columnspacing =1,

handletextpad =0)599 p l t . x l a b e l ( ’ Ze i t $t$ ’ )600 p l t . y l a b e l ( r ’ $\Vert \ mathit {\ Psi }^\ vareps i lon_ {h_1} − \ mathit {\ Psi }^\

vareps i lon_ {h_2} \Vert_{L^2}$ ’ )601 p l t . show ( )602603604 # Vektor aus Ze i t entwick lung be i ver sch i edenen e p s i l o n zu einem festem

Zeitpunkt605 de f mult iple_propagation_eps (x , y , t , e p s i l o n s , ∗potV , e f f e k t i v=True ,

useVBH=True , useM=True , useU=None) :606 " " " Berechnet d i e Ze i t entwick lung bezü g l i c h v e r s c h i e d e n e r e p s i l o n zu

einem f e s t e n Zeitpunkten . " " "607 i f e f f e k t i v i s Fa l se :608 i f l en ( potV ) > 1 :609 p r i n t ( ’Warnung : Zu v i e l e Argumente ! Üb e r f l ü s s i g e werden

i g n o r i e r t . ’ )610 e l i f l en ( potV ) == 0 or potV [ 0 ] i s None :611 r a i s e IOError ( ’ Kein P o t e n t i a l angegeben ! ’ )612 V = potV [ 0 ]613 i f useU i s None :614 useU = False615 typ = " 2D"616 e l s e :617 i f l en ( potV ) > 0 and potV [ 0 ] i s not None :618 p r i n t ( ’Warnung : Zu v i e l e Argumente ! Üb e r f l ü s s i g e s P o t e n t i a l

wird i g n o r i e r t . ’ )619 i f useU i s None :620 i f useM i s True :621 useU = True622 e l s e :623 useU = False624 typ = " 1D"625626 a l l P s i = [ ]627 phi_0 = calc_phi_0 (x , y )628 p r i n t ( ’ \ t [ ’ + typ + ’ ] Todo : ’ + ’ [ ’ + 100 ∗ ’− ’ + ’ ] ’ )629 sys . s tdout . wr i t e ( ’ \ t F o r t s c h r i t t : [ ’ )630 sys . s tdout . f l u s h ( )631 printedMark = 0632 s t a r t = time . time ( )633 f o r i in range ( l en ( e p s i l o n s ) ) :634 psi_0 = calc_psi_0 (x , e p s i l o n s [ i ] )

85

Page 94: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Anhang A. Appendix

635 i f e f f e k t i v i s Fa l se :636 Psi_0 = calc_Psi (x , y , psi_0 , phi_0 , useU , e p s i l o n=e p s i l o n s [ i

] )637 Psi = time_propagation_2d (x , y , Psi_0 , V, t , e p s i l o n=e p s i l o n s [

i ] )638 e l s e :639 p s i = time_propagation_BO (x , psi_0 , t , useVBH , useM , e p s i l o n=

e p s i l o n s [ i ] )640 Psi = calc_Psi (x , y , ps i , phi_0 , useU , e p s i l o n=e p s i l o n s [ i ] )641 a l l P s i . append ( [ e p s i l o n s [ i ] , Ps i ] )642 i f l en ( e p s i l o n s )−1 == 0 :643 currentMark = 100644 e l s e :645 currentMark = 100 ∗ i // ( l en ( e p s i l o n s ) −1)646 i f printedMark < currentMark :647 sys . s tdout . wr i t e ( ’ ∗ ’ ∗ ( currentMark − printedMark ) )648 sys . s tdout . f l u s h ( )649 printedMark = currentMark650 sys . s tdout . wr i t e ( ’ ] \ n ’ )651 end = time . time ( )652 t o t a l = [ i n t ( ( end − s t a r t ) % 60) , ( end − s t a r t ) // 60 ]

# [ s , min , (h) ]653 i f t o t a l [ 1 ] >= 60 :654 t o t a l . append ( t o t a l [ 1 ] // 60)655 t o t a l [ 1 ] = i n t ( t o t a l [ 1 ] % 60)656 p r i n t ( ’ \tBenö t i g t e Ze i t : { : 0 2 . 0 f } : { : 0 2 d } : { : 0 2 d} Stunden\n ’ . format (

t o t a l [ 2 ] , t o t a l [ 1 ] , t o t a l [ 0 ] ) )657 e l s e :658 p r i n t ( ’ \tBenö t i g t e Ze i t : { : 0 2 . 0 f } : { : 0 2 d} Minuten\n ’ . format ( t o t a l

[ 1 ] , t o t a l [ 0 ] ) )659660 return a l l P s i661662663 # Verg l e i che Ze i tentwick lungen be i u n t e r s c h i e d l i c h e n Eps i l ons664 de f compare_error_eps (x , y , t , e p s i l o n s , V, handl ing="new" ) :665 " " " Berechnet d i e Fehler der Born−Oppenheimer−Approximation mit und

ohne Korrekturterm zu gegebenen Eps i l ons . " " "666 i f handl ing [ 0 ] == " load " :667 with open ( handl ing [ 1 ] , ’ rb ’ ) as i n p u t _ f i l e :668 p r i n t ( " \t2D−Daten werden ge laden . . . " )669 a l l P s i , a l lPs iU = p i c k l e . load ( i n p u t _ f i l e )670 e l s e :671 a l l P s i = mult iple_propagation_eps (x , y , t , e p s i l o n s , V, e f f e k t i v=

False )672 a l lPs iU = mult iple_propagation_eps (x , y , t , e p s i l o n s , V, e f f e k t i v=

False , useU=True )673 i f handl ing [ 0 ] == " save " :674 with open ( handl ing [ 1 ] , ’wb ’ ) as i n p u t _ f i l e :675 p r i n t ( " \t2D−Daten werden g e s p e i c h e r t in " + handl ing [ 1 ] + " \n"

)676 p i c k l e . dump ( [ a l l P s i , a l lPs iU ] , i n p u t _ f i l e )677678 a l l P s i _ e f f = mult iple_propagation_eps (x , y , t , e p s i l o n s , useVBH=False ,

useM=False , useU=False )679 a l lP s i_e f fV = mult iple_propagation_eps (x , y , t , e p s i l o n s , useM=False ,

86

Page 95: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

A.2. Quelltext

useU=False )680 al lPsi_effVM = mult iple_propagation_eps (x , y , t , e p s i l o n s , useU=False )681 allPsi_effVMU = mult iple_propagation_eps (x , y , t , e p s i l o n s )682683 normDif f_ef f = np . empty ( l en ( e p s i l o n s ) )684 normDiff_effV = np . empty ( l en ( e p s i l o n s ) )685 normDiff_effVM = np . empty ( l en ( e p s i l o n s ) )686 normDiff_effVMU = np . empty ( l en ( e p s i l o n s ) )687 normFact = s q r t ( (4 ∗ np . abs ( x [ 0 ] ) ∗ np . abs ( y [ 0 ] ) ) / ( l en ( x ) ∗ l en ( y ) ) )688 f o r i in range ( l en ( e p s i l o n s ) ) :689 normDif f_ef f [ i ] = normFact ∗ np . l i n a l g . norm( a l l P s i [ i ] [ 1 ] −

a l l P s i _ e f f [ i ] [ 1 ] )690 normDiff_effV [ i ] = normFact ∗ np . l i n a l g . norm( a l l P s i [ i ] [ 1 ] −

a l lP s i_e f fV [ i ] [ 1 ] )691 normDiff_effVM [ i ] = normFact ∗ np . l i n a l g . norm( a l l P s i [ i ] [ 1 ] −

al lPsi_effVM [ i ] [ 1 ] )692 normDiff_effVMU [ i ] = normFact ∗ np . l i n a l g . norm( a l lPs iU [ i ] [ 1 ] −

allPsi_effVMU [ i ] [ 1 ] )693694 return normDiff_eff , normDiff_effV , normDiff_effVM , normDiff_effVMU695696697 # Plot des Epsi lon−V e r g l e i c h s698 de f plot_error_eps ( normDiff_eff , normDiff_effV , normDiff_effVM ,

normDiff_effVMU , t , e p s i l o n s , show_bound=True ) :699 " " " P l o t t e t d i e Fehler der Born−Oppenheimer−Approximation mit und ohne

Korrekturterm zu gegebenen Eps i l ons . " " "700 p l t . f i g u r e ( )701 i f show_bound i s True :702 c1 , c2 = [ ] , [ ]703 f o r i in range ( l en ( e p s i l o n s ) ) :704 c1 . append ( normDiff_effV [ i ] / ( e p s i l o n s [ i ] ∗ (1 + abs ( t ) ) ) )705 c2 . append ( normDiff_effVMU [ i ] / ( e p s i l o n s [ i ]∗∗2 ∗ (1+abs ( t ) ) ) )706 c1 = np . amax( c1 )707 c2 = np . amax( c2 )708 p l t . p l o t (np . log2 ( e p s i l o n s ) , np . log2 ( e p s i l o n s ) + np . log2 ( c1 ∗ (1 +

abs ( t ) ) ) , c=’ gray ’ , l s= ’ dotted ’ )709 p l t . p l o t (np . log2 ( e p s i l o n s ) , 2 ∗ np . log2 ( e p s i l o n s ) + np . log2 ( c2 ∗

(1 + abs ( t ) ) ) , c=’ gray ’ , l s=’ dotted ’ )710711 p l t . p l o t (np . log2 ( e p s i l o n s ) , np . log2 ( normDif f_ef f ) , c=’ f i r e b r i c k ’ ,

marker=’ v ’ , mec=’ None ’ ,712 l s=’−− ’ , l a b e l=r ’ ohne $V_{\mathrm{BH}}$ und $\ mathcal {M}$ ’ )713 p l t . p l o t (np . log2 ( e p s i l o n s ) , np . log2 ( normDiff_effV ) , ’ bv−− ’ , mec=’ None ’

, l a b e l=r ’ ohne $\ mathcal {M}$ ’ )714 p l t . p l o t (np . log2 ( e p s i l o n s ) , np . log2 ( normDiff_effVM ) , c=’ orange ’ ,

marker=’ v ’ , mec=’ None ’ ,715 l s=’−− ’ , l a b e l=r ’ mit $\ mathcal {M}$ , ohne $U_{(1) }^{\

v a r e p s i l o n \ , ∗}$ ’ )716 p l t . p l o t (np . log2 ( e p s i l o n s ) , np . log2 ( normDiff_effVMU ) , ’ gv−− ’ , mec=’

None ’ ,717 l a b e l=r ’ mit $\ mathcal {M}$ und $U_{(1) }^{\ v a r e p s i l o n \ , ∗}$ ’ )718 p l t . l egend ( l o c=’ bes t ’ )719 p l t . t i t l e ( ’ $t = {}$ ’ . format ( t ) )720 p l t . x l a b e l ( r ’ $\ log_2 (\ v a r e p s i l o n ) $ ’ )721 p l t . y l a b e l ( r ’ $\ log_2 (\ Vert \ mathit {\ Psi }^\ v a r e p s i l o n − \ mathit {\ Psi }^\

87

Page 96: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Anhang A. Appendix

vareps i lon_ {\mathrm{ e f f }} \Vert_{L^2}) $ ’ )722 p l t . show ( )723724725 # Verg l e i che d i e Entwicklung des Feh l e r s über d i e Ze i t726 de f compare_error_time (x , y , times , V, handl ing="new" ) :727 " " " Berechnet den Fehler der Born−Oppenheimer−Approximation mit und

ohne Korrekturterm zu gegebenen Ze i t en . " " "728 i f handl ing [ 0 ] == " load " :729 with open ( handl ing [ 1 ] , ’ rb ’ ) as i n p u t _ f i l e :730 p r i n t ( " \t2D−Daten werden ge laden . . . " )731 a l l P s i , a l lPs iU = p i c k l e . load ( i n p u t _ f i l e )732 e l s e :733 a l l P s i = mult iple_propagation_times (x , y , times , V, e f f e k t i v=False

)734 a l lPs iU = mult iple_propagation_times (x , y , times , V, e f f e k t i v=

False , useU=True )735 i f handl ing [ 0 ] == " save " :736 with open ( handl ing [ 1 ] , ’wb ’ ) as i n p u t _ f i l e :737 p r i n t ( " \t2D−Daten werden g e s p e i c h e r t in " + handl ing [ 1 ] + " \n"

)738 p i c k l e . dump ( [ a l l P s i , a l lPs iU ] , i n p u t _ f i l e )739740 a l l P s i _ e f f = mult iple_propagation_times (x , y , times , useVBH=False ,

useM=False , useU=False )741 a l lP s i_e f fV = mult iple_propagation_times (x , y , times , useM=False , useU

=False )742 al lPsi_effVM = mult iple_propagat ion_times (x , y , times , useU=False )743 allPsi_effVMU = mult iple_propagat ion_times (x , y , t imes )744745 normDif f_ef f = np . empty ( l en ( t imes ) )746 normDiff_effV = np . empty ( l en ( t imes ) )747 normDiff_effVM = np . empty ( l en ( t imes ) )748 normDiff_effVMU = np . empty ( l en ( t imes ) )749 normFact = s q r t ( (4 ∗ np . abs ( x [ 0 ] ) ∗ np . abs ( y [ 0 ] ) ) / ( l en ( x ) ∗ l en ( y ) ) )750 f o r i in range ( l en ( t imes ) ) :751 normDif f_ef f [ i ] = normFact ∗ np . l i n a l g . norm( a l l P s i [ i ] [ 1 ] −

a l l P s i _ e f f [ i ] [ 1 ] )752 normDiff_effV [ i ] = normFact ∗ np . l i n a l g . norm( a l l P s i [ i ] [ 1 ] −

a l lP s i_e f fV [ i ] [ 1 ] )753 normDiff_effVM [ i ] = normFact ∗ np . l i n a l g . norm( a l l P s i [ i ] [ 1 ] −

al lPsi_effVM [ i ] [ 1 ] )754 normDiff_effVMU [ i ] = normFact ∗ np . l i n a l g . norm( a l lPs iU [ i ] [ 1 ] −

allPsi_effVMU [ i ] [ 1 ] )755756 return normDiff_eff , normDiff_effV , normDiff_effVM , normDiff_effVMU757758759 # Plot des Feh l e r s über d i e Ze i t760 de f plot_error_time ( times , normDiff_eff , normDiff_effV , normDiff_effVM ,

normDiff_effVMU ) :761 " " " P l o t t e t d i e Fehler der Born−Oppenheimer−Approximation mit und ohne

Korrekturterm zu gegebenen Ze i t en . " " "762 p l t . f i g u r e ( )763 p l t . p l o t ( times , normDiff_eff , c=’ f i r e b r i c k ’ , marker=’ v ’ , mec=’ None ’ ,764 l s=’−− ’ , l a b e l=r ’ ohne $V_{\mathrm{BH}}$ und $\ mathcal {M}$ ’ )

88

Page 97: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

A.2. Quelltext

765 p l t . p l o t ( times , normDiff_effV , ’ bv−− ’ , mec=’ None ’ , l a b e l=r ’ ohne $\mathcal {M}$ ’ )

766 p l t . p l o t ( times , normDiff_effVM , c=’ orange ’ , marker=’ v ’ , mec=’ None ’ , l s=’−− ’ ,

767 l a b e l=r ’ mit $\ mathcal {M}$ , ohne $U_{(1) }^{\ v a r e p s i l o n \ , ∗}$ ’ )768 p l t . p l o t ( times , normDiff_effVMU , ’ gv−− ’ , mec=’ None ’ , l a b e l=r ’ mit $\

mathcal {M}$ und $U_{(1) }^{\ v a r e p s i l o n \ , ∗}$ ’ )769 p l t . l egend ( l o c=’ bes t ’ )770 p l t . x l a b e l ( ’ Ze i t $t$ ’ )771 p l t . y l a b e l ( r ’ $\Vert \ mathit {\ Psi }^\ v a r e p s i l o n − \ mathit {\ Psi }^\

vareps i lon_ {\mathrm{ e f f }} \Vert_{L^2}$ ’ )772 p l t . show ( )773774775 # Verg l e i che d i e Entwicklung des Feh l e r s der Approximation zwe i t e r Ordnung

über d i e Ze i t776 de f compare_error_time_VMU (x , y , times , e p s i l o n s , V, handl ing="new" ) :777 " " " Berechnet den Fehler der Born−Oppenheimer−Approximation zwe i t e r

Ordnung zu gegebenen Ze i ten . " " "778 i f handl ing [ 0 ] == " load " :779 with open ( handl ing [ 1 ] , ’ rb ’ ) as i n p u t _ f i l e :780 p r i n t ( " \t2D−Daten werden ge laden . . . " )781 a l lPs iU = p i c k l e . load ( i n p u t _ f i l e )782 e l s e :783 a l lPs iU = [ ]784 f o r i in range ( l en ( e p s i l o n s ) ) :785 PsiU = mult iple_propagation_times (x , y , times , V, e f f e k t i v=

False , useU=True , e p s i l o n=e p s i l o n s [ i ] )786 a l lPs iU . append ( [ e p s i l o n s [ i ] , PsiU ] )787 i f handl ing [ 0 ] == " save " :788 with open ( handl ing [ 1 ] , ’wb ’ ) as i n p u t _ f i l e :789 p r i n t ( " \t2D−Daten werden g e s p e i c h e r t in " + handl ing [ 1 ] + " \n"

)790 p i c k l e . dump( al lPs iU , i n p u t _ f i l e )791792 allPsi_effVMU = [ ]793 f o r i in range ( l en ( e p s i l o n s ) ) :794 Psi_effVMU = mult iple_propagation_times (x , y , times , e p s i l o n=

e p s i l o n s [ i ] )795 allPsi_effVMU . append ( [ e p s i l o n s [ i ] , Psi_effVMU ] )796797 normDiff_effVMU_eps = [ ]798 normDiff_effVMU = np . empty ( l en ( t imes ) )799 normFact = s q r t ( (4 ∗ np . abs ( x [ 0 ] ) ∗ np . abs ( y [ 0 ] ) ) / ( l en ( x ) ∗ l en ( y ) ) )800801 f o r j in range ( l en ( e p s i l o n s ) ) :802 f o r i in range ( l en ( t imes ) ) :803 normDiff_effVMU [ i ] = normFact ∗ np . l i n a l g . norm( a l lPs iU [ j ] [ 1 ] [ i

] [ 1 ] − allPsi_effVMU [ j ] [ 1 ] [ i ] [ 1 ] )804 normDiff_effVMU_eps . append ( [ e p s i l o n s [ j ] , np . copy ( normDiff_effVMU )

] )805806 return normDiff_effVMU_eps807808809 # Plot des Feh l e r s der zweiten Ordnung über d i e Ze i t

89

Page 98: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Anhang A. Appendix

810 de f plot_error_time_VMU ( times , normDiff_effVMU_eps , i n t e r p o l a t e=False ) :811 " " " P l o t t e t b i s zu d r e i Fehlerkurven der Born−Oppenheimer−Approximation

zwe i t e r Ordnung . " " "812 i f l en ( normDiff_effVMU_eps ) > 3 or l en ( normDiff_effVMU_eps ) < 1 :813 p r i n t ( "Zu wenig/ v i e l e Argumente angegeben ! ( Zul ä s s i g : 1 b i s 4) " )814 p r i n t ( " Abbruch . " )815 return816817 e p s i l o n s = [ normDiff_effVMU_eps [ i ] [ 0 ] f o r i in range ( l en (

normDiff_effVMU_eps ) ) ]818 eps_frac = [ Fract ion ( e p s i l o n s [ i ] ) . l imit_denominator ( ) f o r i in range (

l en ( e p s i l o n s ) ) ]819820 p l t . f i g u r e ( )821 p l t . p l o t ( times , normDiff_effVMU_eps [ 0 ] [ 1 ] , ’ gs−− ’ , mec=’ None ’ ,822 l a b e l=r ’ $\ v a r e p s i l o n=\f r a c {{{0}}}{{{1}}} $ ’ . format ( eps_frac

[ 0 ] . numerator , eps_frac [ 0 ] . denominator ) )823 i f i n t e r p o l a t e i s True :824 c = np . poly1d (np . p o l y f i t ( t imes [ 1 : ] , normDiff_effVMU_eps [ 0 ] [ 1 ] [ 1 : ] ,

1) )825 p r i n t ( " I n t e r p o l a t i o n s g e r a d e f ür eps = " + s t r ( eps_frac [ 0 ] .

numerator ) + " / " + s t r (826 eps_frac [ 0 ] . denominator ) + " : " )827 p r i n t ( " \ t ( " + s t r ( eps_frac [ 0 ] . numerator ) + " / " + s t r ( eps_frac [ 0 ] .

denominator ) + " ) ^2 ∗ ( " + s t r (828 normDiff_effVMU_eps [ 0 ] [ 0 ] ∗∗ (−2) ∗ c [ 1 ] ) + " ∗ t + " + s t r (

normDiff_effVMU_eps [ 0 ] [ 0 ] ∗∗ (−2) ∗ c [ 0 ] ) + " ) " )829 p l t . p l o t ( t imes [ 1 : ] , c ( t imes [ 1 : ] ) , ’ darkorange ’ )830 i f l en ( normDiff_effVMU_eps ) >= 2 :831 p l t . p l o t ( times , normDiff_effVMU_eps [ 1 ] [ 1 ] , ’ gd−− ’ , mec=’ None ’ ,832 l a b e l=r ’ $\ v a r e p s i l o n=\f r a c {{{0}}}{{{1}}} $ ’ . format (

eps_frac [ 1 ] . numerator , eps_frac [ 1 ] . denominator ) )833 i f i n t e r p o l a t e i s True :834 c = np . poly1d (np . p o l y f i t ( t imes [ 1 : ] , normDiff_effVMU_eps

[ 1 ] [ 1 ] [ 1 : ] , 1) )835 p r i n t ( " I n t e r p o l a t i o n s g e r a d e f ür eps = " + s t r ( eps_frac [ 1 ] .

numerator ) + " / " + s t r (836 eps_frac [ 1 ] . denominator ) + " : " )837 p r i n t ( " \ t ( " + s t r ( eps_frac [ 1 ] . numerator ) + " / " + s t r ( eps_frac

[ 1 ] . denominator ) + " ) ^2 ∗ ( " + s t r (838 normDiff_effVMU_eps [ 1 ] [ 0 ] ∗∗ (−2) ∗ c [ 1 ] ) + " ∗ t + " +

s t r (839 normDiff_effVMU_eps [ 1 ] [ 0 ] ∗∗ (−2) ∗ c [ 0 ] ) + " ) " )840 p l t . p l o t ( t imes [ 1 : ] , c ( t imes [ 1 : ] ) , ’ darkorange ’ )841 i f l en ( normDiff_effVMU_eps ) == 3 :842 p l t . p l o t ( times , normDiff_effVMU_eps [ 2 ] [ 1 ] , ’ gv−− ’ , mec=’ None ’ ,843 l a b e l=r ’ $\ v a r e p s i l o n=\f r a c {{{0}}}{{{1}}} $ ’ . format (

eps_frac [ 2 ] . numerator , eps_frac [ 2 ] . denominator ) )844 i f i n t e r p o l a t e i s True :845 c = np . poly1d (np . p o l y f i t ( t imes [ 1 : ] , normDiff_effVMU_eps

[ 2 ] [ 1 ] [ 1 : ] , 1) )846 p r i n t ( " I n t e r p o l a t i o n s g e r a d e f ür eps = " + s t r ( eps_frac [ 2 ] .

numerator ) + " / " + s t r (847 eps_frac [ 2 ] . denominator ) + " : " )848 p r i n t ( " \ t ( " + s t r ( eps_frac [ 2 ] . numerator ) + " / " + s t r ( eps_frac

[ 2 ] . denominator ) + " ) ^2 ∗ ( " + s t r (

90

Page 99: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

A.2. Quelltext

849 normDiff_effVMU_eps [ 2 ] [ 0 ] ∗∗ (−2) ∗ c [ 1 ] ) + " ∗ t + " +s t r (

850 normDiff_effVMU_eps [ 2 ] [ 0 ] ∗∗ (−2) ∗ c [ 0 ] ) + " ) " )851 p l t . p l o t ( t imes [ 1 : ] , c ( t imes [ 1 : ] ) , ’ darkorange ’ )852853 p l t . l egend ( l o c=’ bes t ’ )854 p l t . x l a b e l ( ’ Ze i t $t$ ’ )855 p l t . y l a b e l ( r ’ $\Vert \ mathit {\ Psi }^\ v a r e p s i l o n − \ mathit {\ Psi }^\

vareps i lon_ {\mathrm{BO}} \Vert_{L^2}$ ’ )856 p l t . show ( )857858859 # Abfrage zur Datenerhebung860 de f ask_data_handling ( ) :861 " " " Frä gt ab , ob d i e nö t i g e n Daten e i n g e l e s e n oder neu berechnet , sowie

a n s c h l i e ßend g e s p e i c h e r t werden s o l l e n . " " "862 r e t r i e s = 3863 whi l e True :864 ask = s t r ( input ( " " "865 S o l l e n d i e Daten neu berechnet oder e i n g e l e s e n werden?866 [ 1 ] 2D−Daten neu berechnen867 [ 2 ] 2D−Daten e i n l e s e n868 −−−−−−−−−−−−−−−−−−−−−−−869 [ Z ] Zurück870 >> " " " ) )871872 i f ask == " 1 " :873 whi le True :874 ask = s t r ( input ( " " "875 S o l l e n d i e 2D−Daten a n s c h l i e ßend g e s p e i c h e r t werden [ J/n ] ? ( [ Z ]

Zurück )876 >> " " " ) )877 i f ask in ( "Z" , " z " ) :878 break879 e l i f ask in ( "n " , "N" ) :880 re turn "new"881 e l s e :882 whi l e True :883 ask = s t r ( input ( " " "884 B i t t e Spe icherp fad angeben : [ . / wave2d . pkl ] ( [ Z ] Zurück )885 >> " " " ) )886 i f ask in ( "Z" , " z " ) :887 break888 e l i f ask == ’ ’ :889 path = " . / wave2d . pkl "890 e l s e :891 path = ask892 i f os . path . e x i s t s ( path ) i s Fa l se :893 re turn [ " save " , path ]894 e l s e :895 whi l e True :896 ask = s t r ( input ( " " "897 Datei e x i s t i e r t b e r e i t s ! Übe r s ch r e iben [ j /N] ? ( [ Z ] Zurück )898 >> " " " ) )899 i f ask in ( " j " , " J " ) :900 re turn [ " save " , path ]

91

Page 100: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Anhang A. Appendix

901 e l s e :902 break903 e l i f ask == " 2 " :904 whi le True :905 ask = s t r ( input ( " " "906 B i t t e Q ue l l da t e i angeben : [ . / wave2d . pkl ] ( [ Z ] Zurück , [A]

Ordner inha l t anze igen )907 S icherhe i t swarnung : Nur vertrauenswü r d i g e Dateien verwenden !908 >> " " " ) )909 i f ask in ( "Z" , " z " ) :910 break911 e l i f ask in ( "A" , " a " ) :912 p r i n t (8 ∗ " " + s t r ( os . l i s t d i r ( " . / " ) ) )913 cont inue914 e l i f ask == ’ ’ :915 path = " . / wave2d . pkl "916 e l s e :917 path = ask918 i f os . path . e x i s t s ( path ) i s Fa l se :919 p r i n t (8 ∗ " " + " Datei e x i s t i e r t n i cht ! " )920 cont inue921 return [ " load " , path ]922 e l i f ask in ( "Z" , " z " ) :923 re turn924 e l s e :925 p r i n t (6 ∗ " " + "Ungü l t i g e Eingabe ! " )926927 r e t r i e s −= 1928 i f r e t r i e s == 0 :929 re turn

92

Page 101: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Symbole und Zeichen

Symbol Bedeutung

~ reduziertes Plancksche WirkungsquantumH Hilbertraum〈·, ·〉H Skalarprodukt auf dem Hilbertraum H‖·‖H Norm auf dem Hilbertraum Hdist(·, ·) Metrik auf dem Hilbertraum HL(H) Menge der linearen beschränkten Operatoren auf HL(H1,H2) Menge der linearen beschränkten Operatoren auf H1 mit Werten

in H2

L2(Λ) Lebesgue-Raum der quadratintegriebaren Funktionen auf Λψ Komplexe Konjugation von ψ[A,B] Kommutator von A und BA∗ Adjungierter Operator von A ∈ L(H1,H2) in L(H2,H1)IdH Identität auf HP ProjektionP⊥ = (Id− P ) Zu P orthogonale Projektionσ(A) Spektrum eines Operators AR ResolventeO(ε) Landau-NotationCkb(Λ,H) Raum der beschränkten k-mal stetig differenzierbaren Funktionen

auf Λ mit Werten in H⊗ Tensorprodukt∇ Gradient∆ Laplace-OperatorχΛ Charakteristische Funktion auf Λδij Kronecker-Delta

93

Page 102: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls
Page 103: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Abbildungsverzeichnis

2.1 Durch eine Lücke lokal getrenntes Spektrum . . . . . . . . . . . . . . . . 112.2 Das Spektrum σ∗(x) ist lokal durch eine Lücke isoliert. . . . . . . . . . . 19

3.1 Einzelne isolierte Energieoberfläche. . . . . . . . . . . . . . . . . . . . . 27

4.1 Kontur- und HSL-Plot des Potentials V . . . . . . . . . . . . . . . . . . . 434.2 Dreidimensionaler Plot des Potentials V . . . . . . . . . . . . . . . . . . . 444.3 Konturplot der Wellenfunktion mit dem Potential (ε = 1/8). . . . . . . . 454.4 HSL-Plot der Wellenfunktion. . . . . . . . . . . . . . . . . . . . . . . . . 464.5 Born-Huang-Potential und Korrekturterm. . . . . . . . . . . . . . . . . . 504.6 Konturplot der Wellenfunktion mit dem Potential (ε = 1/16). . . . . . . 524.7 Vergleich verschiedener Gittergrößen (longitudinale Richtung). . . . . . 534.8 Vergleich verschiedener Gittergrößen (transversale Richtung). . . . . . . 544.9 Vergleich verschiedener Zeitschrittweiten. . . . . . . . . . . . . . . . . . 544.10 Vergleich der Approximationen bezüglich ε (t = 3). . . . . . . . . . . . . 554.11 Vergleich der Approximationen bezüglich ε (t = 12). . . . . . . . . . . . 564.12 Vergleich der Approximationen bezüglich t (ε = 1/8). . . . . . . . . . . . 574.13 Vergleich der Approximationen bezüglich t bei unterschiedlichen An-

fangsimpulsen (ε = 1/8). . . . . . . . . . . . . . . . . . . . . . . . . . . . 584.14 Vergleich der Approximationen bezüglich t (T = 50, ε = 1/4). . . . . . . 604.15 Fehlerkurven der Approximation zweiter Ordnung. . . . . . . . . . . . . 614.16 Fehlerkurven der Approximation zweiter Ordnung (interpoliert). . . . . 62

95

Page 104: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls
Page 105: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Tabellenverzeichnis

4.1 Laufzeiten. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 594.2 Koeffizienten der Interpolation. . . . . . . . . . . . . . . . . . . . . . . . 61

97

Page 106: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls
Page 107: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Literaturverzeichnis

[AE99] Joseph E. Avron und Alexander Elgart. „Adiabatic theorem without a gapcondition“. In: Comm. Math. Phys. 203.2 (1999), S. 445–463. issn: 0010-3616.doi: 10.1007/s002200050620.

[BF28] Max Born und Vladimir Fock. „Beweis des Adiabatensatzes“. In: Zeitschriftfür Physik 51.3 (1928), S. 165–180. issn: 0044-3328. doi: 10.1007/BF01343193.

[BO27] Max Born und Robert Oppenheimer. „Zur Quantentheorie der Molekeln“.In: Annalen der Physik 389.20 (1927), S. 457–484. issn: 1521-3889. doi: 10.1002/andp.19273892002.

[Bor98] Folkmar Bornemann. Homogenization in time of singularly perturbed me-chanical systems. Bd. 1687. Lecture Notes in Mathematics. Springer-Verlag,Berlin, 1998, S. xii+156. isbn: 3-540-64447-4. doi: 10.1007/BFb0092091.

[Böl09] Sebastian Bölzle. „Explicit Error Bounds for the First Order Adiabatic Ap-proximation with an Application to Quantum Wave Guides“. Diplomarbeit.Eberhard Karls Universität Tübingen, 2009.

[Ehr16] Paul Ehrenfest. „Adiabatische Invarianten und Quantentheorie“. In: Annalender Physik 356.19 (1916), S. 327–352. issn: 1521-3889. doi: 10.1002/andp.19163561905.

[HS96] Peter D. Hislop und Israel M. Sigal. Introduction to spectral theory. Bd. 113.Applied Mathematical Sciences. With applications to Schrödinger operators.Springer-Verlag, New York, 1996, S. x+337. isbn: 0-387-94501-6. doi: 10.1007/978-1-4612-0741-2.

[HST01] Frank Hövermann, Herbert Spohn und Stefan Teufel. „Semiclassical limit forthe Schrödinger equation with a short scale periodic potential“. In: Comm.Math. Phys. 215.3 (2001), S. 609–629. issn: 0010-3616. doi: 10.1007/s002200000314.

[Hun07] John D. Hunter. „Matplotlib: A 2D graphics environment“. In: ComputingIn Science & Engineering 9.3 (2007), S. 90–95. doi: 10.1109/MCSE.2007.55.

[Kat50] Tosio Kato. „On the Adiabatic Theorem of Quantum Mechanics“. In: Journalof the Physical Society of Japan 5.6 (1950), S. 435–439. doi: 10.1143/JPSJ.5.435.

99

Page 108: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Literaturverzeichnis

[MS02] André Martinez und Vania Sordoni. „A general reduction scheme for thetime-dependent Born-Oppenheimer approximation“. In: C. R. Math. Acad.Sci. Paris 334.3 (2002), S. 185–188. issn: 1631-073X. doi: 10.1016/S1631-073X(02)02212-4.

[Mát18] Edit Mátyus.Non-adiabatic mass-correction functions and rovibrational statesof 4He+

2 (X 2Σ+u ). Okt. 2018. arXiv: 1810.05584 [physics.chem-ph].

[NS04] Gheorghe Nenciu und Vania Sordoni. „Semiclassical limit for multistate Klein-Gordon systems: almost invariant subspaces, and scattering theory“. In: J.Math. Phys. 45.9 (2004), S. 3676–3696. issn: 0022-2488. doi: 10.1063/1.1782279.

[Oli06] Travis E. Oliphant. A guide to NumPy. USA: Trelgol Publishing. 2006.[PST03] Gianluca Panati, Herbert Spohn und Stefan Teufel. „Space-adiabatic pertur-

bation theory“. In: Adv. Theor. Math. Phys. 7.1 (2003), S. 145–204. issn:1095-0761.

[PST07] Gianluca Panati, Herbert Spohn und Stefan Teufel. „The time-dependentBorn-Oppenheimer approximation“. In: M2AN Math. Model. Numer. Anal.41.2 (2007), S. 297–314. issn: 0764-583X. doi: 10.1051/m2an:2007023.

[Ros95] Guido van Rossum. Python Tutorial. Technical Report CS-R9526. Amster-dam: Centrum voor Wiskunde en Informatica (CWI), Mai 1995.

[RS75] Michael Reed und Barry Simon.Methods of modern mathematical physics. II.Fourier analysis, self-adjointness. Academic Press, New York-London, 1975,S. xv+361.

[ST01] Herbert Spohn und Stefan Teufel. „Adiabatic decoupling and time-dependentBorn-Oppenheimer theory“. In: Comm. Math. Phys. 224.1 (2001). Dedicatedto Joel L. Lebowitz, S. 113–132. issn: 0010-3616. doi: 10.1007/s002200100535.

[Str68] Gilbert Strang. „On the construction and comparison of difference schemes“.In: SIAM J. Numer. Anal. 5 (1968), S. 506–517. issn: 0036-1429. doi: 10.1137/0705041.

[Teu03] Stefan Teufel. Adiabatic perturbation theory in quantum dynamics. Bd. 1821.Lecture Notes in Mathematics. Springer-Verlag, Berlin, 2003, S. vi+236.isbn: 3-540-40723-5. doi: 10.1007/b13355.

[TS02] Stefan Teufel und Herbert Spohn. „Semiclassical motion of dressed elec-trons“. In: Rev. Math. Phys. 14.1 (2002), S. 1–28. issn: 0129-055X. doi:10.1142/S0129055X02001077.

[Wei] Eric W. Weisstein. Hermite Polynomial. From MathWorld–A Wolfram WebResource. url: http://mathworld.wolfram.com/HermitePolynomial.html (besucht am 10. 04. 2018).

[Pyt] Python Software Foundation. Python Language Reference. Version 3.5.2. url:http://www.python.org (besucht am 26. 06. 2018).

100

Page 109: Theorie und Numerik der Born-Oppenheimer-Approximation ... · Theorie und Numerik der Born-Oppenheimer-Approximation zweiter Ordnung Daniel Weber Fachbereich Mathematik Eberhard Karls

Eigenständigkeitserklärung

Hiermit versichere ich, dass ich die vorliegende Arbeit selbstständig verfasst, keine an-deren als die angegebenen Quellen benutzt, alle wörtlich oder sinngemäß aus anderenWerken übernommenen Aussagen als solche gekennzeichnet habe und die Arbeit wedervollständig noch in wesentlichen Teilen Gegenstand eines anderen Prüfungsverfahrensgewesen ist.

Tübingen, den 25. Oktober 2018

Daniel Weber


Recommended