Berechnungen von Molekülschwingungen in der Thermodynamik

Maschinelles Lernen: Bisher zielten Anwendungen vor allem auf das Lernen von Potenzialenergiehyperflächen, derzeit zeigen sich viele neue Forschungsrichtungen. Molekülschwingungen in der Thermodynamik: Berechnung in Theorie und Simulation. Ultrakurze Laserpulse zum Beobachten molekularer Dynamik: Da die Bewegung von Elektronen genauso schnell ist wie die zeitliche Ausdehnung der Pulse, sind die Pulseigenschaften bei Simulationen explizit zu berücksichtigen.

Berechnungen von Molekülschwingungen in der Thermodynamik

Der Brückenschlag von quantenmechanisch berechneten Größen zu solchen, die für das Experiment oder praktische Anwendungen relevant sind, ist oft herausfordernd: So sind thermodynamische Größen makroskopischer Phasen häufig relevanter als reine elektronische Energien für isolierte Moleküle am Nullpunkt im Vakuum, wie sie typischerweise in quantenchemischen Rechnungen bestimmt werden. Molekülschwingungen spielen – neben Rotations- und Translationsbeiträgen – bei Standardbedingungen eine übergeordnete Rolle. Während äußerst genaue Methoden dabei helfen, das Problem der elektronischen Struktur zu lösen, werden die Beiträge der Atomkerne oft in der unzureichenden harmonischen Näherung beschrieben.

Auch darüber hinaus ist die Vorhersage genauer Schwingungsfrequenzen durch theoretische Methoden wichtig, etwa beim Entwickeln von Kraftfeldern für Molekulardynamiksimulationen oder zur Interpretation von Schwingungsspektren.¹⁾

Um Schwingungsfreiheitsgrade bei der Berechnung thermodynamischer Daten besser zu beschreiben, sind zwei Aspekte zu behandeln: zum einen die Schwingungsfrequenzen selbst, die in der quantenchemischen Rechnung aus der Hesse-Matrix bestimmt werden, und zum anderen die Schwingungszustandssumme, mit der diese Frequenzen statistisch-thermodynamisch zu thermodynamischen Größen beitragen.

Effiziente Berechnung anharmonischer Frequenzen

Es existieren mehrere Ansätze zeitunabhängiger Theorien, in denen die Kernschrödingergleichung näherungsweise gelöst wird. Analog zu den elektronischen Strukturmethoden gibt es störungstheoretische Ansätze (vibrational perturbation theory, VPT) sowie variationelle Ansätze (vibrational self-consistent field, VSCF oder vibrational configuration interaction, VCI). Allen Methoden gemein ist die Notwendigkeit, die Potenzialenergiehyperfläche (potential energy surface, PES) darzustellen. Während in der harmonischen Näherung eine Parabel in dem jeweilige Minimum zur Approximation verwendet wird, sind für anharmonische Methoden komplexere Funktionen nötig. Dies erfordert höhere Ableitungen der Energie, die rechenzeitintensiv sind. In der einfachsten Methode, Anharmonizitäten zu berücksichtigen, werden die Schwingungsmoden als voneinander entkoppelt betrachtet. Das Molekül wird entlang der jeweiligen Normalmode aus der Gleichgewichtsgeometrie ausgelenkt, und an Gitterpunkten entlang dieses Vektors werden die elektronische Struktur sowie die Hesse-Matrix bestimmt. Typischerweise sind je Mode 16 Gitterpunkte zu berechnen, der Rechenaufwand skaliert also zwar linear (mit der Zahl der Schwingungsmoden), beträgt aber bereits ein Vielfaches der harmonischen Rechnung. Um genauere Ergebnisse zu erzielen, können Terme höherer Ordnung berücksichtigt werden. Schon Zwei-Moden-Terme, in denen die Auslenkung entlang zweier Moden gleichzeitig betrachtet wird, sind für größere Moleküle nicht mehr durchführbar.

Lokale Moden

Einen Ansatz, die Kernschrödingergleichung für größere Systeme besser zu lösen als es die harmonische Näherung erlaubt, stellten Panek und Jacob vor. Darin führt eine unitäre Transformation der Normalmoden zu lokalisierten Moden. Dadurch werden die Zwei-Moden-Kopplungen klein und sind für bestimmte Paare von Moden vernachlässigbar, die auf verschiedenen Teilen des Moleküls lokalisiert sind.²⁾ Mit diesem Ansatz ließen sich Schwingungsspektren größerer Systeme berechnen, zum Beispiel Polypeptide und weitere Biomoleküle, sowie 2-D-IR-Spektren,³⁾ in denen diese Intermoden-Kopplung zentral ist.

Während diese Methode die Menge der behandelbaren Systeme signifikant erweitert, sind auch die VSCF- oder VCI-Rechnung in diesen lokalen Moden dennoch recht aufwendig. Es wäre wünschenswert, anharmonische Spektren zu bestimmen, ohne Terme höherer Ordnung behandeln zu müssen.

Maschinelles Lernen

In den vergangenen Jahren dienten Techniken des maschinellen Lernens (ML) vermehrt auch dazu, genauere Spektren vorherzusagen.⁴⁾ Dafür gibt es grundsätzlich verschiedene Methoden, die sich darin unterscheiden, welche Größen als Deskriptoren für das Training des Modells verwendet werden und welches die Zielgrößen für das maschinell erlernte Modell sind.

Eine oft verwendete Methode basiert darauf, die PES durch künstliche neuronale Netze (NN) zu bestimmen.⁵⁾ Diese können als Grundlage für Molekulardynamiksimulationen dienen, aus denen durch Fouriertransformation Spektren generiert werden,^1,6) was aber hier nicht im Fokus stehen soll. Als Deskriptoren für die ML-basierte Modellierung dienen häufig Kernkoordinaten, auch unter Einbeziehung von Symmetrie.⁷⁾ Dabei werden die Zielgrößen an bestimmten Punkten des Phasenraums quantenchemisch berechnet und der so erzeugte Datensatz dient dem Training des Modells. Mit Kenntnis der PES lassen sich dann durch numerische Lösung der Kernschrödingergleichung sowohl harmonische als auch anharmonische Frequenzen bestimmen. Durch erweiterte Techniken wie Δ-ML oder Transfer Learning sind die komplexen PES sogar auf genauem Niveau, etwa Coupled Cluster, generierbar.^8,9,10) Diese Ansätze waren in der Lage, Frequenzen mit einer Genauigkeit von nur wenigen cm^–1 verglichen mit explizit anharmonischen Rechnungen oder experimentellen Daten vorherzusagen.

Abgesehen von Deskriptoren basierend auf der Kerngeometrie lässt sich zudem die Elektronenstruktur in den Trainingsdatensatz einbeziehen, zum Beispiel die Elektronendichte⁷⁾ oder die Molekülorbitale.¹¹⁾ Auch können die Zielgrößen des maschinellen Lernens variiert werden, zum Beispiel können Frequenzen direkt – oder deren Abweichung vom Experiment¹²⁾ oder Verschiebung⁷⁾ – modelliert werden. Auch die Größen Energie, Kraft und Dipolmoment als Grundlage für die Spektrenberechnung von protonierten Wasserclustern wurden modelliert.¹³⁾

Insgesamt lassen sich Spektren durch die ML-basierte Modellierung der PES bereits zuverlässig berechnen. Eingeschränkt ist allerdings die Transferierbarkeit, da man sich schon für einzelne Systeme in einem großen Phasenraum bewegt und es zusätzlich eine so große Zahl von Systemen mit unterschiedlichen Eigenschaften gibt. Hilfreich ist dabei, dass in Kombination mit den lokalen Ansätzen sich die größeren molekularen Systeme in einzelne Fragmente zerlegen lassen, die dann separat trainiert werden.⁴⁾

Noch effizienter und somit wünschenswert wäre es, wenn sich anharmonische Frequenzen ohne die Rekonstruktion der PES direkt aus molekularen Deskriptoren vorhersagen ließen. Das Problem hierbei ist, dass wegen der Komplexität der PES die Anharmonizität einer Schwingung von vielen Faktoren abhängt: darunter das System sowie die Umgebung und die jeweils betrachtete Schwingungsmode sowie deren Kopplung mit anderen Moden.¹⁴⁾ In einer aktuellen Arbeit untersuchen wir, ob eine solche Vorhersage möglich ist, und verwenden als Deskriptoren die Molekülgeometrie, die Symmetrie des Systems sowie der jeweiligen Normalmode, die Kraftkonstante der Normalmode sowie die harmonische Frequenz.¹⁵⁾ Der Datensatz für das Training des Modells bestand aus verschieden großen Clustern der Systeme HX, CH₃X und CH₃CH₂X (X=F, Cl, Br) und umfasst 392 Datenpunkte. Diese wurden anharmonisch charakterisiert, wobei sowohl diagonale Frequenzen (ohne Berücksichtigung von Kopplungen zwischen den Moden) als auch VSCF und störungstheoretische Frequenzen (VPT2) unter Berücksichtigung dieser Kopplungen berechnet wurden.

Abbildung 1zeigt die berechneten Infrarotspektren (IR) für das System HF und vergleicht die Methoden. Die Spektren sind aus gleichverteilten HF-Clustern bis zum Hexamer zusammengesetzt. Dabei wird sowohl der Einfluss der Potenzialform einer einzelnen Oszillation deutlich (grüne versus blaue Kurve) als auch die Einbeziehung der intermodalen Kopplungen (blaue versus schwarze und rote Kurve).

Die Zielgröße des Gradient-Boosting-ML-Algorithmus war die Differenz zwischen anharmonischer und harmonischer Frequenz. Das erhaltene Modell liefert (nach fünffacher Kreuzvalidierung) gute Ergebnisse für die diagonalen Frequenzen (MAE = 23,10 cm^–1), unter Berücksichtigung der intermodalen Kopplungen sind die Ergebnisse weniger genau (VSCF: MAE = 25,54 cm^–1, VPT2: MAE = 37,99 cm^–1). Die Ursache dafür ist, dass die Deskriptoren lediglich die isolierten Normalmoden enthalten. Eine deutliche Verbesserung ist zu erwarten, wenn die Deskriptoren erweitert werden oder ein komplexerer ML-Algorithmus mit einem erweiterten Datensatz kombiniert wird.

Anharmonische Zustandssumme

Um die erhaltenen anharmonischen Frequenzen in genauere thermodynamische Größen zu übertragen, muss auch die Schwingungszustandssumme um anharmonische Terme korrigiert werden. In der harmonischen Näherung vereinfacht sich die Summation über alle (unendlich viele) Schwingungszustände des harmonischen Oszillators in einen einfachen analytischen Ausdruck. Für die anharmonischen Energieeigenwerte ist das nicht mehr ohne weiteres möglich. Gründe dafür sind die zusätzlichen Terme und die Tatsache, dass die Energieniveaus nur bis zur Dissoziation existieren.

Als Lösung lässt sich die Summation direkt ausführen. Das ist zwar aufwendiger, erlaubt es aber, beliebige anharmonische Potenziale zu verwenden. Eine weitere Möglichkeit besteht darin, genäherte anharmonische Potenziale anzuwenden, für die eine analytisch geschlossene und ähnlich einfache Form der Zustandssumme bestimmt wird. Ein Beispiel hierfür ist das Morse-Potenzial mit der Zustandssumme bei der Temperatur T, der Anharmonizitätskonstante χ_e und der Schwingungstemperatur Θ^vib des Oszillators (Gleichung rechts).¹⁶⁾

Die Kenntnis der Anharmonizitätskonstanten ist somit essenziell für die Berücksichtigung von Anharmonizitäten in der Zustandssumme. Dabei ist anzumerken, dass die obige Gleichung lediglich isolierte Morseoszillatoren betrachtet und Kopplungen zu anderen Schwingungen nicht enthalten sind.

In unserer aktuellen Studie wählen wir deshalb diagonale, anharmonische Frequenzen und die daraus bestimmten Anharmonizitätskonstanten, um aus diesen Größen in einem statistischen Ansatz – der Quantum-Cluster-Equilibrium-Methode (QCE) – thermodynamische Daten makroskopischer Phasen zu bestimmen.^17,18) Wir betrachten erneut HF, wobei wir die flüssige Phase durch verschieden große Cluster im Gleichgewicht modellieren. Im Ergebnis der QCE-Methode erhalten wir temperaturabhängige Populationen der einzelnen Cluster sowie thermodynamische Größen der flüssigen Phase und der Gasphase. Abbildung 2 zeigt das molare Volumen der flüssigen und gasförmigen Phase, wobei Clusterstrukturen bis zum Dodekamer berücksichtigt wurden. Die Schwingungsfrequenzen wurden mit den oben erwähnten Methoden bestimmt, wobei die Frequenzen der Cluster, die größer als das Pentamer sind, durch das ML-Modell bestimmt wurden. Das molare Phasenvolumen wird deutlich durch die Anharmonizitäten beeinflusst, wobei der allgemeine Kurvenverlauf erhalten bleibt.

Spannend ist die Frage nach den vorherrschenden Strukturmotiven in Fluorwasserstoff: Abbildung 3 zeigt die Häufigkeit der beiden dominanten Strukturen, dem cyclischen Pentamer und Hexamer. Während das Level der Anharmonizität das Volumen nur leicht beeinflusst, zeigen sich starke Unterschiede in den Populationen. In der flüssigen Phase dominiert das Hexamer. Unter Berücksichtigung von Anharmonizität nimmt dessen Population allerdings stark zugunsten des Pentamers ab. Ähnliche Effekte sind auch in der Gasphase zu sehen, wobei hier sogar die Kurvenverläufe variieren.

Offene Herausforderung

Anharmonische Effekte zu berücksichtigen, ist ein wichtiger Schritt hin zu genaueren thermodynamischen Daten basierend auf zuverlässigen quantenchemischen Berechnungen.¹⁹⁾ Auf molekularer Ebene gibt es dafür Verfahren für die Lösung der Kernschrödingergleichung (für kleine Systeme), Kombination dieser mit lokalen Ansätzen (für mittelgroße Systeme) sowie Techniken des maschinellen Lernens (für große Systeme). Die so gewonnenen Daten lassen sich derzeit nur eingeschränkt in thermochemische Berechnungen einbringen, da es sehr komplex ist, eine Zustandssumme zu formulieren, in der Kopplungen zwischen den einzelnen Schwingungsmoden berücksichtigt werden. Analytisch ist dies nicht ohne weiteres effizient berechenbar.

Die Verwendung der Morse-Zustandssumme ist ein erster Schritt auf dem Weg zu genaueren thermochemischen Daten.¹⁹⁾ Werden die damit konsistenten, diagonalen anharmonischen Frequenzen verwendet, sind signifikante Einflüsse insbesondere auf die vorherrschenden Strukturmotive zu beobachten. Für bessere Resultate müssen intermodale Kopplungen zwingend berücksichtigt werden. ML-basierte Ansätze, die auf der Modellierung der PES beruhen, werden hier schon häufig und äußerst erfolgreich angewandt. Die Modellierung der anharmonischen Frequenzen mit Berücksichtigung der Kopplungen direkt aus den Geometrien und harmonischen Frequenzen ist möglich und erfordert

einen umfangreicheren und diverseren Datensatz,ähnlich zu Lit. 20),fortgeschrittene Techniken des maschinellen Lernens undmehrere Schwingungsmoden als Deskriptoren einzubeziehen, statt nur die jeweils einzelne Schwingung zu betrachten.

Die Autorinnen

Eva Perlt, Jahrgang 1987, ist seit dem Jahr 2021 Nachwuchsgruppenleiterin an der Universität Jena. Sie erforscht die Rolle der Kerne in Elektronenstrukturberechnungen von Materialien, also Kernquanteneffekte sowie Anharmonizitäten. Nach Chemiestudium und Promotion in Leipzig war sie zu Forschungsaufenthalten an der Uni Bonn, der UC Irvine sowie an der FU Berlin. eva.von.domaros@uni-jena.de

Jamoliddin Khanifaev, Jahrgang 1993, ist Doktorand in der Untergruppe von Eva Perlt in der Computational Materials Science Gruppe von Marek Sierka an der Universität Jena. Er untersucht anharmonische Schwingungen mit maschinellem Lernen.

Tim Schrader, Jahrgang 1996, ist ebenfalls Doktorand in der Unterarbeitsgruppe von Eva Perlt. Er beschäftigt sich mit der Simulation von Kernquanteneffekten, insbesondere für optische Anregungen.

• ¹ M. Meuwly, CHIMIA 2022, 76(6), 589
• ² P. T. Panek, C.R. Jacob, ChemPhysChem 2014, 15,3365
• ³ J. Brüggemann, M. Wolter, C. R. Jacob, J. Chem. Phys. 2022, 157, 244107
• ⁴ R. Han, R. Ketkaew, S. Luber, J. Phys. Chem. A 2022, 126, 801
• ⁵ J. Behler, M. Parrinello, Phys. Rev. Lett. 2007, 98, 146401
• ⁶ M. Meuwly, J. Phys. Chem. B 2022, 126, 2155
• ⁷ K. Kwac, H. Freedman, M. Cho, J. Chem. Theory Comput. 2021, 17, 6353
• ⁸ R. Conte, A. Nandi, C. Qu, Q. Yu, P. L. Houston, J. M. Bowman, J. Phys. Chem. A 2022, 126, 7709
• ⁹ D. S. Rasheeda, A. M. Santa Daría, B. Schröder, E. Mátyus, J. Behler, Phys. Chem. Chem. Phys. 2022, 24, 29381
• ¹⁰ S. Käser, E. D. Boittier, M. Upadhyay, M. Meuwly, J. Chem. Theory Comput. 2021, 17, 3687
• ¹¹ F. Lu, L. Cheng, R. J. DiRisio et al., J. Phys. Chem. A 2022, 126, 4013
• ¹² J. Fu, Z. Wan, Z. Yang, L. Liu, Q. Fan, F. Xie, Y. Zhang, J. Ma, Int. J. Quant. Chem 2022, 122, e26953
• ¹³ R. Beckmann, F. Brieuc, C. Schran, D. Marx, J. Chem. Theory Comput. 2022, 18, 5492
• ¹⁴ E. Perlt, S. A. Berger, A.-M. Kelterer, B. Kirchner, J. Chem. Theory Comput. 2019, 15, 2535
• ¹⁵ J. Khanifaev, T. Schrader, E. Perlt, unveröffentlicht; eingereicht und in Begutachtung bei J. Chem. Theory Comput. 2023
• ¹⁶ M. L. Strekalov, Chem. Phys. Lett. 2004, 393, 192
• ¹⁷ B. Kirchner, F. Weinhold, J. Friedrich, E. Perlt, S. B. C. Lehmann, Quantum Cluster Equilibrium. In: Many-Electron Approaches in Physics, Chemistry and Mathematics, Cham, Springer, 2014, S. 77–96
• ¹⁸ M. von Domaros, E. Perlt, J. Ingenmey, G. Marchelli, B. Kirchner, SoftwareX 2018, 7, 356
• ¹⁹ M. von Domaros, E. Perlt, J. Chem. Phys. 2017, 146, 124114
• ²⁰ J. C. Z. Trujillo, L. K. McKemmish, J. Phys. Chem. A 2022, 126, 4100

Wissenschaft + ForschungTrendberichte