Molekulardynamische Simulation des Buckminster...

Universität RostockMathematisch-Naturwissenschaftliche FakultätInstitut für ChemieAbteilung Physikalische und Theoretische Chemie

MolekulardynamischeSimulation des Buckminster

Fullerens C60 in unterschiedlichsubstituierten aromatischen

Solventien

Arbeit zur Erlangung des akademischen Grades Bachelor of Science30. Juli 2014

Autor: Erstgutachter:Jan Neumann Prof. Dr. Ralf LudwigWiggersstraße 818057 Rostock Hansaviertel Zweitgutachter:[email protected] Dr. Dietmar Paschek

Danksagung

Ich möchte mich bei Herrn Prof. Dr. Ralph Ludwig bedanken, der mir die Möglich-keit gab, in seinem Arbeitskreis, im Bereich der Physikalischen und TheoretischenChemie, diese Arbeit anzufertigen.Außerdem gilt mein Dank in besonderem Maße Herrn Dr. Dietmar Paschek fürdas Vermitteln der Grundlagen der molekulardynamischen Simulation und die Ar-beitsweise des Arbeitskreises, für das Bereitstellen geeigneter Literatur und für dasErstellen der Simulation. Des Weiteren danke ich für das Korrekturlesen dieser Ar-beit und die daraus resultierenden Verbesserungsvorschläge und für die intensiveBetreuung während der gesamten Arbeitszeit.Abschließend danke ich meiner Familie und meinen Freunden, ohne deren Unter-stützung dieses Studium nicht möglich gewesen wäre.

2

Selbstständigkeitserklärung

Hiermit versichere ich, Jan Neumann, geboren am 30.10.1992 in Pritzwalk, dass dieseArbeit von mir selbstständig verfasst und ausgearbeitet wurde und keine anderenals die angegebenen Quellen und Hilfsmittel benutzt, sowie Zitate kenntlich gemachtwurden.

3

Inhaltsverzeichnis

Tabellenverzeichnis 5

Abbildungsverzeichnis 6

1 Einleitung 7

2 Grundlagen 8

3 Methodik 113.1 Kraftfeld für das Fulleren C60 . . . . . . . . . . . . . . . . . . . . . . 113.2 Solvens-Kraftfeld . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 113.3 Aufbau der Simulation . . . . . . . . . . . . . . . . . . . . . . . . . . 123.4 Monte-Carlo-Bestimmung des Eigenvolumens der Moleküle . . . . . . 13

4 Auswertung 144.1 Thermodynamische und kinetische Eigenschaften . . . . . . . . . . . 14

4.1.1 Dichte . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 144.1.2 Diffusionskoeffizient . . . . . . . . . . . . . . . . . . . . . . . . 15

4.2 Struktur der Solvatsphäre des C60 . . . . . . . . . . . . . . . . . . . . 174.2.1 Paarverteilungen . . . . . . . . . . . . . . . . . . . . . . . . . 174.2.2 Winkelverteilungen . . . . . . . . . . . . . . . . . . . . . . . . 234.2.3 Abstandsabhängigkeit der Winkelverteilungen . . . . . . . . . 28

4.3 Struktur des Lösungsmittels . . . . . . . . . . . . . . . . . . . . . . . 304.3.1 Paarverteilungen . . . . . . . . . . . . . . . . . . . . . . . . . 304.3.2 Winkelverteilungen . . . . . . . . . . . . . . . . . . . . . . . . 354.3.3 Abstandsabhängigkeit der Winkelverteilungen . . . . . . . . . 37

5 Zusammenfassung und Ausblick 39

Literatur 41

4

Tabellenverzeichnis

3.1 Verwendete Potenzialparameter für die intramolekularen Wechselwir-kungen des Buckminster-Fullerens C60 . . . . . . . . . . . . . . . . . 11

3.2 Parameter der nichtbindenden Wechselwirkungen der Solvensmoleküle 123.3 Verwendete Potenzialparameter für die intramolekularen Wechselwir-

kungen der Solvensmoleküle . . . . . . . . . . . . . . . . . . . . . . . 12

4.1 Temperaturabhängigkeit der Dichte der Lösungsmittel . . . . . . . . . 144.2 Vergleich der ermittelten Dichten bei 300 K mit Literaturwerten der

reinen Solvenzien bei 298 K . . . . . . . . . . . . . . . . . . . . . . . 144.3 Diffusionskoeffizienten der Solvenzien und des Fullerens bei 300 K

verglichen mit den Literaturwerten . . . . . . . . . . . . . . . . . . . 164.4 Steigerung des freien Volumens mit der Größe der Solvensmoleküle . . 174.5 Verwendete Cutoff-Radien für die Einteilung der ersten Solvatspähre

in zwei Bereiche . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 184.6 Anzahl der Nachbarmoleküle in den untersuchten Abstandsbereichen

aus Tabelle 4.5 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 214.7 Anteil der Moleküle im Bereich 1 an der Gesamtzahl der Moleküle in

der ersten Solvatsphäre des C60 . . . . . . . . . . . . . . . . . . . . . 224.8 Einteilung der ersten Solvatspähre der Solvensmoleküle in zwei Bereiche 304.9 Anzahl der Nachbarmoleküle um die Solvensmoleküle in den Berei-

chen aus Tabelle 4.8 . . . . . . . . . . . . . . . . . . . . . . . . . . . 334.10 Anteil der Moleküle im Bereich 1 an der Gesamtzahl der Moleküle in

der ersten Solvatsphäre . . . . . . . . . . . . . . . . . . . . . . . . . . 34

5

Abbildungsverzeichnis

4.1 Quadratische Verschiebung gegen die Zeit am Beispiel von C60 inBenzen bei 300 K . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15

4.2 Diffusionskoeffizienten des Solvenz (schwarz) und des Fullerens imSolvenz (rot) bei 300 K . . . . . . . . . . . . . . . . . . . . . . . . . . 16

4.3 Paarverteilungen der Solvensmoleküle um das C60 . . . . . . . . . . . 194.4 Temperaturabhängigkeit der Nachbarzahl um das C60 . . . . . . . . . 204.5 Anzahl der Moleküle in den Bereichen der verschiedenen Solvenzien

(T=300 K) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 224.6 Illustration der orientierungsdefinierenden Vektoren . . . . . . . . . . 234.7 Verteilung der Winkel in den zwei Bereichen für Benzen und Toluen . 244.8 Verteilung der Winkel in den beiden Abstandsbereichen für m-Xylen

und Mesitylen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 264.9 Fulleren mit den acht Solvensmolekülen mit dem jeweils kleinsten

Schwerpunkt-Schwerpunkt-Abstand aus zwei Perspektiven . . . . . . 274.10 Paarverteilungsfunktion der Solvensmoleküle um das C60 und Ab-

standsabhängigkeit von 〈P2〉 bei 300 K . . . . . . . . . . . . . . . . . 294.11 Solvens-Solvens-Paarverteilungen . . . . . . . . . . . . . . . . . . . . 314.12 Temperaturabhängigkeit der Nachbarzahl um die Solvensmoleküle . . 324.13 Anzahl der Moleküle in den beiden Abstandsbereichen der ersten Sol-

vathülle als Funktion der Anzahl der Substituenten . . . . . . . . . . 344.14 Verteilung der Winkel in den beiden Abstandsbereichen für Benzen

und Toluen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 354.15 Verteilung der Winkel in den beiden Abstandsbereichen für m-Xylen

und Mesitylen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 364.16 Solvens-Solvens-Paarverteilungsfunktion und Abstandsabhängigkeit

von 〈P2〉 bei 300 K . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38

5.1 Gesamte erste Solvatsphäre von Toluen um das C60 . . . . . . . . . . 40

6

1 Einleitung

Heutzutage ist die Computersimulation aus den Naturwissenschaften nicht mehrwegzudenken. Sie erlaubt Beobachtungen an Systemen, die außerhalb unserer expe-rimentellen Reichweite liegen. Weder hohe Temperaturen noch extreme Drücke sindHindernisse. Sie ermöglicht uns sowohl Einblicke in die mikroskopische Welt derAtome und Moleküle als auch ins Reich der Bewegungen der Galaxien. Außerdemkönnen Theorien durch Simulationen auf ihre Gültigkeit überprüft werden.In dieser Arbeit wird die Computersimulation in Form der Molekulardynamik zurUntersuchung des Lösungsverhaltens des Buckminsterfullerens C60 in unterschied-lichen aromatischen Lösungsmitteln eingesetzt. Fullerene sind Allotrope des Koh-lenstoffs und existieren als Nanopartikel in vielen Formen, wie zum Beispiel kugel-oder stäbchenförmig. Ihre Entdeckung im Jahr 1985 [15] wurde 1996 mit dem No-belpreis für Chemie geehrt. Sie sind seit dem durch ihre vielen Anwendungen undihre ungewöhnliche Chemie für viele Wissenschaftler von Interesse. Zum Beispieldokumentiert eine Studie die lebensverlängernde Wirkung bei Ratten durch oraleAufnahme von C60 [5]. Fullerene wurden nach ihrer Entdeckung auch in der Naturund selbst im Weltall nachgewiesen [7].Das Lösungsverhalten von C60 ist hier von besonderem Interesse, da IR-spektroskopische Untersuchungen in der Arbeitsgruppe von Prof. Dr. Ralf Ludwigdarauf hindeuten, dass planar-orientierte, stark-wechselwirkende Spezies in der ers-ten Solvatsphäre von nicht-planar-orientierten Spezies spektroskopisch unterschie-den werden können. In dieser Arbeit soll deshalb insbesondere das Orientierungs-verhalten der unterschiedlich substituierten Benzen-Derivate (Benzen, Toluen, m-Xylen, Mesitylen) in der unmittelbar benachbarten Solvathülle des C60 als Funktionder systematisch varierten Konstitution als auch als Funktion der Temperatur un-tersucht werden.

7

2 Grundlagen

Die molekulardynamische (MD) Simulation ist in den letzten Jahrzehnten zu einemwichtigen Instrument in der Untersuchung molekularer Systeme herangewachsen.Trotz der enormen Leistungsfähigkeit moderner Computer kann dabei jedoch meistnicht eine makroskopische Anzahl von Teilchen, zum Beispiel 6,022 ·1023 (1 Mol),simuliert werden. Daher wird nur ein kleiner Ausschnitt aus dem System betrachtet,der je nach Problemstellung zwischen 102 und 107 Teilchen enthält. [11]Zur Vermeidung von Oberflächeneffekten, die bei einer solchen kleinen Anzahl vonMolekülen auftreten würden, simuliert man die betrachtete Zelle unter periodischenRandbedingungen. Das heißt, dass die Zelle in allen Raumrichtungen von ihrenAbbildern umgeben ist. Dadurch treten Teilchen, die an einer Seite aus der Zelleaustreten, auf der gegenüberliegenden Seite wieder ein. Auch die Wechselwirkungder Teilchen untereinander wird durch die periodischen Randbedingungen berück-sichtigt. Die Teilchen können über die Grenzen der eigenen Zelle hinaus von Teilchenbeeinflusst werden, die sich im Wechselwirkungsradius befinden. Somit kann auchdurch eine kleine Anzahl von Molekülen ein Bulk-System simuliert werden.Des Weiteren wird angenommen, dass die Bewegung der Teilchen der klassischenMechanik gehorchen. Diese Vernachlässigung der quantenmechanischen Effekte isteine sehr gute Näherung für die meisten Probleme. Lediglich bei der Simulationsehr kleiner Teilchen, wie He oder H2, können diese Effekte meist nicht vernachläs-sigt werden. [9]Die Wechselwirkung der Teilchen wird durch Potenzialmodelle beschrieben, welchedie Grundlage der Simulation darstellen. Aus diesen Potenzialmodellen ergibt sichdie Kraft ~Fj, die auf ein Teilchen j an dem Ort ~rj wirkt als negativer Gradient desPotenzials E:

~Fj = −∇j E(~rj...~rn)

Aus diesen Kräften kann man nach Lösen der Newton’schen Bewegungsgleichungendie zeitliche Entwicklung des Systems berechnen. Dazu werden die Teilchen zu Be-ginn so in der Simulationsbox platziert, dass sie nicht überlappen. Dies wird meistdurch eine Anordnung in einem kubischen Gitter erreicht. Des Weiteren werden denTeilchen zufällig Anfangsgeschwindigkeiten zugewiesen. Dabei wird darauf geachtet,dass der Gesamtimpuls des Systems null ist.Da gegenwärtig sämtliche auf die Teilchen wirkenden Kräfte, sowie ihre Orte undGeschwindigkeiten bekannt sind, können mit Hilfe des Velocity-Verlet-Algorithmus

8

die Positionen und Geschwindigkeiten der Teilchen zu einem Zeitpunkt t + ∆t be-rechnet werden [11]:

~rj(t+ ∆t) = ~rj(t) + ~vj(t) ·∆t+1

2·~Fj(t)

mj

·∆t2

~vj(t+ ∆t) = ~vj(t) +1

2·[~Fj(t) + ~Fj(t+ ∆t)

] ∆t

mj

Dabei wird das Zeitintervall ∆t üblicherweise im Femtosekundenbereich gewählt,da sich die Stärke der Wechselwirkung zwischen den Teilchen sehr schnell ändert.Dieser Schritt wird so lange wiederholt, bis ein ausreichend langer Zeitraum unter-sucht worden ist. Aus der zeitlichen Entwicklung der Teilchenkoordinaten ergebensich die Trajektorien. Da die zu Beginn gewählten Orte und Geschwindigkeitenmeist nicht in der gewünschten Dichte und Temperatur des Systems resultieren,wird das System durch einen Equilibrierungslauf, bei dem die thermodynamischenGrößen durch Skalierung auf die gewünschten Werte gebracht werden, in dasthermodynamische Gleichgewicht gebracht. Nach dessen Ende wird die letzteKonfiguration des Systems herausgeschrieben und als Startkonfiguration für dieeigentliche Simulation verwendet.

Um die Wechselwirkung zwischen den Teilchen gut beschreiben zu können,wird sie in einen intramolekularen und einen intermolekularen Anteil unterteilt. DieGesamtwechselwirkung des Systems ergibt sich aus der Summe dieser beiden Teile.Die Beschreibung der intermolekularen Wechselwirkung erfolgt durch ein Coulomb-potenzial und durch das Lennard-Jones-(12,6)-Potenzial. Das Potenzial besitztzwei Parameter, wobei ε die Potenzialtiefe am Ort ~rε = 2

16 · σ beschreibt, während

durch σ der Stoßdurchmesser (Nullstelle der Potenzialfunktion) des Teilchenscharakterisiert wird. Die intermolekulare Wechselwirkungsenergie zwischen denMolekülen a und b mit den Wechselwirkungszentren i auf a und j auf b ergibt sichwie folgt:[12]

Eab =∑i

∑j

[1

4πε0· qiqjrij

+ 4εij

((σijrij

)12

−(σijrij

)6)]

rij = |~ri − ~rj|

Die im folgenden angegebenen Parameter beziehen sich immer auf Wechselwir-kungen mit Teilchen gleichen Atomtyps. Um die Parameter für Wechselwirkun-gen mit unterschiedlichen Atomtypen zu ermitteln, wurden die Lorentz-Berthelot-Mischungsregeln verwendet [1]:

σij =1

2(σii + σjj)

εij = (εii · εjj)12

9

2 Grundlagen

Bei der intramolekularen Wechselwirkung muss zusätzlich zur Coulomb- undLennard-Jones-Wechselwirkung noch die Bindungs-, Bindungswinkel- und die Die-derwechselwirkung berücksichtigt werden. Die Bindugswechselwirkung wird über ein"mit einem harmonischen Potenzial gekoppeltes Bindungspotenzial mit der Kraft-konstanten kb

ng und der Gleichgewichstbindungslänge r0ng" [19] beschrieben:

EBindung =∑

Bindungen

1

2kb(rij − r0ij)2

Dabei wurde die Bindungslänge auf die Gleichgewichtsbindungslängen beschränkt.Um dies zu erreichen, gibt es verschiedene Algorithmen, wie zum Beispiel SETTLE[17], SHAKE [20] und LINCS [10]. In dieser Arbeit wurde der LINCS -Alogrithmusverwendet.Ein harmonisches Bindungswinkelpotenzial mit der Kraftkonstanten ka und demGleichgewichtswinkel φ0

ijk charakterisiert die Bindungswinkelwechselwirkung:

EBindungwinkel =∑

Winkel

1

2ka(φijk − φ0

ijk)2

Das Torsionspotential um eine zentrale Bindung wird durch die nichtbindendenWechselwirkungen der beteiligten Gruppen als auch durch Diederwinkelpotentialebeschrieben. Ein Diederwinkelpotential wird hierbei in Form einer Fourierentwick-lung beschrieben [19, 1]:

EDieder =∑

Diederwinkel

kdp [1 + cos(mψijkl − ψ0ijkl)

]Hierbei skaliert kdp die Stärke der Wechselwirkung, m beschreibt die Multiplizitätund ψ0

ijkl den Phasenwinkel.Die nichtbindende, sogenannte 1-4 Wechselwirkung innerhalb eines Moleküls wird,wie die intermolekulare Wechselwirkung, durch das oben beschriebene Lennard-Jones- und Coulomb-Potenzial beschrieben. Jedoch wird diese Wechselwirkunghierbei jeweils mit dem Faktor 1

2skaliert.[12]

10

3 Methodik

Der Simulationsaufbau, der Equilibrierungslauf und die Auswertung erfolgte mitProgrammen aus dem Programmpaket MOSCITO 4.180 [18], während die Simulati-on durch das Programm GROMACS 4.6.3 [1] durchgeführt wurde. Die Auswertungder Simulationsdaten erfolgte außerdem mit Hilfe von selbstgeschriebenen Shell-Skripten.

3.1 Kraftfeld für das Fulleren C60

Zur Beschreibung der inter- und intramolekularen Wechselwirkungen des Fullerenswurde das OPLS All Atom Forcefield [12, 16, 13] von Jorgensen et al. verwendet.Dabei wurde auf jedem Kohlenstoff des C60 ein nicht geladenes Wechselwirkungszen-trum (Site) mit einer molaren Masse 12.011 g·mol−1 gelegt. Als Potenzialfunktionwurde ein Lennard-Jones-(12,6)-Potenzial mit den Parametern σ = 3, 3997 Å undε = 43.28 K verwendet.In Tabelle 3.1 sind die Bindungslänge und Bindungswinkel sowie der Diederwinkelzwischen den Kohlenstoffatomen des C60 und deren Kraftkonstanten zusammenge-fasst.

r0 [Å] kb [kJ mol−1 Å−2] φ0 ka [kJ mol−1 rad−2] ψ0 kdp [kJ mol−1] m1.4 3924.6 120◦ 527.2 180◦ 15.167 2

Tabelle 3.1: Verwendete Potenzialparameter für die intramolekularen Wechselwir-kungen des Buckminster-Fullerens C60

3.2 Solvens-Kraftfeld

Die Sites der Lösungsmittelmoleküle wurden ebenfalls auf den Atompositionenplatziert. Die betrachteten aromatischen Solvenzien sind einheitlich aufgebaut,wobei das Kohlenstoff- und die drei Wasserstoffatome der Methylgruppen alseine Teilchensorte zusammengefasst wurden. Die Parameter des Lennard-Jones-Potenzials sind in Tabelle 3.2 dargestellt.

11

3 Methodik

Teilchenart σ [Å] ε [K] M [g mol−1] q [e]C 3.55 35.23 12.0110 -0.115H 2.42 15.10 1.0079 0.115

CH3 3.80 85.55 15.035 0.115

Tabelle 3.2: Parameter der nichtbindenden Wechselwirkungen der Solvensmoleküle

Die Parameter, welche den Bindungsanteil der innermolekularen Wechselwirkungenbeschreiben, sind in Tabelle 3.3 zusammengefasst.

Bindung r0 [Å] kb [kJ mol−1 Å−2] Winkel φ0 ka [kJ mol−1 rad−2]C-C 1.4 3925 120◦ C-C-C 527.2C-H 1.08 3071 120◦ C-C-H 292.9

C-CH3 1.51 3925 120◦ C-C-CH3 527.2

Tabelle 3.3: Verwendete Potenzialparameter für die intramolekularen Wechselwir-kungen der Solvensmoleküle

Die Gleichgewichtsdiederwinkel sind, wie im C60, 180◦ mit einer Kraftkonstante vonkdi = 15.178 kJ mol−1 und einer Multiplizität m von 2.

3.3 Aufbau der Simulation

Um das System aus C60 und Solvens aufzubauen, wurde zunächst ein reines Solvens-System erzeugt und durch eine kurze MD-Simulation equilibriert. Anschließend wur-de das C60-Molekül in der Box platziert und sämtliche überlappenden Solvensmo-leküle entfernt. Hierbei muss darauf geachtet werden, dass kein Solvensmolekül indas Innere des Fullerens gelangt. Dies wurde dadurch sichergestellt, das im Innerendes Fullerens ein Argonatom platziert und nach Beendigung der Zusammenführungsofort wieder entfernt wurde.Eine Simulationsbox besteht aus einem C60 und 429 Solvensmolekülen. Der Equili-brierungslauf wurde durchgeführt und die Endkoordinaten als Anfangskoordinatendes eigentlichen Produktionslaufes verwendet. Das System wurde in einem Tempe-raturbereich von 260 K bis 360 K mit Abständen von jeweils 20 K simuliert. AlsZeitintervall ∆t wurden 2 fs gewählt. Nach 2 ns erfolgte das Erstellen der Output-dateien, in denen Energien und die Konfigurationen des Systems zusammengefasstsind. Die Gesamtlänge der Simulation betrug im Falle des Toluens 60 ns und beiallen anderen Solvenzien 40 ns. Die Temperaturkontrolle erfolgte über das Noose-Hoover-Thermostat und die Druckkontrolle über das Parrinello-Rahman-Barostat.

12

3.4 Monte-Carlo-Bestimmung des Eigenvolumens der Moleküle

3.4 Monte-Carlo-Bestimmung des Eigenvolumensder Moleküle

Um das Eigenvolumen jeweils einzelner Moleküle zu ermitteln, wurde ein Monte-Carlo-Verfahren angewendet, welches im folgenden kurz beschrieben werden soll:Zunächst wurde ein Referenzvolumen definiert, welches die von uns untersuchtenMoleküle vollständig einschließt. Hierzu wurde jeweils die maximale beziehungsweiseminimale Koordinate in jeder Raumrichtung aus dem Satz der Atomkoordinaten desReferenzmoleküls ermittelt. Zu dieser Koordinate wurde der entsprechende halbeStoßdurchmesser σii

2zuzüglich eines Offsets von 0.1 Å addiert (maximaler Wert)

beziehungsweise abgezogen (minimaler Wert).Dieses Referenzvolumen wurde danach mit Hilfe von insgesamt 108 zufällig erzeugtenKoordinaten abgetastet. Eine ausgewählte Koordinate wird als sich innerhalb desMolekülvolumens befindend definiert, wenn der Abstand zu einem der Atome kleinerist als der entsprechende halbe Stoßdurchmesser.Aus dem Verhältnis der Koordinaten innerhalb des Molekülvolumens bezogen aufsämtliche Koordinaten lässt sich das Molekülvolumen ermitteln:

Vmol =nmol

nges · Vref

Die Zufallszahlen wurden mit Hilfe des in gFortran eingebauten Pseudozufallszahlen-generators ermittelt. Durch Variation des Startwertes (seed-value) wurde der relativestatistische Fehler der Volumenbestimmung als <0.1% ermittelt.Um im Fall des Fullerens zu verhindern, dass das Volumen innerhalb des Käfigsmitgezählt wird, wurde ein zusätzliches Wechselwirkungszentrum im Schwerpunktdes Moleküls platziert, dessen Stoßdurchmesser dem Durchmesser des C60-Molekülsentspricht.

13

4 Auswertung

4.1 Thermodynamische und kinetischeEigenschaften

4.1.1 Dichte

Die Temperaturabhängigkeit der Dichte ist in Tabelle 4.1 zu sehen. Wie erwartetnimmt die Dichte mit steigender Temperatur ab. Des Weiteren nimmt sie ebenfallsmit der Zahl der Substituenten ab.

T [K] ρ [kg m−3]Benzen Toluen m-Xylen Mesitylen

260 922.40 892.01 881.18 875.71280 897.85 871.85 863.73 858.51300 875.01 852.31 845.37 841.64320 850.83 830.53 826.86 824.78340 826.76 811.05 809.13 805.53360 801.69 787.91 790.32 786.13

Tabelle 4.1: Temperaturabhängigkeit der Dichte der Lösungsmittel

In Tabelle 4.2 wurde die ermittelten Dichten mit Literaturdaten der reinen Solven-zien verglichen.

Solvens ρ [kg m−3]Simulation Literatur

Benzen 875.01 873.61Toluen 852.31 862.19m-Xylen 845.37 860.07Mesitylen 841.64 861.12

Tabelle 4.2: Vergleich der ermittelten Dichten bei 300 K mit Literaturwerten derreinen Solvenzien bei 298 K

Die Daten des m-Xylen und des Mesitylen zeigen hierbei größere Abweichungen,welche allerdings durch die Anwesenheit des C60 erklärt werden können.

14

4.1 Thermodynamische und kinetische Eigenschaften

4.1.2 Diffusionskoeffizient

Der Diffusionskoeffizient ist definiert durch die Stokes-Einsteinbeziehung:

D =kBT

6πηR

Dabei ist η die dynamische Viskosität des Solvens und R der hydrodynamischeRadius des diffundierenden Teilchens.Durch die mittlere quadratische Verschiebung der Molekülschwerpunkte〈|~ri(t)− ~ri(t+ τ)|2〉 lässt sich der Diffusionskoeffizient aus den Simulationsda-ten ermitteln:

D =1

6limτ→∞

∂

∂τ

⟨|~ri(t)− ~ri(t+ τ)|2

⟩

0 100 200 300 400t [ps]

0

1

2

3

4

<|r

i(t)-

r i(t+

τ)|

2>

Benzen

C60

Abbildung 4.1: Quadratische Verschiebung gegen die Zeit am Beispiel von C60 inBenzen bei 300 K

In Abbildung 4.1 ist die quadratische Verschiebung der Molekülschwerpunkte gegendie Zeit aufgetragen. Der minimale Zeitschritt (Zeit zwischen zwei abgespeichertenKonfigurationen) betrug dabei 0.2 ps. Um die Diffusionskoeffizienten zu bestimmen,wurde nun ein Polynom erster Ordnung an die letzten 1500 Punkte des Graphengefittet und aus dem Anstieg der Diffusionskoeffizient bestimmt.In Abbidung 4.2 und Tabelle 4.3 sind Diffusionskoeffizienten der Solvenzienund des C60 dargestellt. Es ist auffällig, dass die Werte zum Teil stark von denLiteraturwerten abweichen. Diese Abweichungen von experimentellen Werten sindjedoch keine Seltenheit für MD-Simulationen, da die Diffusionskoeffizienten sehrempfindlich auf Details der intermolekularen Wechselwirkungen reagieren. [6, 4].

15

4 Auswertung

Des Weiteren würde es den Erwartungen entsprechen, wenn der Diffusionskoeffizientdes Solvenzmoleküls abnehmen würde, je größer das Solvenzmolekül ist. Da diesoffensichtlich nicht der Fall ist wirkt diesem Trend offenbar noch ein anderer Trendentgegen.

0

0.5

1

1.5

2

2.5D

[1

0-9

m2 s

-1]

Solvens

C60

Benzen Toluen m-Xylen Mesitylen

Abbildung 4.2: Diffusionskoeffizienten des Solvenz (schwarz) und des Fullerens imSolvenz (rot) bei 300 K. Die gestrichelte Linie dient der Illustra-tion und soll in keiner Weise einen physikalischen Zusammenhangdarstellen.

Solvenz DSolvenz [10−9m2s−1] DSolvenz,Lit. [10−9m2s−1] DC60 [10−9m2s−1]Benzen 1.31 2.02 (293 K) [2] 0.28Toluen 2.43 2.27 (298 K) [2] 0.62m-Xylen 2.39 2.60 (313 K) [21] 0.77Mesitylen 2.40 1.58 (298 K) [8] 1.00

Tabelle 4.3: Diffusionskoeffizienten der Solvenzien und des Fullerens bei 300 K ver-glichen mit den Literaturwerten

Tatsächlich nimmt das freie Volumen in der Simulationsbox mit der Größe der Sol-vensmoleküle zu (Tabelle 4.4).

16

4.2 Struktur der Solvatsphäre des C60

j VBox, 300 K [nm3] Vj [Å3] Vf [nm3]C60 - 525.24 -

Benzen 64.96 89.04 26.24Toluen 78.42 100.93 34.60m-Xylen 90.89 112.83 41.96Mesitylen 103.16 124.74 49.12

Tabelle 4.4: Steigerung des freien Volumens mit der Größe der Solvensmoleküle

Das freie Volumen wurde mit folgender Gleichung bestimmt:

Vf = VBox −∑j

Nj · Vj

Wobei VBox das Volumen der Simulationsbox, Nj die Anzahl der Moleküle vomTypen j und Vj das von uns ermittelte Eigenvolumen des Moleküls j ist.Das Eigenvolumen der Moleküle wurde mittels eines einfachen Monte-Carlo-Verfahrens berechnet. Das Volumen der Box wurde als Mittelwert aus denSimulationsdaten entnommen. Jede Simulation besteht aus 429 Solvenzmolekülenund einem Fulleren.

Durch die Zunahme des freien Volumens haben die Moleküle mehr Platz,was zu einem größeren Diffusionskoeffizienten führt. Dieser Effekt kompensiert dieAbnahme der Diffusionskoeffizienten mit der Größe der Moleküle.


4.2.1 Paarverteilungen

Die radiale Paarverteilungsfunktion beschreibt die lokalen Dichteunterschiede desSystems als Funktion des Abstandes zu einem betrachteten Teilchen. Sie ist definiertals lokale Dichte im Abstand r und normiert auf die mittlere Dichte des Systems.

g(r) =ρi(r, δr)

< ρi >

Beziehungsweise [3]:

g(r) =V

N2

⟨∑i

∑j 6=i

δ(|~r − ~ri|)

⟩In Abbildung 4.11 sind die Paarverteilung der Solvensmoleküle um das Fullerenbei verschiedenen Temperaturen dargestellt. Um den Temperaturverlauf zu ver-deutlichen, sind im unteren Teil der Graphen außerdem Differenzpaarverteilungen

17

4 Auswertung

enthalten. Dabei wurden die Funktionswerte der Paarverteilung bei 260 K vondenen bei 300 K abgezogen.Man kann deutlich erkennen, dass sich die erste Solvatsphäre (bis ca. 1 nm) beiallen Solvenzien in zwei Bereiche aufteilen lässt, die jeweils durch die gestricheltenLinien gekennzeichnet sind und die im Folgenden als Bereich 1 und Bereich 2bezeichnet werden. Der Abstand, bei welchem der Übergang zwischen den beidenBereichen stattfindet, ist leicht temperaturabhängig. Dies wird jedoch aus Gründender besseren Vergleichbarkeit im Folgenden vernachlässigt (Tabelle 4.5).

Die zwei getrennten Maxima in der ersten Solvatsphäre lassen auf zwei un-terschiedliche Anordnungen der Solvensmoleküle um das C60 schließen.Es wird weiterhin beobachtet, dass mit steigender Temperatur das Maximum imBereich 2 auf Kosten des Maximums im Bereich 1 größer wird. Dies lässt auf einenmöglichen unterschiedlichen Entropieeinfluss in der Population beider Bereicheschließen: Je näher ein Solvensmolekül dem C60 kommt, desto eingeschränkter ist esin seinen Orientierungsmöglichkeiten und desto kleiner wird offenbar seine Entropie.

Bereich r [nm] (Benzen) r [nm] (Toluen) r [nm] (m-Xylen) r [nm] (Mesitylen)1 bis 0.73 bis 0.74 bis 0.74 bis 0.742 bis 0.96 bis 0.98 bis 1.01 bis 1.01

Tabelle 4.5: Verwendete Cutoff-Radien für die Einteilung der ersten Solvatspähre inzwei Bereiche

In Abbildung 4.4 und Tabelle 4.6 ist die Temperaturabhängigkeit der Anzahl derNachbarmoleküle in den Bereichen dargestellt. Auch hier wird beobachtet, dass mitsteigender Temperatur die Gesamtzahl der Moleküle in der ersten Solvatsphäre ab-nimmt. Dieser Effekt ist grundsätzlich auf die Verringerung der Dichte der Lösungs-mittel mit steigender Temperatur zurückzuführen.Bei der Betrachtung des Anteils der Moleküle im Bereich 1 der ersten Solvatsphärewird der oben vermutete Entropieeffekt besonders deutlich. Wie in Tabelle 4.7 dar-gestellt ist, nimmt dieser Anteil über das betrachtete Temperaturintervall von 260 Kbis 360 K etwa zwischen 5% und 10% ab. Hierbei fällt auf, dass die Temperatur-änderung im Falle der symmetrischen Moleküle Benzen und Mesitylen (etwa 10%)stärker ausfällt als im Falle von Toluen und m-Xylen (etwa 5%).

18


0 0.5 1 1.5 2r [nm]

-0.2

-0.1

0

0.1

g(r

) 300 K

- g

(r) 2

60 K 0

1

2

3

g(r

)

260 K

280 K

300 K

320 K

340 K

360 K

(a) Benzen

0 0.5 1 1.5 2r [nm]

-0.1

-0.05

0

0.05g

(r) 3

00 K

- g

(r) 2

60 K 0

0.5

1

1.5

2

2.5

g(r

)

260 K

280 K

300 K

320 K

340 K

360 K

(b) Toluen

0 0.5 1 1.5 2r [nm]

-0.1

-0.05

0

0.05

g(r

) 300 K

- g

(r) 2

60 K 0

0.5

1

1.5

2

2.5

g(r

)

260 K

280 K

300 K

320 K

340 K

360 K

(c) m-Xylen

0 0.5 1 1.5 2r [nm]

-0.1

-0.05

0

0.05

g(r

) 300 K

- g

(r) 2

60 K 0

0.5

1

1.5

2

2.5

g(r

)

260 K

280 K

300 K

320 K

340 K

360 K

(d) Mesitylen

Abbildung 4.3: Paarverteilungen der Solvensmoleküle um das C60

19

4 Auswertung

260 280 300 320 340 360T [K]

0

5

10

15

20

25N

Bereich 1 Bereich 2

Gesamt

(a) Benzen

260 280 300 320 340 360T [K]

0

5

10

15

20

25

N

Bereich 1 Bereich 2

Gesamt

(b) Toluen

260 280 300 320 340 360T [K]

0

5

10

15

20

25

N

Bereich 1 Bereich 2

Gesamt

(c) m-Xylen

260 280 300 320 340 360T [K]

0

5

10

15

20

25

N

Bereich 1 Bereich 2

Gesamt

(d) Mesitylen

Abbildung 4.4: Temperaturabhängigkeit der Nachbarzahl um das C60 aus Tabelle4.6

20


(a) Benzen

T [K] Bereich 1 Bereich 2 Gesamt260 4.59 19.34 23.93280 4.48 19.00 23.49300 4.24 18.92 23.17320 4.12 18.59 22.71340 3.96 18.23 22.18360 3.76 18.07 21.83

(b) Toluen


(c) m-Xylen


(d) Mesitylen


Tabelle 4.6: Anzahl der Nachbarmoleküle in den untersuchten Abstandsbereichenaus Tabelle 4.5

21

4 Auswertung

T [K] Benzen Toluen m-Xylen Mesitylen260 0.192 0.243 0.253 0.277280 0.191 0.243 0.258 0.265300 0.183 0.239 0.245 0.272320 0.181 0.235 0.239 0.257340 0.179 0.232 0.237 0.258360 0.172 0.231 0.240 0.244

Tabelle 4.7: Anteil der Moleküle im Bereich 1 an der Gesamtzahl der Moleküle inder ersten Solvatsphäre des C60

Bei dieser Betrachtung wird außerdem ein zweiter lösungsmittelabhängiger Effektdeutlich. Zwar sinkt die Gesamtzahl der Moleküle in der ersten Solvatsphäre mitzunehmender Größe des Solvensmoleküls, jedoch nimmt der Anteil der Moleküleim Bereich 1 zu. Dies führt dazu, dass die Anzahl der Moleküle im Bereich 1 beieiner Temperatur für alle Lösungsmittel nahezu konstant bleibt (Abbildung 4.5),während die Gesamtzahl der Moleküle in der ersten Solvatsphäre deutlich abfällt.Dieser Abfall kann möglicherweise durch die Anordnung größer werdender Moleküleum ein Zentralteilchen konstanter Größe erklärt werden. Der steigende Anteil derMoleküle im Bereich 1 lässt darauf schließen, dass mit zunehmender Größe derSolvensmoleküle die Wechselwirkung mit dem Fulleren stärker wird.

0

5

10

15

20

25

30

N

Bereich 1 Bereich 2

Gesamt


Abbildung 4.5: Anzahl der Moleküle in den Bereichen der verschie-denen Solvenzien (T=300 K). Die gestrichelte Liniedient nur der Illustration und keinen physikalischenZusammenhang darstellen.

22


4.2.2 Winkelverteilungen

Um die Orientierung der Solvensmoleküle um das C60 zu charakterisieren, wurdendrei Vektoren definiert (Abbildung 4.6). Der erste Vektor liegt innerhalb derMolekülebene und verbindet das vierte mit dem ersten Kohlenstoffatom innerhalbdes aromatischen Ringsystems, an dem im Falle des Toluens die Methylgruppegebunden ist. Dieser Vektor wird im Weiteren als innerer Vektor bezeichnet.Der zweite Vektor ist der Normalenvektor der Molekülebene. Er wird durch dasKreuzprodukt zweier innermolekularer Vektoren definiert. Der eine Vektor verläuftvom C1 zum C3 und der andere vom C1 zum C5. Das Kreuzprodukt dieser beidenVektoren ergibt einen Vektor, der zu beiden senkrecht steht, also genau denNormalenvektor der Molekülebene definiert. Anschließend wurde die Verteilungdes Kosinus des Winkels θ ausgewertet, welchen diese Vektoren mit dem drittenVektor, dem Schwerpunktverbindungsvektor zwischen Solvensmolekül und C60,einschließen. (Abbildung 4.7).

Abbildung 4.6: Illustration der orientierungsdefinierenden Vektoren. Hier am Bei-spiel von Benzen. Zu sehen sind der innere Vektor (blau), der Norma-lenvektor (violett) und der Schwerpunktverbindungsvektor (grün).

Im Falle einer statistischen Verteilung der betrachteten Vektoren wäre die Vertei-lung des Kosinus des Winkels eine Konstante. Dies würde bedeuten, dass es keineVorzugsorientierung gibt.Bei der berechneten Verteilung der Winkel für das Benzen (dargestellt in Abbildung4.7) ist erkennbar, dass die Verteilung beider Vektoren im Bereich 1 sehr geordnetsind. Der innere Vektor bevorzugt eine Konfiguration senkrecht zum Schwer-punktverbindungsvektor, während der Normalenvektor parallel orientiert ist. DieMolekülebene des Benzens ist also parallel zur Oberfläche des Fullerens angeordnet.Diese Anordnung wird auch bei allen anderen Solvensmolekülen gefunden.Im Bereich 2 wird im Gegensatz dazu eine senkrechte Anordnung bevorzugt, bei derdie Winkel jedoch über einen deutlich größeren Bereich variieren. Die wahrschein-lichste Orientierung des Normalenvektors ist allerdings nicht parallel oder senkrecht

23

4 Auswertung

zum Verbindungsvektor, sondern stellt offenbar eine gekippte Anordnung dar. Derinnere Vektor bevorzugt hierbei eine Ausrichtung senkrecht (mittleres Maximum)und rund 35◦ bzw. 145◦ (äußere Maxima) zum Schwerpunktverbindungsvektor. DieBenzenmoleküle zeigen mit den C-H-Bindungen auf einer Seite in Richtung auf dasFulleren.

Für die Verteilung der Winkel für Toluen finden wir ein ähnliches Verhalten.Im Gegensatz zum Benzen ist die Verteilung für den inneren Vektor im Bereich 1leicht zu positiven Werten verschoben. Zwar richtet sich auch das Toluenmolekülparallel zur Oberfläche des Fullerens aus, aber die Methylgruppe orientiert sichgeringfügig zum C60. Dies ermöglicht offenbar eine stärkere Wechselwirkung mitdem Fulleren.Auch im Bereich 2 können Abweichungen von der Benzenverteilung beobachtetwerden. Der innere Vektor ist hier zu stärker negativen Werten verschoben undweist weniger ausgeprägte Maxima bei cos θ ≈ −0.09 (θ = 95◦) und cos θ ≈ 0.63(θ = 51◦) auf. Das Toluenmolekül zeigt also ebenfalls mit einer Bindung in Richtungauf das Fulleren. Jedoch sind Konfigurationen, bei denen die Methylgruppe vomC60 wegzeigt, bevorzugt.

-1 -0.5 0 0.5 1cosθ

0

1

2

3

4

P(c

os

θ)

0

0.5

1

1.5

2

P(c

os

θ)

Benzen

Toluen

(a) Bereich 1

-1 -0.5 0 0.5 1cosθ

0

0.2

0.4

0.6

P(c

os

θ)

0.2

0.4

0.6

0.8

P(c

os

θ)

Benzen

Toluen

(b) Bereich 2

Abbildung 4.7: Verteilung der Winkel in den zwei Bereichen für Benzen und Toluen.Im oberen Teil der Abbildung ist die Verteilung des inneren Vektorsund im unteren Teil die Verteilung des Normalenvektors dargestellt.

24


Bei den Verteilungen für den Bereich 1 setzt sich der im Falle des Toluens beobachte-te Trend für m-Xylen und Mesitylen weiter fort (Abbildung 4.8): Der Normalenvek-tor ist weiterhin parallel zum Schwerpunktverbindungsvektor orientiert, während derinnere Vektor bevorzugt einen Winkel einnimmt, der geringfügig kleiner als 90◦ ist.Hierbei versuchen die Methylgruppen analog zum Toluen sich zum Fulleren hin zuorientieren.Auch im Bereich 2 setzt sich der Trend für die Verteilung der Orientierung desNormalenvektors fort. Er richtet sich weniger senkrecht zum Schwerpunktverbin-dungsvektor aus sondern bevorzugt stattdessen einen 45◦-Winkel (oder 135◦). Beider Verteilung des inneren Vektors gibt es jedoch nur noch zwei Maxima.Für m-Xylen liegt das erste Maximum bei cos θ ≈ −0, 8 (θ = 143◦) und das zweiteMaximum bei cos θ ≈ 0.43 (θ = 65◦). Das Molekül steht im Gegensatz zum Benzenund Toluen nicht mehr auf einer Bindung, sondern auf einem Atom aus dem Ring.Dabei ordnet es sich entweder so an, dass beide Molekülgruppen in Richtung aufdas Fulleren zeigen (θ = 65◦) oder beide vom Fulleren weg zeigen (θ = 143◦).Die Anordnung der Mesitylenmoleküle kann dabei wie folgt erklärt werden: Auf-grund der Symmetrie des Mesitylen-Moleküls gibt es keine unterscheidbaren ge-kippten Anordnungen, daher existiert nur noch ein Maximum. Das Molekül ordnetsich so an, dass es auf einem Ringkohlenstoffatom steht, welches sich zwischen zweiMethylgruppen befindet. Damit zeigen die zwei benachbarten Methylgruppen etwamit einem Winkel von rund 60◦ in Richtung auf das Fulleren, während die dritteMethylgruppe in die entgegengesetzte Richtung ausgerichtet ist (180◦).

25

4 Auswertung

-1 -0.5 0 0.5 1cosθ

0

1

2

3

4

5

6

P(c

os

θ)

0

0.5

1

1.5

2

P(c

os

θ)

m-Xylen

Mesitylen

(a) Bereich 1

-1 -0.5 0 0.5 1cosθ

0.1

0.2

0.3

0.4

0.5

0.6

0.7

P(c

os

θ)

0.4

0.6

0.8

P(c

os

θ)

m-Xylen

Mesitylen

(b) Bereich 2

Abbildung 4.8: Verteilung der Winkel in den beiden Abstandsbereichen für m-Xylenund Mesitylen. Im oberen Teil der Abbildung ist die Verteilung desinneren Vektors und im unteren Teil die Verteilung des Normalen-vektors dargestellt.

Zusammenfassend lässt sich für alle Solvenzien sagen, dass die erste Solvatsphäreaus zwei Bereichen mit hoher Orientierungsordnung besteht. Im Bereich 1 bevorzu-gen die Moleküle eine Anordnung parallel zur Oberfläche des Fullerens, während dieMoleküle im Bereich 2 sich dazu senkrecht ausrichten.Um dies zu verdeutlichen, haben wir ausgewählte Konfigurationen des C60 inklusiveeiniger Solvensmoleküle in Abbildung 5.1 dargestellt. Die Abbildungen wurden mitHilfe des Programmes VMD aus den Simulationsdaten erstellt. Deutlich zu erkennenist, dass die vier Moleküle, welche den kürzesten Abstand zum Fulleren aufweisen(hier rot dargestellt), parallel zur Oberfläche des Fullerens ausgerichtet sind, wäh-rend die vier in der Abstandshierarchie folgenden Moleküle eine stärker senkrech-te Anordnung bevorzugen. Auch die getroffene Aussage über die Ausrichtung derMesitylenmoleküle im Bereich 2 lässt sich beobachten: Keines der vier hellblauenSolvenzmoleküle aus Bereich 2 zeigt mit einer Methylgruppe in Richtung auf dasFulleren.

26


(a) Benzen

(b) Toluen

(c) m-Xylen

(d) Mesitylen

Abbildung 4.9: Fulleren mit den acht Solvensmolekülen mit dem jeweils kleinstenSchwerpunkt-Schwerpunkt-Abstand aus zwei Perspektiven. Die vierroten Moleküle kommen dem C60 dabei jeweils am nächsten.

27

4 Auswertung

4.2.3 Abstandsabhängigkeit der Winkelverteilungen

Um die Abstandsabhängigkeit der Winkelverteilungen zu charakterisieren, wurdedas zweite Legendre-Polynom auf den Kosinus des Winkels θ zwischen den Vektorenangewendet:

P2 =1

2

(3 cos2 θ − 1

)Durch eine Mittelwertbildung über die Orientierung aller Vektoren lässt sich durchdiese Funktion parallele und senkrechte Orientierungen der Vektoren von einerGleichverteilung unterscheiden. Der Mittelwert 〈P2〉 nimmt den Wert 1 für par-allele Anordnungen (cos θ = ±1) und den Wert −1

2für senkrechte Anordnungen

(cos θ = 0) an, wobei leichte Abweichungen von diesen Werten auf die Mittelwert-bildung über viele Moleküle zurückzuführen sind. Bei einer Gleichverteilung derWinkel nimmt der Mittelwert der Funktion den Wert 0 ein.Der Funktionswert an der Stelle r0 ergibt sich aus dem Mittelwert der Funktion überalle Moleküle, die sich in einem Abstandsbereich von r0 bis r0 + 0.05 nm befinden.Wie in Abbildung 4.10 zu erkennen ist, stimmt der Verlauf der Funktion gut mitdem Verlauf der Paarverteilungsfunktionen überein. Während für kleine Abständeeine hohe Ordnung sichtbar ist, geht die Funktion für große Abstände gegen null.Besonders bei Mesitylen ist erkennbar, dass sich die Orientierungsverteilung offenbarauch nach dem Ende der ersten Solvatsphäre im geringen Maße fortsetzt.Der Funktionsverlauf bestätigt die bisherigen Vermutungen. Während die Moleküleim Bereich 1 parallel zum Fulleren orientiert sind, ordnen sie sich im Bereich 2 ehersenkrecht an. Es wurde keine Temperaturabhängigkeit des Orientierungsverhaltensgefunden.

28

0 0.5 1 1.5 2r [nm]

-0.4

0

0.4

0.8



innerer Vektor

Normalenvektor

0

0.5

1

1.5

2

2.5

3

3.5g

(r)

260 K

280 K

300 K

320 K

340 K

360 K

(a) Benzen

0 0.5 1 1.5 2r [nm]

-0.4

0

0.4

0.8<

P2>

innerer Vektor

Normalenvektor

0

0.5

1

1.5

2

2.5

g(r

)

260 K

280 K

300 K

320 K

340 K

360 K

(b) Toluen

0 0.5 1 1.5 2r [nm]

-0.4

0

0.4

0.8



innerer Vektor

Normalenvektor

0

0.5

1

1.5

2

2.5

g(r

)

260 K

280 K

300 K

320 K

340 K

360 K

(c) m-Xylen

0

0.5

1

1.5

2

2.5

g(r

)

260 K

280 K

300 K

320 K

340 K

360 K

0 0.5 1 1.5 2r [nm]

-0.4

0

0.4

0.8



innerer Vektor

Normalenvektor

(d) Mesitylen

Abbildung 4.10: Paarverteilungsfunktion der Solvensmoleküle um das C60 und Ab-standsabhängigkeit von 〈P2〉 bei 300 K

29

4 Auswertung

4.3 Struktur des Lösungsmittels

4.3.1 Paarverteilungen

Im Gegensatz zu den Paarverteilungen um das Fulleren zeigen die Solvens-Solvens-Paarverteilungen (Abbildung 4.11) keine zwei voneinander getrennten Bereiche inder ersten Solvatsphäre auf. Da jedoch davon auszugehen ist, dass es trotzdem einähnliches Orientierungsverhalten wie im Falle des C60 existiert, wurden auch hierzwei Abstandsbereiche willkürlich definiert. Der Bereich 1 geht hierbei bis zum Ma-ximum der Paarverteilung und der Bereich 2 geht von dort aus bis zum Ende derersten Solvatsphäre (Minimum der Paarverteilungsfunktion) (Tabelle 4.8).Nur bei der Paarverteilung des Mesitylens ist bei niedrigen Temperaturen eine Schul-ter zu erkennen. Der Bereich 1 für Mesitylen wurde daher nicht an der Lage desMaximums, sondern an der Lage dieser Schulter festgesetzt.

Bereich r [nm] (Benzen) r [nm] (Toluen) r [nm] (m-Xylen) r [nm] (Mesitylen)1 bis 0.57 bis 0.60 bis 0.62 bis 0.58 (0.64)2 bis 0.77 bis 0.83 bis 0.88 bis 0.95

Tabelle 4.8: Einteilung der ersten Solvatspähre der Solvensmoleküle in zwei Bereiche.Für Mesitylen wurde die Lage des Maximums in Klammern hinzugefügt.

Bei der temperaturabhängigen Betrachtung der Nachbarzahl (Abbildung 4.12) kön-nen genau die gleichen Beobachtungen wie beim C60 angestellt werden. Mit stei-gender Temperatur nimmt die Gesamtzahl der Moleküle in der ersten Solvatsphäreaufgrund der Dichteänderung ab.

30


0 0.5 1 1.5 2r [nm]

-0.2

-0.1

0

0.1

g(r

) 300 K

- g

(r) 2

60 K 0

0.5

1

1.5

2

2.5g

(r)

260 K

280 K

300 K

320 K

340 K

360 K

(a) Benzen

0 0.5 1 1.5 2r [nm]

-0.1

-0.05

0

0.05g

(r) 3

00 K

- g

(r) 2

60 K 0

0.5

1

1.5

2

g(r

)

260 K

280 K

300 K

320 K

340 K

360 K

(b) Toluen

0 0.5 1 1.5 2r [nm]

-0.08

-0.04

0

0.04

g(r

) 300 K

- g

(r) 2

60 K 0

0.5

1

1.5

g(r

) 260 K

280 K

300 K

320 K

340 K

360 K

(c) m-Xylen

0 0.5 1 1.5 2r [nm]

-0.08

-0.04

0

g(r

) 300 K

- g

(r) 2

60 K 0

0.5

1

1.5

g(r

) 260 K

280 K

300 K

320 K

340 K

360 K

(d) Mesitylen

Abbildung 4.11: Solvens-Solvens-Paarverteilungen

31

4 Auswertung

260 280 300 320 340 360T [K]

0

5

10

15N

Bereich 1 Bereich 2

Gesamt

(a) Benzen

260 280 300 320 340 360T [K]

0

5

10

15

N

Bereich 1 Bereich 2

Gesamt

(b) Toluen

260 280 300 320 340 360T [K]

0

5

10

15

N

Bereich 1 Bereich 2

Gesamt

(c) m-Xylen

260 280 300 320 340 360T [K]

0

5

10

15

N

Bereich 1 Bereich 2

Gesamt

(d) Mesitylen

Abbildung 4.12: Temperaturabhängigkeit der Anzahl der nächsten Nachbarn um dieSolvensmoleküle aus Tabelle 4.9

32


(a) Benzen


(b) Toluen


(c) m-Xylen


(d) Mesitylen


Tabelle 4.9: Anzahl der Nachbarmoleküle um die Solvensmoleküle in den Bereichenaus Tabelle 4.8

Beim Anteil der Moleküle im Bereich 1 an der Gesamtzahl der Moleküle in der erstenSolvatsphäre (Tabelle 4.10) kann jedoch ein Unterschied zum beim C60 beobachtetenVerhalten festgestellt werden. Im Gegensatz zur Situation beim Fulleren nimmt derAnteil der Moleküle in der ersten Solvatsphäre mit zunehmender Größe der Moleküle

33

4 Auswertung

ab. Dies lässt vermuten, dass die Solvens-Solvens-Wechselwirkung mit zunehmenderMolekülgröße offenbar schwächer wird.

T [K] Benzen Toluen m-Xylen Mesitylen260 0.358 0.321 0.273 0.143280 0.343 0.312 0.267 0.140300 0.331 0.304 0.261 0.136320 0.320 0.296 0.257 0.133340 0.309 0.291 0.259 0.130360 0.300 0.284 0.248 0.128

Tabelle 4.10: Anteil der Moleküle im Bereich 1 an der Gesamtzahl der Moleküle inder ersten Solvatsphäre

Bei der Betrachtung der Anzahl der Nachbarmoleküle als Funktion der Anzahl derSubstituenten (Abbildung 4.13) wird ein weiterer Unterschied zu den Auftragungenbeim C60 deutlich. Während dort die Gesamtzahl stark abnimmt, bleibt sie hiernahezu konstant. Dies erscheint auch plausibel, da man hier im Gegensatz zum Ful-leren auch jeweils die Größe des Zentralteilchens ändert.Der starke Abfall der Anzahl der Moleküle im Bereich 1 lässt sich durch die oben er-wähnte, abweichende Definition des Bereiches 1 erklären. Aus diesem Grund wurdein Abbildung 4.13 die Zahl der Nachbarmoleküle ebenfalls für die bei den ande-ren Molekülen verwendete Definition berechnet und in Form der Quadratsymboledargestellt.

0

5

10

15

20

N

Bereich 1 Bereich 2

Gesamt


Abbildung 4.13: Anzahl der Moleküle in den beiden Abstandsbereichen derersten Solvathülle als Funktion der Anzahl der Substituenten(T=300 K). Die gestrichelte Linie dient nur der Illustration undsoll keinen physikalischen Zusammenhang darstellen.

34


4.3.2 Winkelverteilungen

Auch für die Solvensmoleküle lohnt sich ein Blick in die Orientierungsverteilungder Moleküle der ersten Solvatsphäre. Dazu wurden ebenfalls die Verteilungen derKosinusse der Winkel zwischen den analog zum Falle des C60 definierten Vektorenbetrachtet. Da es jedoch nicht möglich war, die zwei Bereiche so genau voneinanderzu trennen, wird besonders die Struktur des Bereiches 1 von Molekülen aus demBereich 2 überlagert. Die Verteilung von Benzen ist in allen Bereichen sehr symme-

-1 -0.5 0 0.5 1cosθ

0

0.5

1

P(c

os

θ)

0

0.2

0.4

0.6

0.8

P(c

os

θ)

Benzen

Toluen

-1 -0.5 0 0.5 1cosθ

0.2

0.3

0.4

0.5

0.6

P(c

os

θ)

0.3

0.4

0.5

0.6

0.7

P(c

os

θ)

Benzen

Toluen

(a) Bereich 2

Abbildung 4.14: Verteilung der Winkel in den beiden Abstandsbereichen für Benzenund Toluen. Im oberen Teil der Abbildung ist die Verteilung desinneren Vektors und im unteren Teil die Verteilung des Normalen-vektors dargestellt.

trisch. Im Bereich 1 bevorzugt der innere Vektor eine Ausrichtung senkrecht zumSchwerpunktverbindungsvektor, während der Normalenvektor vorzugsweise parallelorientiert ist.Im Bereich 2 ist der Normalenvektor zwischen ungefähr 60◦ und 120◦ nahezu gleich-verteilt. Der innere Vektor bevorzugt eine Ausrichtung von rund 30◦, 90◦ und 150◦.Es kann also eine ähnliche Struktur wie beim Fulleren beobachtet werden. Im erstenBereich stehen die Molekülebenen parallel zueinander, während im zweiten Bereicheine T-Struktur bevorzugt wird. Diese beiden Strukturen wurden bereits in der

35

4 Auswertung

Literatur diskutiert [14]. Weiterhin werden Strukturen diskutiert, bei denen die Mo-leküle parallelverschoben sind. Dies kann hier ebenfalls beobachtet werden, da auchdie Winkel zwischen den oben Genannten häufig auftreten.Auch beim Toluen gibt es Ähnlichkeiten zu den Verteilungen vom C60. Der innereVektor im Bereich 1 nimmt einen Winkel größer als 90◦ ein, während der Norma-lenvektor bevorzugt parallel zum Schwerpunktverbindungsvektor steht. Dies könnteein Hinweis auf eine verschobene parallele Anordnung sein, welche in der Literaturhäufig beschrieben wird [22].Im Bereich 2 gilt für den Normalenvektor das Gleiche wie beim Benzen. Die Maxi-ma der Verteilung für den inneren Vektor beschreiben jedoch eine Situation, bei derdie Methylgruppe bevorzugt auf den Schwerpunkt des anderen Moleküls oder genaudavon wegzeigt. Eine solche Struktur wird, neben anderen, auch von Schauer undBernstein [22] beschrieben.

-1 -0.5 0 0.5 1cosθ

0

0.5

1

1.5

P(c

os

θ)

0

0.2

0.4

0.6

0.8

P(c

os

θ)

m-Xylen

Mesitylen

(a) Bereich 1

-1 -0.5 0 0.5 1cosθ

0.3

0.35

0.4

0.45

0.5

0.55

P(c

os

θ)

0.4

0.45

0.5

0.55

0.6

P(c

os

θ)

m-Xylen

Mesitylen

(b) Bereich 2

Abbildung 4.15: Verteilung der Winkel in den zwei Bereichen für m-Xylen und Mesi-tylen. Im oberen Teil der Abbildung ist die Verteilung des innerenVektors und im unteren Teil die Verteilung des Normalenvektorsdargestellt.

Für m-Xylen und Mesitylen fälllt im Bereich 1 auf, dass der innere Vektor deutlichstärker zu Winkeln größer als 90◦ verschoben ist, während der Normalenvektoreine ähnliche Verteilung wie beim Toluen zeigt. Dies lässt vermuten, dass die oben

36


beschriebene, verschobene und parallele Anordnung für Mesitylen und m-Xylennoch stabiler ist und dadurch häufiger auftritt.Für den Bereich 2 lassen sich für m-Xylen die gleichen Aussagen wie für Toluentreffen. Es wird die Anordnung bevorzugt, bei der die Moleküle senkrecht zueinan-der stehen und eine Methylgruppe direkt auf das andere Molekül zeigt.Für Mesitylen lassen sich keine genauen Aussagen über die Struktur im Bereich 2treffen, da die Winkel nahezu gleichverteilt sind.

4.3.3 Abstandsabhängigkeit der Winkelverteilungen

Die Abstandsabhängigkeit der Winkelverteilungen um die Solvensmoleküle lässt sichanalog zum C60 betrachten (Abbildung 4.16). Auch hier wurde wieder über einenAbstandsbereich von 0.05 nm gemittelt.Der Funktionsverlauf stimmt in der ersten Solvatsphäre gut mit den Paarverteilungs-funktionen überein. Während beim C60 jedoch auch nach Ende der ersten Solvat-sphäre ein geringes Orientierungsverhalten gezeigt werden konnte, nähern sich dieFunktionswerte für die Orientierung der Solvensmoleküle schon nach kleinen Ab-ständen einem Wert von null an.Die Funktionen bestätigen die bei den Winkelverteilungen getroffenen Vermutun-gen. Ähnlich wie beim Fulleren ordnen sich die Moleküle um ein Solvensmolekülbei geringen Abständen zunächst parallel und dann senkrecht an. Auch hier wurdenahezu keine Temperaturabhängigkeit des Orientierungsverhaltens gefunden.

37

4 Auswertung

0 0.5 1 1.5 2r [nm]

-0.4

0

0.4

0.8



innerer Vektor

Normalenvektor

0

0.5

1

1.5

2

2.5g

(r)

260 K

280 K

300 K

320 K

340 K

360 K

(a) Benzen

0 0.5 1 1.5 2r [nm]

-0.4

0

0.4

0.8



innerer Vektor

Normalenvektor

0

0.5

1

1.5

2

g(r

)

260 K

280 K

300 K

320 K

340 K

360 K

(b) Toluen

0 0.5 1 1.5 2r [nm]

-0.4

0

0.4

0.8



innerer Vektor

Normalenvektor

0

0.5

1

1.5

g(r

)

260 K

280 K

300 K

320 K

340 K

360 K

(c) m-Xylen

0

0.5

1

1.5

g(r

) 260 K

280 K

300 K

320 K

340 K

360 K

0 0.5 1 1.5 2r [nm]

-0.4

0

0.4

0.8



innerer Vektor

Normalenvektor

(d) Mesitylen

Abbildung 4.16: Solvens-Solvens-Paarverteilungsfunktion und Abstandsabhängig-keit von 〈P2〉 bei 300 K

38

5 Zusammenfassung und Ausblick

Durch das Berechnen und Auswerten von Paarverteilungsfunktionen, Dichtenund Winkelverteilungen konnte in dieser Arbeit die erste Solvatsphäre des Buck-minsterfullerens C60 hinsichtlich ihrer Struktur in unterschiedlich substituiertenaromatischen Lösungsmitteln charakterisiert und verglichen werden.Bei den betrachteten Solvenzien lässt sich die erste Solvatsphäre um das Fullerenin zwei Abstandsbereiche aufteilen, in denen sich die Moleküle unterschiedlichorientieren. Im ersten Bereich, der je nach Solvensmolekül und Temperatur aus vierbis fünf Molekülen besteht, richten sich die Moleküle parallel zur Oberfläche desC60 aus, während im zweiten Bereich eine dazu senkrechte Anordnung bevorzugtwird. Die Größe der gesamten ersten Solvatsphäre ist stark solvens- und tempera-turabhängig. Sie besteht aus 15 (Mesitylen) bis 24 (Benzen) Molekülen, die in ihrerAusrichtung stark geordnet sind. In Abbildung 5.1 ist diese Ordnung gut zu sehen.Dort wurde mit Hilfe des Programmes VMD die gesamte erste Solvatsphäre desFullerens in Toluen dargestellt. Die vier roten Moleküle stellen den ersten Bereichda und die restlichen 16 den zweiten Bereich. Es ist gut zu erkennen, dass dieMoleküle des ersten Bereiches parallel zur Oberfläche ausgerichtet sind, währenddie hellblauen Moleküle meist senkrecht zur Oberfläche orientiert sind.Des Weiteren wurde versucht, die erste Solvatsphäre der Solvensmoleküle inähnlicher Weise zu beschreiben. Zwar ist eine Einteilung in zwei Bereiche anhandder Schwerpunkts-Schwerpunkts-Paarverteilungsfunktion nicht möglich, jedochkonnte auch hier grundsätzlich dasselbe Verhalten wie in der Solvatsphäre des C60

gefunden werden. Die Zahl der Moleküle in der ersten Solvatsphäre ist allerdingsweit weniger temperatur- und solvensabhängig. Sie liegt zwischen 13 und 15Molekülen. Es konnte gezeigt hier werden, dass sich die molekulardynamischeSimulation hervorrangend eignet, um einen Einblick in das Verhalten des Fullerensin Lösung zu erhalten.

Es gibt noch viel an diesem System zu untersuchen. Durch die Einführung von weite-ren Vektoren an den Solvensmolekülen könnte sowohl ihre Anordnung um das Fulle-ren als auch um andere Solvensmoleküle besser beschrieben und verstanden werden.Außerdem wäre eine Betrachtung, bei der zwei oder mehr C60-Moleküle in Lösungmiteinander wechselwirken können, sehr interessant. Die geplanten Experimente imArbeitskreis von Prof. Dr. Ludwig werden zeigen, inwieweit sich die verschiedenenStrukturen um das Fulleren durch IR-Spektroskopie nachweisen lassen.

39

5 Zusammenfassung und Ausblick

Abbildung 5.1: Gesamte erste Solvatsphäre von Toluen um das C60. Rot dargestelltsind dabei die vier Moleküle aus dem ersten Bereich.

40

Literatur

[1] Mark Abraham et al. GROMACS - Groningen Machine for Chemical Simu-lations User Manual Version 4.6.1.

[2] Edge A Vernon J Albright John G und Reginald Mills. “Relations betweenintradiffusion and viscosity coefficients. Tests of the albright equation”. In:Journal of the Chemical Society, Faraday Transactions 1: Physical Chemistryin Condensed Phases 79.6 (1983), S. 1327–1334.

[3] Mike P Allen, Dominic J Tildesley et al. “Computer simulation of liquids”. In:(1987).

[4] Jon Anderson, John J Ullo und Sidney Yip. “Molecular dynamics simulationof hydrocarbon molecules in condensed phases. II. Benzene”. In: The Journalof chemical physics 86.7 (1987), S. 4078–4089.

[5] Tarek Baati et al. “The prolongation of the lifespan of rats by repeated oraladministration of [60] fullerene”. In: Biomaterials 33.19 (2012), S. 4936–4946.

[6] M Isabel Cabaço et al. “Neutron diffraction and molecular dynamics study ofliquid benzene and its fluorinated derivatives as a function of temperature”.In: The Journal of Physical Chemistry B 101.35 (1997), S. 6977–6987.

[7] Jan Cami et al. “Detection of C60 and C70 in a young planetary nebula”. In:Science 329.5996 (2010), S. 1180–1182.

[8] Jens Fischer und Alarich Weiss. “Transport Properties of Liquids. V. Self Dif-fusion, Viscosity, and Mass Density of Ellipsoidal Shaped Molecules in the PureLiquid Phase”. In: Berichte der Bunsengesellschaft für physikalische Chemie90.10 (1986), S. 896–905.

[9] Daan Frenkel und Berend Smit. Understanding molecular simulation: fromalgorithms to applications. Bd. 1. Academic press, 2001. url: http://www.sciencedirect.com/science/book/9780122673511.

[10] Berk Hess et al. “LINCS: a linear constraint solver for molecular simulations”.In: Journal of computational chemistry 18.12 (1997), S. 1463–1472.

[11] J.K.Lehmann und E. Bich. Physikalische Chemie - Hauptpraktikum (Univer-sität Rostock). 2014.

41

Literatur

[12] William L Jorgensen, David S Maxwell und Julian Tirado-Rives. “Develop-ment and testing of the OPLS all-atom force field on conformational energe-tics and properties of organic liquids”. In: Journal of the American ChemicalSociety 118.45 (1996), S. 11225–11236.

[13] William L Jorgensen und Nora A McDonald. “Development of an all-atomforce field for heterocycles. Properties of liquid pyridine and diazenes”. In:Journal of Molecular Structure: THEOCHEM 424.1 (1998), S. 145–155.

[14] William L Jorgensen und Daniel L Severance. “Aromatic-aromatic interac-tions: free energy profiles for the benzene dimer in water, chloroform, andliquid benzene”. In: Journal of the American Chemical Society 112.12 (1990),S. 4768–4774.

[15] Harold W Kroto et al. “C 60: buckminsterfullerene”. In: Nature 318.6042(1985), S. 162–163.

[16] Nora A McDonald und William L Jorgensen. “Development of an all-atomforce field for heterocycles. Properties of liquid pyrrole, furan, diazoles, andoxazoles”. In: The Journal of Physical Chemistry B 102.41 (1998), S. 8049–8059.

[17] Shuichi Miyamoto und Peter A Kollman. “SETTLE: an analytical version ofthe SHAKE and RATTLE algorithm for rigid water models”. In: Journal ofcomputational chemistry 13.8 (1992), S. 952–962.

[18] Dietmar Paschek und Alfons Geiger. “MOSCITO 4”. In: Department of Phy-sical Chemistry, University of Dortmund (2002).

[19] Dietmar Helmut Paschek. Molekulardynamik Simulation der hydrophoben Hy-dratation nichtionischer Tenside. Mainz, 1999.

[20] Jean-Paul Ryckaert, Giovanni Ciccotti und Herman JC Berendsen. “Nume-rical integration of the cartesian equations of motion of a system with cons-traints: molecular dynamics of n-alkanes”. In: Journal of Computa-tional Physics 23.3 (1977), S. 327–341.

[21] FM Samigullin. “Study of translational self-diffusion of molecules in liquids”.In: Journal of Structural Chemistry 14.4 (1974), S. 569–574.

[22] Mark Schauer und ER Bernstein. “Calculations of the geometry and bindingenergy of aromatic dimers: benzene, toluene, and toluene–benzene”. In: TheJournal of chemical physics 82.8 (1985), S. 3722–3727.

42

Molekulardynamische Simulation des Buckminster...

Documents

Transcript of Molekulardynamische Simulation des Buckminster...