15. Vorlesung SS 2006 Computational Chemistry V5: Klassische statistische Mechanik –...

1 5. Vorlesung SS 2006

Computational Chemistry

V5: Klassische statistische Mechanik – thermodynamische Ensembles

Die statistische Mechanik behandelt Systeme mit vielen (im Grunde unendlich

vielen) Freiheitsgraden. Diese sollen durch wenige Makrovariablen beschrieben

werden. Eine besondere Rolle spielt hierbei die Makrovariable Entropie.

Die Freiheitsgrade sind für Biomoleküle speziell die Positionen q und Impulse

p = m v der Atome (mit deren Masse m und ihren Geschwindigkeiten v).

In klassischen Systemen (d.h. ohne Berücksichtigung von Quanteneffekten)

kommutieren alle Freiheitsgrade miteinander. [q,p] = 0.

Den Raum, der durch die Freiheitsgrade aufgespannt wird, nennt man den

Phasenraum. Jeder Punkt im Phasenraum bezeichnet eine bestimmte

Konfiguration des Systems.

Im zeitlichen Verlauf beschreibt das System einen Pfad im Phasenraum, der durch

die Bewegungsgleichungen bestimmt wird.



Wiederholung: Boltzmann-Verteilung(1) In einem System mit N Teilchen sei die Teilchenzahl konstant.

(2) Gesamtenergie des Systems sei konstant.

D.h. es gibt Energieaustausch zwischen den Teilchen, aber nicht mit der Umgebung.

(3) Wenn solch ein System im Gleichgewicht ist, ist die Energie der Teilchen E

entsprechend einer Boltzmann-Verteilung populiert:

00 exp

E

ENEN

Boltzmann-verteilte Systeme findet man in vielen Bereichen der Physikalischen Chemie.

Welche Systeme besitzen eine Boltzmann-Verteilung?



Wiederholung: PhasenraumdichteDie Wahrscheinlichkeitsdichte im Phasenraum (= kurz die “Phasenraumdichte”) ist

im kanonischen Ensemble proportional zum Boltzmann-Faktor:

XEXNVT exp

wobei E die Gesamtenergie des Systems ist und = kBT.

Für zwei Zustände des Systems X und X’ lautet das Verhältnis ihrer Wahrscheinlich-

keiten:

'

,exp'

XEXEE

EX

X

NVT

NVT

Was ist daskanonischeEnsembleeigentlich?



Wiederholung: Phasenraumdichte

Der Normalisierungsfaktor der ersten Gleichung ist die Zustandssumme des

gesamten Phasenraums (Raum der 3N Koordinaten und 3N Geschwindigkeiten):

dxxENh

Q

Q

XE

NhX

NNVT

NVTNNVT

exp!

1

exp

!

1

3

3

Der Erwartungswert einer Observablen A des Systems lässt sich darstellen als:

dxxAxxA NVTNVT

Im Metropolis-Algorithmus erzeugt man eine geeignete Markov-Kette von

Konfigurationen, so dass der Erwartungswert von A als einfacher Mittelwert folgt:

M

ii

MNVT

xAM

xA1

1lim



Wiederholung: Markov-Kette

Betrachte Markov-Kette von N molekularen Zuständen {X1, X2, X3, ...} mit einer

Verteilung NVT(X) für N .

In einer Markov-Kette gehört jeder Zustand zu einer endlichen Menge an

Zuständen aus dem Zustandsraum D0 D.

Für die konditionelle Verteilung jedes Zustands P bezüglich aller vorherigen

Zustände gilt: nnnn XDXPXXDXP 01001 ,...,

d.h. das Ergebnis Xn+1 hängt nur von Xn ab.

Der Metropolis-Algorithmus erzeugt eine stochastische und ergodische

Übergangsmatrix für die Markovkette, so dass die Verteilung für jeden Zustand Xi

im Limit i = NVT (Xi) ist.

So wird eine Phasenraumtrajektorie im kanonischen Ensemble erzeugt.



Wiederholung: mikroskopische Reversibilität (“detailed balance”)

Lege Übergangsmatrix fest durch Definition einer Übergangswahrscheinlichkeit für

jeden Übergang von Xi nach Xj, so dass mikroskopische Umkehrbarkeit erfüllt ist:

jijiji kk

Das Verhältnis der Übergangswahrscheinlichkeiten hängt damit nur vom

Energieunterschied zwischen den Zuständen i und j ab:

jiij

ijij

ji

j

i

XEXEE

Ek

k

exp

Was für ein Ensemble bekommen wir mit diesen Übergangswahrscheinlichkeiten?



Wiederholung: Metropolis AlgorithmusDie am häufigsten verwendete Technik zur Auswahl von Konformeren („importance

sampling“) mittels Monte-Carlo-Methoden ist der Metropolis Algorithmus:

(1) konstruiere Anfangskonfiguration des Moleküls

(2) führe zufällige Änderung eines Freiheitsgrades (z.B. eines Torsionswinkel)

durch.

(3) berechne Änderung der Energie E aufgrund dieser Änderung der

Konformation.

(4) falls E < 0 akzeptiere die neue Konfiguration

falls E > 0 berechne die Wahrscheinlichkeit

erzeuge Zufallszahl r im Intervall [0,1]

akzeptiere die neue Konfiguration, falls w r, sonst verwerfe sie.

Da die Boltzmann-gewichtete Energiedifferenz mit einer Zufallszahl verglichen

wird, werden auch vereinzelt Konformere hoher Energie akzeptiert. Daher erhält

man ein Ensemble (Menge) von Konformationen mit einer Energieverteilung

entsprechend einer Boltzmann-Verteilung.

Tk

Ew

B

exp



abgeschlossene und nicht abgeschlossene Systeme

Wenn ein System thermisch und mechanisch abgeschlossen, also von seiner

Umgebung isoliert, ist, dann bleibt seine interne Energie zeitlich konstant.

Warum?

Wenn ein System nicht abgeschlossen ist, wird es mit der Zeit die Temperatur

seiner Umgebung annehmen (die wir als konstant annehmen).

Physikalische Grössen, die entweder konstant sind oder deren Durchschnittswert

durch die Umgebung bestimmt wird, bezeichnet man als Systemparameter.

Unterschiedliche experimentelle Umstände sorgen dafür, dass unterschiedliche

Parameter konstant gehalten werden.

In der Theorie der statistischen Physik entsprechen diese Fälle unterschiedlichen

Ensembles.



Ensembles Simulationen

Ein Ensemble ist eine gedankliche Konstruktion: wir denken uns eine sehr große Anzahl M

physikalisch identischer Kopien des isolierten Originalsystems, die sich in beliebigen

Bewegungszuständen befinden können, dargestellt durch Ensemblepunkte (q,p) im

Phasenraum. Jeder Ensemblepunkt charakterisiert den Bewegungszustand eines

Ensemblemitglieds. Wir wollen annehmen, daß wir das Ensemble durch eine

Ensembledichte (q,p,t) beschreiben können: (q,p,t)d soll der Bruchteil der M

Phasenraumpunkte sein, der sich zur Zeit t im Phasenraumelement d = dp dq bei (q,p)

aufhält. Damit muss auch die Ensembledichte auf den Wert 1 normiert sein:

Diese Normierung erlaubt auch die Interpretation, daß (q,p,t)d die Wahrscheinlichkeit ist,

den Ensemblepunkt eines beliebig herausgegriffenen Ensemblemitgliedes im Element d

bei (q,p) zu finden.

Zur Simulation von klassischen Vielteilchensystemen gibt es im Wesentlichen zwei

Methoden: Monte Carlo und Moleküldynamik.

Um diese Simulationstechniken an die experimentellen Situationen anzupassen ist in jedem

Fall eine sorgfältige Betrachtung notwendig, welches Ensemble dazu geeignet ist, d.h. in

welchem Ensemble man simulieren muss.

1),,( tpqd



Ergoden-Hypothese.

Das grundlegende Postulat der statistischen Mechanik gilt für Systeme mit

konstanter Energie E, Volumen V und Teilchenzahl N.

Das Postulat besagt, dass alle für das System zugänglichen Zustände (die eine

bestimmte Energie, Volumen und Teilchenzahl besitzen) im Verlauf der Zeit mit

gleicher Häufigkeit = Wahrscheinlichkeit eingenommen werden.

Dies ist die Ergoden-Hypothese.

Daher kann der zeitliche Mittelwert einer physikalischen Messgrösse A (also

einer „Observablen“) durch den Mittelwert über alle zugänglichen Zustände ersetzt

werden A.

Mit den Zuständen X des Systems gilt:

A



Ensemble-Mittelwerte

Hier ist H(X) der Hamiltonian, der die Energie des Systems an einem Punkt X im

Phasenraum angibt.

Damit gilt H = T + V, also die Summe aus den kinetischen und potentiellen

Energien T und V.

Der Nenner bewirkt die entsprechende Normierung. Die Summe

bezeichnet eine Summe über alle Zustände X mit einer bestimmten

Energie E. In der rechten Summe wird diese Einschränkung durch die

Delta-Funktion erreicht. (Zusätzlich ist das System auf ein bestimmtes Volumen

und auf eine bestimmte Teilchenzahl eingeschränkt.)

Der Mittelwert heisst der Ensemble-Mittelwert.

Das Ensemble NVE = konstant heisst das mikrokanonische Ensemble.

A

EXH

EXHXAXAA

X

X

EX

EX

EX



Entropie

Die Anzahl an Zuständen mit Energie E ist

X

EXHEVN ,,

Die Entropie wird mittels (N, V, E) und der Boltzmann-Konstante kB definiert als:

EVNkEVNS B ,,ln,,

Die thermodynamischen Grössen Temperatur T, chemisches Potential und

Druck P ergeben sich als Ableitungen der Entropie nach den Systemparametern:

NEVEVN V

STp

N

ST

E

ST

,,

1

,

,,

Dies folgt aus dem ersten Gesetz der Thermodynamik:

dNpdVTdSdE



Phasenraum

Ein System aus N Teilchen besitzt 3N Koordinaten und 3N Geschwindigkeiten.

Daher besetzt es zu jedem Moment einen Punkt X im 6N-dimensionalen

Phasenraum. Dessen Koordinaten sind:

X = (x1, y1, z1, px,1,py,1, pz,1, x2, y2, z2, px,2, py,2, pz,2, ... , xN, yN, zN, px,N,py,N, pz,N )

Dabei kann man separat die Orts- und Impulskoordinaten definieren als

q = (x1, y1, z1, x2, y2, z2, ..., xN, yN, zN )

p = (px,1,py,1, pz,1, px,2, py,2, pz,2, ..., px,N, py,N, pz,N )

und damit diesen Punkt im Phasenraum als X = (q,p) darstellen.



Beispiel: der harmonische Oszillator

Der Phasenraum dieses harmoni-

schen Oszillators hat nur 2

Dimensionen, Position und Impuls.

Wir definieren den Ursprung dieses

Phasenraums als den Punkt, an dem

der Ball der Masse m in Ruhe ist

(d.h. Impuls = 0) wobei die Feder in

ihrer Ruhelage ist.

Dieser Punkt entspricht einem

stationären Zustand des Systems.

Nun verfolgen wir das dynamische

Verhalten des Systems falls wir von

einem anderen Punkt als dem

Ursprung starten.

Die Gesamtenergie = potentielle

Energie + kinetische Energie ist

konstant.

Aus: C. Cramer, Computational Chemistry, Wiley




Betrachten wir nun die Phasenraumtrajektorie des Oszillators.

Zu jedem beliebigen Zeitpunkt gilt: 2

1

12

t

t

dtm

tptqtq

wobei wir die Beziehung zwischen Geschwindigkeit und Impuls verwendet haben:

m

pv

Zwei Impulsvektoren zu verschiedenen

Zeitpunkten sind analog gemäß 2

1

12

t

t

dttamtptp

verknüpft, wobei a die Beschleunigung ist. Diese beiden Gleichungen sind die

Newtonschen Bewegungsgleichungen. Weiterhin benutzen wir das zweite

Newtonsche Gesetz:

m

Fa

mit der Kraft F und den Zusammenhang zwischen Kraft und der Ableitung der

potentiellen Energie nach der Ortskoordinate:

r

EF

(1)

(2)




Für kompliziertere Systeme ist es fast nie möglich, eine analytische Lösung der

Bewegungsgleichungen zu finden.

Wenn wir die Gleichungen (1) und (2) durch

t

m

tttt

pqq

und

tmttt tapp

annähern, können wir die exakte Lösung der Bewegungsgleichung numerisch

annähern.

Dann erhalten wir als analytische Lösung der Bewegungsgleichungen:

t

m

kmkbtpt

m

kbtq sin,cos

Sie erzeugen eine Bewegung auf der elliptischen Bahn im Phasenraum.




Aus: C. Cramer, Computational Chemistry, Wiley



Das kanonische NVT-Ensemble

In experimentellen Situationen muss oft die Temperatur konstant gehalten

werden und nicht die Energie E (da man das System gar nicht perfekt isolieren

kann bzw. möchte).

Dies stellt man sich durch Kopplung des Systems an ein externes Wärmebad vor,

das viel grösser als das System ist (und dessen Temperatur daher konstant ist).

Dadurch ergibt sich, dass ein zeitlicher Mittelwert des Systems gleich einem

gewichteten Mittelwert über alle Zustände mit festem Volumen und Teilchenzahl

ist. Den Gewichtungsfaktor bezeichnet man als Boltzmann-Faktor

XHTk

XH

B

expexp



Das kanonische NVT-Ensemble

Damit gilt nun für den Ensemble-Mittelwert einer Observablen A

X

XNVT

XHTVNZ

XHXAZ

A

exp,,

exp1

Der Faktor Z bewirkt die korrekte Normierung. Er heisst Zustandssumme.

Z hängt mit der freien Energie F (free enthalpy) zusammen:

TVNZTkF B ,,ln

Bezüglich der thermodynamischen Grössen gilt auch:

TSEF

Im Gleichgewicht nimmt die freie Energie ihr Minimum unter der Bedingung festes

Volumen und feste Teilchenzahl an.

Der Mittelwert ANVT heisst der kanonische Ensemble-Mittelwert.



Das kanonische Ensemble

Aus diesen Betrachtungen ergibt sich sofort, dass die Newtonsche

Bewegungsgleichungen keine geeignete Methode sind um die Zustände in diesem

Ensemble zu samplen.

Warum?

Weil die Newtonschen Bewegungsgleichungen ja die Energie erhalten.

Jede Trajektorie generiert die Dynamik auf einer Hyperfläche des Phasenraums

bei konstanter Energie.

Damit dies einem System mit Boltzmann-Verteilung entspricht, benötigt man

Übergänge zwischen Oberflächen konstanter Energie, die jeweils mit dem

Boltzmann-Faktor gewichtet werden sollten.



Beziehung zwischen den Ensembles

Die Zustandssumme kann als Summe über die Zustände fester Energie

dargestellt werden: E

EVNETVNZ ,,exp,,

Die Anzahl der Zustände (N,V,E) nimmt sehr rasch mit E zu, wogegen die

Boltzmann-Verteilung sehr rasch mit E abnimmt.

Das Produkt der beiden Funktionen hat daher ein scharfes Maximum bei einem

Wert und das System wird meist eine Energie sehr eng bei diesem Wert haben.

Daher besteht in der Praxis meist kein grosser Unterschied zwischen dem

kanonischen und dem mikrokanonischen System.

Aufgrund des Gesetzes der grossen Zahlen, haben messbare physikalische

Grössen sehr kleine Fluktuationen.

Die Abweichungen zwischen den Ensembles werden jedoch grösser, je kleiner die

Systeme werden Vorsicht also bei Simulationssystemen!

E




Mit der Definition der Entropie

kann man

umschreiben in:

wobei FE die Freie Energie E – TS ist, wobei S im mikrokanonischen Ensemble für

die Energie berechnet wurde.

Man sieht, dass die Summe tatsächlich durch die Zustände dominiert wird, für die

die freie Enthalpie minimal ist.

EVNkEVNS B ,,ln,,

E

EVNETVNZ ,,exp,,

E E

EFTSETVNZ expexp,,



Vollständiges Differential

Definition: Sei f(x) eine Funktion einer Veränderlichen, dann ist

dxxfdf das totale oder vollständige Differential.

Das Differential beschreibt die Änderung des Funktionswerts

auf der in (x, f(x)) errichteten Tangente an die Kurve von f(x).

Definition: Sei f(x) eine Funktion von zwei

Veränderlichen, dann ist

dyy

fdx

x

fdf

das vollständige Differential.




Der erste Satz der Thermodynamik

dNpdVTdSdE

lässt sich mit dem Differential

NVTNTV T

FS

V

Fp

N

F

,,,

,,

SdTTdSdEdF

TSEF

Damit lassen sich aus der freien Energie F folgende thermodynamische Grössen

ableiten:

umformen in: dNpdVSdTdF



Das isotherm-isobare NpT-Ensemble

In experimentellen Situationen muss zusätzlich zur Temperatur oft ebenfalls der

Druck als konstant angenommen werden. Dann ist aber das Volumen V nicht

mehr konstant.

Man erhält für die Mittelwerte von Observablen im isotherm-isobaren NpT-

Ensemble:

TVNZpVdVXHpVdVTpNQ

XAXHpVdVQ

A

X

XNPT

,,expexpexp,,

expexp1

Q heisst wieder die Zustandssumme.

Q verhält sich zur kanonischen Zustandssumme Z ähnlich wie Z mit der Funktion

im mikrokanonischen Ensemble zusammenhing.



Das isotherm-isobare NpT-Ensemble

Q hängt mit der Gibbs‘schen Freien Enthalpie G (englisch: Gibbs free energy)

folgendermassen zusammen: TpNQTkG B ,,ln

G kann wieder mit den thermodynamischen Grössen ausgedrückt werden:

pVTSEG Analog zu oben können folgende thermodynamische Grössen aus G(N,P,T)

abgeleitet werden:

NpTNTp T

GS

p

GV

N

G

,,,

,,



Das grosskanonische VT-Ensemble

Wenn nun wieder das Volumen konstant ist, aber die Zahl der Teilchen variiert

werden kann, erhalten wir den gross-kanonischen Mittelwert:

TVNZNXHNTVZ

XAXHNZ

A

X NNG

XNGVT

,,expexpexp,,

expexp1

ist das chemische Potential für die Addition oder Entfernung eines Teilchens,

d.h.

1

NNN

G

ZG sollte nicht mit der kanonischen Zustandssumme Z(N,V,T) verwechselt werden.

(Ihr Zusammenhang ist oben dargestellt.)

ZG definiert das grosskanonische Potential G:

FNTVG ,,



Das grosskanonische VT-Ensemble

Die interne Energie kann mit den Variablen S, V und N dargestellt werden

(Gibbs-Duhem Gleichung)

NpVTSNVSE ,,

Analog zu oben können folgende thermodynamische Grössen aus dem

grosskanonischen Potential abgeleitet werden:

,,,

,,V

G

T

G

TV

G

TS

VpN

so dass gilt: pVTVG ,,



Zusammenfassung: gebräuchliche statistische Ensembles

Das mikrokanonische Ensemble NVE = konstant

Das kanonische Ensemble NVT = konstant

Das isotherme-isobare Ensemble NPT = konstant

Das grosskanonische Ensemble VT = konstant



Anwendungen

Das mikrokanonische Ensemble NVE = konstant

überprüfe Stabilität von Integrationsalgorithmen,

Erzeugung von „wahrer“ Dynamik

Das kanonische Ensemble NVT = konstant

simuliere Vorgänge unter Druckschwankungen

Das isotherme-isobare Ensemble NPT = konstant

übliches Ensemble für Simulationen von Biomolekülen

Das grosskanonische Ensemble VT = konstant

Stuchebrukhov-Paper über Zahl der Wassermoleküle in Cytochrom c Oxidase



Anwendung: identifiziere Wassermoleküle in Proteinen

Kristallstrukturen von Proteinen zeigen

oft nicht die Position interner Wasser-

moleküle, da diese zu mobil sind.

Verwende Computersimulationen um

die „wahre“ Hydratation zu finden.

Welches ist das beste Ensemble um

Wassermoleküle in ein Protein zu

positionieren?

Das gross-kanonische Ensemble!

Tashiro, Stuchebrukhov, J Phys Chem B 109, 1015 (2005)



Anwendung: identifiziere Wassermoleküle in Proteinen

Betrachte Transfer aus dem Lösungsmittel

Wasser ins Proteininnere.

In Lösung besitzt ein Wassermolekül eine

freie Lösungsenthalpie von -12.3 kcal/mol.

Wenn es im Proteininneren eine günstigere

Position findet, dann wird diese im zeitlichen

Mittel mit einem Wasser besetzt sein.

Führe Simulation im semi-grosskanonischen

Ensemble durch, wobei während der

Simulation Wassermoleüle in das Protein

hinzugefügt bzw. daraus entfernt werden.

In diesem Fall ist es gerade so günstig, die

Kavität mit 2 Wassermolekülen zu füllen.

Tashiro, Stuchebrukhov, J Phys Chem B 109, 1015 (2005)

15. Vorlesung SS 2006 Computational Chemistry V5: Klassische statistische Mechanik –...

Documents

Transcript of 15. Vorlesung SS 2006 Computational Chemistry V5: Klassische statistische Mechanik –...