Beschreibung dynamischer molekularer Systeme mit effizienten linear-skalierenden QM- und QM/MM-Methoden

der Fakultät für Chemie und Pharmazie

der Ludwig-Maximilians-Universität Mün hen

Bes hreibung dynamis her molekularer Systeme

mit ezienten linear-skalierenden

QM- und QM/MM-Methoden

Denis Benjamin Flaig

aus

Stuttgart

DieseDissertation wurde im Sinne von Ÿ 7 der Promotionsordnung vom 28. November 2011

vonHerrn Prof.Dr. C. O hsenfeld betreut.

Eidesstattli he Versi herung

DieseDissertationwurde eigenständig und ohneunerlaubte Hilfe erarbeitet.

Mün hen, 24.09.2012

(Denis Flaig)

Dissertation eingerei ht am24.09.2012

1.Guta hter: Prof. Dr. C.O hsenfeld

2.Guta hter: Prof. Dr. H.Ebert

Herrn Prof. O hsenfeld danke i h herzli h für die Mögli hkeit, in seiner Arbeitsgruppe die

Dissertation anfertigenzu können. VielenDankfür dieBereitstellung desFors hungsthemas,

die freundli he Arbeitsatmosphäreund beständige Unterstützung.

Prof. Ebertdanke i h herzli h fürdie AnfertigungdesZweitguta htens.

Auÿerdem bedanke i h mi h beiden Kolleginnen und Kollegen inder Arbeitsgruppe für die

nette, aufbauende Zusammenarbeit in den letzten Jahren und insbesondere den vielfältigen

HilfenvonallenSeiten.Danke,Dr.M.Beer,Dr.T.Berezniak,I.Blank,Dr.L.Clin,Dr.B.

Do-ser, Dr. M. Hanni, T. Kessel, J. Kluge, Dr. J. Kussmann,Dr. D.Lambre ht, Dr. M. Löer,

A.Lünser,Dr.B.Maryasin,M.Maurer,S.Maurer,S.Roÿba h,Dr.K.Sadeghian,J.S häer,

1 Einleitung 7

2 Theoretis he Grundlagen 11

2.1 S hrödinger-Glei hung undKern-Elektron-Separation . . . 11

2.2 Molekulardynamis he Methoden . . . 14

2.3 Quanten hemis he Elektronenstruktur-Methoden . . . 16

2.4 Bere hnung statis her Moleküleigens haften . . . 24

2.5 Molekularme hanis he Methoden . . . 30

2.6 QM/MM-Methoden . . . 32

3 Eziente linear-skalierende SCF-Verfahren 37 3.1 Einführung . . . 37

3.2 Linear-skalierendes quasi-indirektes bzw. semidirektes SCF-Verfahren . . . 38

3.3 Abstands-abhängiges S hwarz-S reening . . . 51

3.4 Zusammenfassung undAusbli k . . . 53

4 QM-Ansätze zur ezienten Bes hreibung der Molekulardynamik 55 4.1 Einführung . . . 55

4.2 Extrapolationsmethode für CPSCF-Anfangss hätzungen . . . 57

4.3 Integralwiederverwendung fürveränderte Konformationen . . . 62

4.4 QM/MM-Implementierung für Molekulardynamik . . . 70

4.5 Zusammenfassung undAusbli k . . . 71

5 QM- und QM/MM-Bere hnung von NMR-Vers hiebungen 73 5.1 Einführung . . . 73

5.3 IntermediäreReferenzen beiderBere hnung von NMR-Vers hiebungen . . . 82

5.4 Konvergenz von QM/MM-Abs hirmungen mit derGröÿe derQM-Region . . . . 88

5.5 Zusammenfassung undAusbli k . . . 105

6 Enzymatis he Reparatur oxidativer DNA-S häden 109 6.1 Einführung . . . 109

6.2 Einzelheiten zurStrukturauswahl und-vorbereitung . . . 111

6.3 Gesamte S hadensstabilisierung . . . 112

6.4 Individuelle Beiträge derEnzymresiduen zurS hadensstabilisierung . . . 116

6.5 Zusammenfassung undAusbli k . . . 121

7 Allgemeine Zusammenfassung 123

Anhang A Q-Chem-ChemShell-Interfa e 127

Anhang B Molekularer Testsatz 129

Anhang C Details zur QM-Gröÿenkonvergenz von NMR-Abs hirmungen 131

Anhang D Details zur Untersu hung der DNA-Reparatur 143

Tabellen- und Abbildungsverzei hnis 145

Literaturverzei hnis 151

Einleitung

DieTheoretis he Chemie kann heute auf eineVielzahl anMethoden zurü kgreifen, um

ther-modynamis he und kinetis he Eigens haften oder die Dynamik molekularer Systeme zu

be-s hreiben. Die Bandbreitetheoretis her Methoden erstre kt si hbeispielsweise von einfa hen

molekularme hanis hen Verfahren über semi-empiris he bis hin zu quanten hemis hen

Me-thoden (füreine Übersi ht siehe z.B.Referenz [1℄). Die vers hiedenen methodis hen Ansätze gründen dabei auf sehr unters hiedli h genauen physikalis hen Bes hreibungen der Materie.

MitsteigenderGenauigkeitderzugrundeliegendenBes hreibungerrei hendieMethoden

typi-s herweise einehöhereVorhersagekraft,allerdingszumeist aufKosteneines maÿgebli h

höhe-renRe henaufwands.UmmolekulareSystemederaktuellen hemis henFors hungsehrgenau

theoretis hzu bes hreiben, wirdderRe henaufwand trotz dererhebli hen methodis hen und

te hnologis henFortentwi klungen derletztenJahrzehntes hnellzumeigentli h

bes hränken-den Faktor.Mitdenimmens gestiegenenMögli hkeiten hat si hnämli hglei hzeitigau hder

Anspru h an dieTheoretis he Chemie massiv erweitert: Galt dieAufmerksamkeit früher oft

besonderseinzelnenund eherkleinenMolekülen,hat si hdertheoretis he Bli kpunkt

spätes-tensmit denForts hritten dersupramolekularen Chemieinden 80erund 90erJahrendarauf

verlagert,weitausgröÿere molekulareEinheiten zuuntersu hen.Heute strebtdieTheoriean,

sowohlsynthetis hesupra-undmakromolekularealsau hgroÿebiomolekulareSysteme

zuver-lässig zubes hreiben.Dabeispieltzuglei h einegroÿeRolle,weitrei hende Umgebungseekte

(z.B. Lösungsmitteleekte) mögli hst genau zuberü ksi htigen.

Eineents heidendeKenngröÿedertheoretis henMethodenistindiesemHinbli kdasSk

a-lenverhaltendesRe henaufwandsmitderSystemgröÿe

M

.Sosteigtbeispielsweisefür

Molekülgröÿean, sondern mit höheren Potenzen von

M

.Umdas Skalenverhalten wesentli h

zureduzieren,konntenallerdingsindenletztenJahrzehntenwi htigemethodis heForts hritte

erzieltwerden.Diemethodis henFortentwi klungensenkenbeispielsweisedasSkalenverhalten

für viele quanten hemis he Re hnungen auf Basis der Hartree-Fo k- oder

Di htefunktional-Theorie mit zunehmender Systemgröÿe auf linear und au h für Post-Hartree-Fo k-Theorien

konnteninzwis hen linear-skalierende Methoden entwi kelt werden(füreinen Überbli ksiehe

z.B.Referenzen [2,3℄).

NebendemSkalenverhaltenistfürdenRe henaufwanddertheoretis henMethodenzudem

sehrbedeutsam,wiegroÿdergegenüberdemAnstiegmitderSystemgröÿekonstanteVorfaktor

derRe henzeit ist. Sobegründet beispielsweise derkleinereVorfaktorder

molekularme hani-s hen Methoden deren typis hen Ges hwindigkeitsvorteil au h gegenüber linear-skalierenden

quanten hemis hen Methoden allerdings umden Preis einer Bes hreibung auf niedrigerem

theoretis hen Niveau und daher einer geringeren Anwendbarkeit und Aussagekraft. Um die

Vorzüge derlinear-skalierenden quanten hemis hen Methoden optimal ausnutzen zu können,

besteht ein wi htiges Ziel darin, den Vorfaktor für linear-skalierende QM-Methoden ohne

si-gnikante Genauigkeitseinbuÿen weiter zu senken.

Diesem Ziel entspre hend bes hreibt ein wesentli her Teil derDissertationW

eiterentwi k-lungen, die zu einer Reduktion des Vorfaktors für linear-skalierende quanten hemis he

Me-thoden auf Basis der Hartree-Fo k- oder Kohn-Sham-Di htefunktional-Theorie führen. Das

Kapitel 3 geht dabei auf methodis he Fortentwi klungen zur Reduktion des Vorfaktors für

eine gegebene Molekülgeometrie ein, während das Kapitel 4 ergänzende Ansätze für

Re h-nungen mehrerer Molekülgeometrien im Rahmen der Born-Oppenheimer-Molekulardynamik

umfasst.

So stellt dasKapitel 3 ein linear-skalierendes Integralspei her-Verfahren dar,das die

Vor-teiledertraditionell indirektenodersemidirektenAnsätzemitdenen direkter SCF-Methoden

kombiniert. Darüber hinaus geht das Kapitel auf die wesentli hen Beiträge für ein neues

Abstands-abhängiges S hwarz-S reening ein, das zu einer weiteren Ezienzsteigerung von

SCF-Methoden führt.

DasKapitel4bes hreibtzweimethodis heAnsätze,umbeidynamis henmolekularen

Sys-temendieBorn-Oppenheimer-Hyperä heezientabzutasten:einenExtrapolationsansatzfür

Wiederverwen-Molekulardynamik. Auÿerdem werden die E kpunkte einer neu implementierten

Programm-s hnittstelle fürgemis hte QM/MM-Energie-und Gradientenre hnungen dargestellt.

Das Kapitel 5 widmet si h ans hlieÿend der Kombination und Untersu hung bestehender

MethodenzurgenauenBere hnungkernmagnetis herVers hiebungen.Einbesonderer

S hwer-punktliegthierbeidarauf,QM/MM-Ansätze[4℄mitlinear-skalierender NMR-Methoden[57℄ zu kombinieren und zu untersu hen, wie groÿ dieQM-Regionen gewählt werden müssen,um

NMR-Vers hiebungenvonzentralenMolekülregionenzuverlässigzubere hnen.Auÿerdem

um-fasst das Kapitel eine breite Untersu hung der Methodengenauigkeit von GIAO-SCF-, und

-MP2 -MethodenfürNMR-Vers hiebungsre hnungen.ZudemwirdeinweiterführenderAnsatz

zurSteigerungderMethodengenauigkeiten amBeispielvon Polypeptiden untersu ht.

Im letzten Teil der Arbeit (Kapitel 6) werden s hlieÿli h Ergebnisse einer

quanten hemi-s hen Untersu hung zur enzymatis hen Reparatur oxidativer DNA-S häden dargestellt. Die

Untersu hungquantiziert dieunters hiedli he Stabilisierung ges hädigter DNAimVerglei h

zu intakter DNA dur h das Reparaturenzym und identiziert wi htige Proteinresiduen bei

der S hadensstabilisierung. Ein besonderes Augenmerk der Untersu hung liegt ebenfalls auf

Aspektendertheoretis hen Bes hreibung,wie Basissatz-undMethodenfehler sowie der

Theoretis he Grundlagen

Das folgende Kapitel umfasst die theoretis hen Grundlagen der quantenme hanis hen (QM)

und molekularme hanis hen (MM) Methoden, die im Rahmen der Dissertation angewendet,

untersu ht undweiterentwi kelt werden.DieDarstellung beginnt ausgehendvonden

Postula-ten der(ni ht-relativistis hen) Quantenme hanik und entwi kelt vondiesem Standpunkt aus

die grundlegendeTheorie der

ˆ molekulardynamis hen Methoden (Abs hnitt 2.2)

ˆ quanten hemis hen Elektronenstruktur-Methoden (Abs hnitt 2.3)

ˆ Bere hnung statis herMoleküleigens haften(Abs hnitt 2.4)

ˆ molekularme hanis hen Me hanik (Abs hnitt 2.5)

ˆ kombiniertenquanten-/molekularme hanis hen (QM/MM) Methoden (Abs hnitt 2.6)

2.1 S hrödinger-Glei hung und Kern-Elektron-Separation

DieQuantenme hanikpostuliert,dasseinmolekularesSystemvollständigdur heineorts-und

zeitabhängigeWellenfunktion

Φ(R, τ , t)

bes hriebenwerden kann.DieWellenfunktionerfüllt

(im ni ht-relativistis hen Fall) diefür dieQuanten hemie zentrale Glei hung,

b

HΦ(R, τ , t) = i~

∂

∂t

Φ(R, τ , t),

(2.1)

diezeitabhängigeS hrödinger-Glei hung,mitderZeit

t

,demSatzder Koordinaten der

Atom-kerne

R

= {R

1

, R

2

· · · R

N

K

}

und Elektronen

τ

= {τ

1

, τ

2

, · · · τ

N

b

H

bestimmt dieFormderS hrödinger-Glei hung undist fürein abges hlossenes molekulares Systemgegeben dur h

b

H =

b

H

e

z

}|

{

−

1

₂

X

i

∇

2

i

−

X

i,A

Z

A

r

iA

+

X

i,j>i

1

r

ij

+

X

A,B>A

Z

A

Z

B

R

AB

|

{z

}

V

NN

b

T

N

z

}|

{

−

1

₂

X

A

1

M

A

∇

2

A

(2.2)

inatomaren Einheiten mit den Kernladungen

Z

A

,Kernmassen

M

A

sowie den

Elektronenin-dizes

i, j

und Kernindizes

A, B

.Lösungen derzeitabhängigen S hrödinger-Glei hung basieren

imAllgemeinenauf Separatationansätzen für dieWellenfunktion. Dabeiwerdenim

Wesentli- henzweisi hergänzende Ansätzeverfolgt,diesi himWegderSeparationunters heidenund

in der Literatur als Ehrenfest- bzw. Born-Oppenheimer-Separationsweg bespro hen werden

(siehe Abbildung 2.1und Referenz[8℄).

Φ(R, τ , t)

Φ

N

(R, t), Ψ

e

(τ )

Φ

N

(R, t), Φ

e

(τ , t)

Born-Oppenheimer-Separationsweg

Ehrenfest-Separationsweg

Abb.2.1: S hematis he Darstellungder Ehrenfest- bzw. Born-Oppenheimer-Wege zurOrt-Zeit- und

Kern-Elektron-Separation.

ImFolgendenwirdderBorn-Oppenheimer-Separationsweg alsGrundlageder

Born-Oppen-heimer-Moleküldynamik und den quanten hemis hen Elektronenstruktur-Methoden

ausführ-li her dargestellt. Der Hamilton-Operator hängt für ein abges hlossenes molekulares System

ni ht explizit vonderZeit t ab,sodasssi hmit demSeparationsansatz

Φ(R, τ , t) = Ψ(R, τ ) · f (t)

(2.3)

diezeitabhängige S hrödinger-Glei hung ineine nur zeitabhängige Glei hung

i~

∂

∂t

f (t) = Ef (t)

(2.4)

undnurortsabhängigeGlei hungseparierenlässt,diezeitunabhängigeS hrödinger-Glei hung,

b

mit dem stationären Energieeigenwert

E

. Verwendet manden Born's hen Separationsansatz

[9,10℄

Ψ(R, τ ) =

X

n

Ψ

N,n

(R) · Ψ

e,n

(R, τ )

(2.6)

übervers hiedene Zustände

n

,ergibt si hfür diezeitabhängige Gesamtwellenfunktion,

Φ(R, τ , t) = f (t)

X

n

Ψ

N,n

(R)

|

_P

{z

}

n

Φ

N,n

(R,t)

·Ψ

e,n

(R, τ ),

(2.7)

Dabeizeigensi hdiezeitabhängigen Kernwellenfunktionen

Φ

N,n

(R, t)

als zeitabhängige

Ent-wi klungskoezientenderzeitunabhängigenelektronis henWellenfunktion

Ψ

e,n

(R, τ )

.Indem

mandenSeparationsansatz indiezeitabhängige S hrödinger-Glei hung einsetzt,mit einer zu

Ψ

e,n

(R, τ )

orthonormalen elektronis hen Wellenfunktion

Ψ

e,m

(R, τ )

multipliziert und über dieelektronis henKoordinatenintegriert,erhältmandieni ht-adiabatis he

Kern-S hrödinger-Glei hung:

h

b

T

N

+ E

m

(R)

i

Φ

N,m

(R, t) +

X

n

b

T

mn

(R, τ )Φ

N,n

(R, t) = i~

∂

∂t

Φ

N,m

(R, t).

(2.8)

Darin sind die Terme

E

m

(R)

die Energieeigenwerte des elektronis hen Hamilton-Operators

b

H

e

,d.h.Lösungen derelektronis hen S hrödinger-Glei hung

b

H

e

Ψ

e,m

(R, τ ) = E

m

(R)Ψ

e,m

(R, τ ).

(2.9)

Dieni ht-adiabatis hen Kopplungsterme

T

b

mn

(R, τ )

hängen explizitvon denKern- und

Elek-tronenkoordinaten

R

und

τ

ab(siehez.B.Referenz[8℄).UmeineEntkopplungderKern-und Elektronenkoordinaten zu errei hen,werden beideradiabatis hen Näherungdie

Auÿerdiago-nalelemente von

T

b

mn

(R, τ )

verna hlässigt:

h

b

T

N

+ E

m

(R)

i

Φ

N,m

(R, t) + b

T

mm

(R)Φ

N,m

(R, t) = i~

∂

∂t

Φ

N,m

(R, t)

(2.10)

Die adiabatis he Näherung impliziert, dass si h während der Bewegung der Kerne der

elek-tronis heZustand ni ht ändertund istgere htfertigt,solange dieEnergiedierenzen der

elek-tronis hen Zustände für alle relevanten Molekülkonformationen

R

groÿ sind und keine

elek-tronis henAnregungen währendderKernbewegung bes hrieben werdensollen(fürMethoden

Bei der Born-Oppenheimer Näherung, die auf dem sehr groÿen Verhältnis von Kern- und

Elektronenmassebasiert, verna hlässigt manzudemdie Diagonalelemente

T

b

mm

(R)

[13℄:

h

b

T

N

+ E

m

(R)

i

Φ

N,m

(R, t) = i~

∂

∂t

Φ

N,m

(R, t)

(2.11)

Inderverbliebenen Kern-S hrödinger-Glei hung bilden dieEigenwerte

E

m

(R)

des

elektroni-s henHamilton-OperatorsalsodenpotentiellenBeitragzurBes hreibung derKernbewegung.

2.2 Molekulardynamis he Methoden

DieTabelle 2.1gibteineÜbersi htübervers hiedenegrundlegendemolekulardynamis he

Me-thoden, die auf unters hiedli hen theoretis hen Niveaus aufbauen (für weiterführende

MD-Te hniken siehez.B. Referenz [8℄). Als Grundlage derDissertationwird imFolgenden ausge-hend vom vorangehenden Abs hnitt diegrundlegende Theorie der

Born-Oppenheimer-Mole-kulardynamik (BOMD)dargestellt. AuÿerdemwirdderAnsatz derklassis hen MDerläutert.

Die zugrundeliegenden quanten hemis hen Elektronenstrukturmethoden bzw.

molekularme- hanis hen Methoden werdeninden darauffolgenden Abs hnitten bes hrieben.

Tab.2.1:Bezei hnungundÜbersi htübervers hiedenemolekulardynamis heMethoden.

Bezei hnung Behandlung Bes hreibung

der Atomkerne

Quanten-MD quantenme h. Quantenme h. Re hnung jenseits der

klass.Behandlung derKernbewegung

Ab-initio-MD klassis h

Ab-initio-Elektronenstruktur-Re hnung währendderMD

a)Ehrenfest-MD basiertauf Ehrenfest-Separation

b)Born-Oppenheimer-MD basiertauf Born-Oppenheimer-Separation

)Car-Parinello-MD Propagationderelektr.Wellenfunktion

klassis he MD klassis h vorparametrisierte Kraftfelder

Bei der Born-Oppenheimer-Molekulardynamik (BOMD) geht man von der

quantenme- hanis hen Bes hreibung der Kerndynamik (siehe Glei hung 2.11) über zu einer klassis hen

N

K

Glei hungen

M

A

∂

2

∂t

2

R

A

= −∇

_|

A

_{z

E

0

(R)

_}

F

A

(R)

,

(2.12)

wobei

E

0

(R)

derEigenwert des elektronis hen Hamiltonoperatorsfür den Grundzustand ist,

also

E

0

(R) = hΨ

e,0

(R, τ )| b

H

e

|Ψ

e,0

(R, τ )i,

(2.13)

fürdienormierte elektronis he Wellenfunktion

Ψ

e,0

(R, τ )

.Dergroÿe QuotientausKern-und

Elektronenmasse bildetdieGrundlagedafür,dieKerneklassis hzubes hreiben.ImZuge(on

the y) der MD werden der Eigenwert des elektronis hen Hamiltonoperators

E

0

(R)

und die

Kräfte

F

A

(R)

andenKernpositionenfür jedeMolekülkonformation

R

mittelseiner

quanten- hemis hen Elektronenstruktur-Re hnung neu bestimmt (siehe Abs hnitte 2.3 und 2.4). Die

Glei hung 2.12 bildet zuglei h die Grundlage zu den einfa heren klassis hen MD-Methoden,

in denen die Kräfte

F

A

(R)

ni ht aus quanten hemis hen Elektronenstruktur-Re hnungen,

sondernauseinfa henmolekularme hanis henEnergieausdrü kenbere hnetwerden(siehe

Ab-s hnitt2.5).EineweitereMögli hkeitbestehtdarin,dieKräfte

F

A

(R)

ausQM/MM-Methoden

zu bere hnen (sieheAbs hnitt2.6).

Dieklassis henKernbewegungs-Glei hungen werdenimAllgemeinen dur hnumeris he

In-tegrationsverfahren gelöst, wie die Glei hung 2.14 am Beispiel eines Verlet-Integrators [14℄ zeigt:

R

A

(t

0

+ h) ≈ R

A

(t

0

) + h

_∂t

∂

R

A





_t

₀

+

1

₂

h

2 F

_M

A

(t

_A

0

)

R

A

(t

0

− h) ≈ R

A

(t

0

) − h

_∂t

∂

R

A





t

0

+

1

2

h

2 F

A

(t

0

)

M

A

R

A

(t

0

+ h) ≈ 2R

A

(t

0

) − R

A

(t

0

− h) + h

2 F

_M

A

(t

_A

0

)

(2.14)

DerZeits hritthmussbeidennumeris henIntegrationsverfahrenmögli hstkleingewählt

wer-den,umdennumeris henFehlerbeiderIntegrationderKernbewegungs-Glei hungenmögli hst

klein zu halten. Glei hzeitig wä hst natürli h mit kleinem

h

die Anzahl benötigter S hritte

umeinebestimmteSimulationsdauerzuerrei hen.Übli herweise bemisstman

h

antypis hen

Zeitenders hnellstenMoleküls hwingung imSystem(im Berei h von1 fs

≈

41.3a.u.[15,16℄). Für dasBeispiel eines Verlet-Integrators [14℄ lässt si h für

t

0

= 0

die benötigte Molekülkon-formation

R

A

(−h)

z.B.bestimmenna h

R

A

(−h) = R

A

(0) − h

∂

∂t

R

A









0

| {z }

v

A

(0)

+

1

2

h

2

F

A

(0)

M

A

.

(2.15)

Die Ges hwindigkeiten

v

A

(0)

werden typis herweise als Maxwell-Boltzmann-Verteilung bei

einer bestimmten Temperatur erzeugt und legen zusammen mit der Molekülkonformation

R

A

(0)

die Anfangsbedingungen einer einzelnenTrajektorie fest.

Molekulardynamik-Simulationenerlaubendendur hdieOrtsvektorenallerAtomkerne

auf-gespanntenKonformationsraum abzutasten,umbeispielsweise

ˆ vers hiedenelokaleMinimumstrukturenoderprinzipiellau hdieglobale

Minimumstruk-tur zunden

ˆ hemis he Umwandlungen (Strukturänderungen, hemis he Reaktionen) zu verfolgen

ˆ thermodynamis her Gröÿen mittelsderstatistis hen Thermodynamik zu bere hnen.

2.3 Quanten hemis he Elektronenstruktur-Methoden

Im Allgemeinen ist die elektronis he S hrödingerglei hung (Glei hung 2.9) für ein

Mehrelek-tronensystem ni ht ges hlossen analytis h lösbar. Allerdings steht eine Hierar hie an

nume-ris hen, Wellenfunktions-basierten ab-initio Methoden zur Verfügung, mit denen im Prinzip

beliebige Genauigkeiten errei ht werdenkönnen. DasProblem derWellenfunktions-basierten

ab-initio Methoden besteht daher maÿgebli h in ihrem Re henaufwand, der von der Gröÿe

desmolekularenSystemsabhängt: Selbst die Hartree-Fo k-(HF)-Methode, alseinfa hste

Nä-herung innerhalb der Methodenhierar hie, skaliert konventionell mit der dritten Potenz der

Systemgröÿe (für Details siehe Abs hnitt 2.3.1). Die Abbildung 2.2 zeigt das

Skalenverhal-ten für einige weitere Wellenfunktions-basierte ab-initio Methoden. Als Ausgangspunkt der

weiterführenden Kapitel der Dissertation wird die Hartree-Fo k-Methode (siehe Abs hnitt

2.3.1)vorgestellt. AuÿerdemwirdinAbs hnitt2.3.3 derAnsatzderDi htefunktional-Theorie

erläutert,die ni ht auf einer elektronis hen Wellenfunktion aufbaut, sondernauf der

Bere h-nung der elektronis hen Grundzustandsenergie als Funktional der Elektronendi hte basiert.

Die Darstellung geht dabei au h auf die methodis hen S hlüssels hritte ein, mit denen man

einlineares SkalenverhaltenmitderMolekülgröÿeerzielenkann.DerAbs hnitt2.3.2 skizziert

HF MP2 CCSD CCSD(T) . . .

O(M

3

)

O(M

5

₎

O(M

6

)

O(M

7

₎

Wellenfunktions-basierte Methoden

O(M

3

₎

Di htefunktional-theorie Exakte Lösung

b

H

e

Ψ

e,0

= E

0

Ψ

e,0

Abb.2.2: Skalenverhalten

O(M

x

₎

mit der Gröÿe

M

des molekularen Systems

Wellenfunktions-basierter Elektronenstruktur-Methoden (HF: Hartree-Fo k, MP2:

Møller-Plesset-Störungstheorie2. Ordnung, CCSD: Coupled-Cluster-Singles-Doubles, CCSD(T):

Coupled-Cluster-Singles-Doublesund störungstheoretis heKorrekturder Triples) sowieEinordnung

derDi htefunktionaltheorie.

2.3.1 Die Hartree-Fo k-Methode

Umden Satzan elektronis hen Koordinaten

τ

= {τ

1

, τ

2

, · · · τ

N

e

}

derelektronis hen W

ellen-funktionzuseparieren, setzt manbeiderHF-Methode eine einzelneSlater-Determinanteder

Form

Ψ

e,0

(τ ) =

1

√

N

e

!

·

























ϕ

1

(τ

1

)

ϕ

2

(τ

1

)

· · ·

ϕ

N

e

(τ

1

)

ϕ

1

(τ

2

)

ϕ

2

(τ

2

)

· · ·

ϕ

N

e

(τ

2

)

. . . . . . . . . . . .

ϕ

1

(τ

N

e

) ϕ

2

(τ

N

e

) · · · ϕ

N

e

(τ

N

e

)

























≡ |ϕ

1

ϕ

2

· · · ϕ

N

e

i

(2.16)

an.Die Spinfunktionen

ϕ

i

hängen jeweils nur von den elektronis hen Koordinaten

τ

j

(kom-binierteRaum-/Spinkoordinaten) einesElektrons ab(

j = 1, 2, · · · , N

e

)undkönnen

orthonor-miert gewählt werden. Für ges hlossens halige Systeme kann manPaarevon Spinfunktionen

β(σ

j

)

erhalten(Restri ted-Hartree-Fo k, RHF):

ϕ

2i−1

(τ

j

) = φ

i

(r

j

) · α(σ

j

)

i = 1, 2, · · · ,

1

2

N

e

(2.17)

ϕ

2i

(τ

j

) = φ

i

(r

j

) · β(σ

j

)

InnerhalbdesSlater-Determinantenansatzes ergibt si h derEnergieerwartungswert als

E

0

= hΨ

e,0

| b

H

e

|Ψ

e,0

i =

X

i

hϕ

i

(τ

j

)|bh

j

|ϕ

i

(τ

j

)i +

1

2

X

i,j6=i

[hϕ

i

ϕ

j

|ϕ

i

ϕ

j

i − hϕ

i

ϕ

j

|ϕ

j

ϕ

i

i] + V

N N

(2.18) mit

bh

j

= −

1

2

∇

2

j

−

X

A

Z

A

r

jA

(2.19) und

hϕ

i

ϕ

j

|ϕ

i

ϕ

j

i =

Z

dτ

1

dτ

2

ϕ

∗

i

(τ

1

)ϕ

∗

j

(τ

2

)

1

r

12

ϕ

i

(τ

1

)ϕ

j

(τ

2

).

(2.20) DieMinimierungdesEnergie-Erwartungswerts auf derGrundlage desVariationsprinzips und

unter der Nebenbedingung der Orthonormalität der Molekülorbitale liefert die allgemeinen

Hartree-Fo k-Glei hungen.

b

F ϕ

i

=

X

j

ǫ

ji

ϕ

j

(2.21)

mitdenLagrange-Multiplikatoren

ǫ

ji

.Sielassensi hunitär transformieren,wodur hmandie

kanonis hen Hartree-Fo k-Glei hungen erhält (siehez.B. Referenz[17℄):

b

F ϕ

i

= ǫ

i

ϕ

i

(2.22)

DerFo k-Operator

F

b

repräsentiertalseektiverEinelektronen-Operatordiekinetis heEnergie

eines Elektrons, dessen We hselwirkung mit den Kernen (Einelektronen-Teil) und ein

gemit-teltesPotential,daseinElektrondur hdieAnwesenheitderanderenElektronenerfährt

(Zwei-elektronen-Teil). Die kanonis hen Hartree-Fo k-Glei hungen können imAllgemeinen ni ht in

einem S hritt gelöst werden, vielmehr hängt

F

b

von allen (zunä hst unbekannten)

Molekülor-bitalen

ϕ

i

ab,sodassdieGlei hungenna hdemVerfahrendesselbst-konsistentenFeldes

(self- onsistent eld, SCF) iterativ gelöst werden müssen. Für das algebrais h-iterative Vorgehen

werden die den Spinorbitalen

ϕ

i

zugrunde liegenden Raumorbitale

φ

i

(r

j

)

(siehe Glei hung

2.17) übli herweise als Linearkombination eines endli hen Satzes

{χ

µ

}

an atomzentrierten

(ni ht-orthogonalen) Basisfunktionen(Atomorbitale, AO)genähert:

φ

i

(r

j

) =

N

_X

bas

µ

FürdieAtomorbitaleverwendetmanbeimolekularenBere hnungenheuteüberwiegend

Gauÿ-typ-Funktionen derArt

χ

µp

(r

j

, R

µ

) = (x − R

µ,x

)

l

µ,x

(y − R

µ,y

)

l

µ,y

(z − R

l,z

)

l

µ,z

e

−ζ

p

(r

j

−R

µ

)

2

(2.24)

bzw.feste Linearkombinationen mehrerer primitiverGauÿtyp-Funktionen,

χ

µ

(r

j

, R

µ

) =

K

µ

X

p=1

D

µp

· χ

µp

(r

j

, R

µ

),

(2.25)

mitdemKontraktionskoezient

D

µp

,demAtomzentrum

R

µ

,demExponent

ζ

p

,den

Drehim-pulsen

l

µ,x

,

l

µ,y

und

l

µ,z

, und dem Kontraktionsgrad

K

µ

(hier gezeigt für kartesis he

Basis-funktionen).Dur hdieBasissatzentwi klunggehendiekanonis henHartree-Fo k-Glei hungen

fürden RHF-Fallindiealgebrais he Form derRoothaan-Hall-Glei hung über:

F C

= SCǫ

(2.26)

DieMatrixglei hungverknüpftdieFo k-Matrix

F

,dieMO-Koezienten-Matrix

C

,die

Über-lappungsmatrix

S

und die Matrixder Orbitalenergien

ǫ

miteinander. Die Fo k-Matrix setzt

si h ausdemEinelektronen-Teil

h

und demZweielektronen-Teil

G

zusammen,

F

= h + G(P )

(2.27)

FüreinElement desZweielektronen-Teils

G

µν

inder AO-Basis gilt fürden RHF-Fall

G

µν

(P ) =

X

κλ

2P

κλ

(µν|κλ)

|

{z

}

J

µν

−

X

κλ

P

κλ

(µλ|κν)

|

{z

}

K

µν

,

(2.28)

mitden Matrixelementen

P

κλ

derEinteil hen-Di htematrix

P

P

κλ

=

1

2

N

e

X

i

C

κi

C

λi

∗

,

(2.29)

den Elementen

J

µν

des Coulombteils bzw.

K

µν

des Austaus hteils und den

Zweielektronen-Integralen

(µν|κλ)

(ele tron repulsionintegral, ERI).DerCoulombteillässtsi halsklassis he

We hselwirkung derLadungsverteilungen

Ω

µν

(r

1

) = χ

µ

(r

1

, R

µ

) · χ

ν

(r

1

, R

ν

)

(2.30)

und

interpretieren, während der Austaus hteil als Folge des quantenme hanis hen Ansatzes

ent-steht und im klassis hen Bild allein ni ht gedeutet werden kann. Jeder Iterationss hritt des

SCF umfasst die Bere hnung der Fo k-Matrix, die Orthogonalisierung der Basis und einen

Diagonalisierungss hritt, der eine neue MO-Koezienten-Matrix

C

bzw.

Einteil hen-Di hte-matrix

P

erzeugt. Ist die Konvergenz gemäÿ eines festgelegten Grenzwerts

10

−ϑ

SCF

errei ht,

wird der Zyklus abgebro hen, ansonsten erneut dur hlaufen. Als Konvergenzkriterium kann

beispielsweisedieNorm(bzw.dasMaximalelement)derDIIS-Fehlermatrixbere hnetwerden,

ε

_{≡ F P S − SP F ,}

(2.32)

dieimFallederKonvergenz

0

wird.DIIS(dire tinversion of theiterative subspa e,[18℄) be-zei hnetdabeieinenhäugverwendetenAlgorithmuszurBes hleunigungderSCF-Konvergenz.

DieersteIterationbeginntnotwendigerweisemiteinerAnfangss hätzung,z.B.der

Superpositi-onatomarerDi hten(SAD),zudermandur hdieSummationsphäris hgemittelter,atomarer

Di htengelangt. InAbhängigkeit desEinelektronen-Teils

h

, derEinteil hen-Di htematrix

P

undderFo k-Matrix

F

lässt si h derRHF-Energieerwartungswert bere hnen als

E

0

= tr{P [2 · h + G(P )

|

{z

}

F

]} + V

NN

.

(2.33)

DerAufwand, umdie Fo k-Matrix

F

zu bere hnen,steigt formalmit dervierten Potenzder

Anzahl an Basisfunktionen

N

bas

an, da zur Bildung von

F

in einem konventionellen

SCF-Verfahren

N

4

bas

Zweielektronen-Integrale

(µν|κλ)

bere hnet und mit der Einteil hen-Di hte-matrix

P

kontrahiertwerdenmüssen.Jedo hklingendieatomzentriertenBasisfunktionenmit

zunehmendem Abstand vom jeweiligen Zentrum

R

µ

exponentiell ab.Zu jeder Basisfunktion

χ

µ

(r

j

, R

µ

)

gibt es dahernur eine konstante Zahlan Basisfunktionen

χ

ν

(r

j

, R

ν

)

,die mit ihr einenumeris hsignikanteLadungsverteilung

Ω

µν

bilden.DasIntegrals reening,basierendauf

derCau hy-S hwarz's hen Unglei hung, wurde 1989 von Häser und Ahlri hs [19℄ eingeführt und berü ksi htigtdieses Verhalten:

|(µν|κλ)| ≤ (µν|µν)

1

2

| {z }

Q

µν

(κλ|κλ)

1

2

| {z }

Q

κλ

(2.34)

Zweielektronen-Integrale werden verna hlässigt, wenn die mit

O(N

2

bas

)

-Aufwand bere hnete obere S hwarz-S hranke

Q

µν

Q

κλ

kleiner als ein festgelegter Grenzwert

10

−ϑ

Int

ist, wodur h

insgesamt das Skalenverhalten zur Bildung des Zweielektronen-Teils

G

auf

O(N

2

Darüberhinausistfüreinmolekulares,ni ht-metallis hesSystemdieEinteil hen-Di htematrix

typis herweisedünnbesetzt(sparse),d.h.dieElemente

P

κλ

klingenfürzunehmendeAbstände

der Zentren

R

κ

und

R

λ

der Basisfunktion

χ

κ

und

χ

λ

ras hab [20℄.Deshalb istes besonders vorteilhaft,dieEinteil hen-Di htematrix

P

n

deraktuellenIteration

n

indasIntegrals reening

miteinzubeziehen,

|P

κλ

n

(µκ|λν) | ≤ |P

κλ

n

|Q

µκ

Q

λν

,

(2.35)

bzw.die aktuelle Dierenzdi htematrix

∆P

n

_{= P}

n+1

_{− P}

n

|∆P

κλ

n

(µκ|λν) | ≤ |∆P

κλ

n

|Q

µκ

Q

λν

,

(2.36)

fürden Fall einer inkrementellen BildungderFo k-Matrix[21℄gemäÿ

F

n+1

= F

n

+ G(∆P

n

).

(2.37)

FürdenAustaus hteil, indenen dieEinteil hen-Di htematrix den Bra-und Ketteilder

Zwei-elektronen-Integrale miteinander koppelt, kann man dann in einem weiten Berei h an

Mo-lekülgröÿenein lineares Skalenverhalten erzielen, mittels eines linear-skalierenden S reenings

(LinK,[22,23℄). DasSkalenverhaltendesCoulombteils lässtsi hmithilfedervonWhiteetal. entwi kelten Continuous fast multipole method(CFMM, [24℄) auflinear reduzieren.

Asymptotis h wird das Skalenverhalten eines konventionellen SCF-Verfahrens dur h den

AufwandzurDiagonalisierungderFo k-Matrixbestimmt(

O(N

3

bas

)

).Jedo histesinzwis hen mögli h mit diagonalisierungsfreien Algorithmen (siehe beispielsweise Referenz [2℄ und darin enthaltene Referenzen) den

O(N

3

bas

)

-S hritt zu umgehen. Dabei sind die diagonalisierungs-freien Algorithmeninsbesonderefür sehrgroÿemolekulareSysteme (mehreretausendAtome)

wi htig, dader

O(N

3

bas

)

-S hritt typis herweise einen kleinen Vorfaktor hat und daher vergli- hen mit derBildungderFo k-Matrixerst spätzumzeitbestimmenden S hrittwird.

2.3.2 Post-Hartree-Fo k-Methoden

In der Hartree-Fo k-Theorie wird die Mehrelektronen-Wellenfunktion als eine einzelne, aus

Spinorbitalen zusammengesetzte Slater-Determinante (siehe Glei hung 2.16) genähert. Als

Konsequenz bes hreibt die Hartree-Fo k-Theorie die We hselwirkung eines jeden Elektrons

näherungsweise als We hselwirkung mit einem gemittelten Feld der anderen Elektronen; die

einem Ort

r

1

zu nden und einanderes glei hzeitig amOrt

r

2

ergibt si h als einfa hes

Pro-dukt derjeweils separatenAufenthaltswahrs heinli hkeiten. DerAusdru k

Korrelationsener-gie istdeniertals Dierenzzwis hen derexaktenEnergie undderHartree-Fo k-Energie im

Basissatz-Limit.

Bei der HF-Methode sind typis herweise mehr Molekülorbitale zugängli h als für die

Be-setzung mit Elektronen benötigt. Während in die HF-Slater-Determinante nur die

Spinor-bitale mit den niedrigsten Orbitalenergien eingehen, können weitere Slater-Determinanten

dur h andere Orbital-Kombinationen (Kongurationen) erzeugt werden. Prinzipiell lässtsi h

die im Rahmen der Basis optimale Wellenfunktion als Linearkombination aller mögli hen

Slater-Determinanten ansetzen (Full onguration intera tion, FCI). Aufgrund des

Re hen-aufwands muss man si h für die allermeisten molekularen Systeme auf bestimmte

Kombi-nationen eins hränken, sei es im Sinne von abgebro henen CI- oder

Coupled-Cluster-(CC)-Methoden(siehez.B.Referenzen[2527℄)mitunters hiedli hen Vor-undNa hteilen.Eine wei-tere Elektronenkorrelations-Methode, die Møller-Plesset-Störungstheorie, geht ebenfalls von

derHartree-Fo k-Methode ausundleiteteineKorrekturfürdieElektronenkorrelationauf

Ba-sis der zeitunabhängigen Störungstheorie her[17,28℄. Die CCSD(T)-Methode basiert sowohl aufeinerabgebro henen CC-Entwi klung undalsau heiner störungstheoretis henKorrektur

(siehez.B.Referenzen [29,30℄).

2.3.3 Di htefunktionaltheorie

Die Grundlage der Di htefunktionaltheorie (DFT) ist der Beweis von Hohenberg und Kohn

[31℄, dass die elektronis he Energie des Grundzustands

E

0

ni ht nur als Funktional der elek-tronis hen Wellenfunktion (siehe Glei hung 2.13) bestimmt ist, sondern au h vollständig als

Funktional

E[ρ]

der Elektronendi hte

ρ(x, y, z)

. Anstelle die elektronis he

S hrödingerglei- hung(näherungsweise)dur h Wellenfunktions-basierte Methodenzulösen, bietetsi hdamit

prinzipiell dieMögli hkeit,eine gültigeBes hreibung des imAllgemeinen unbekannten

Ener-giefunktionals

E[ρ]

herzuleiten. Der formale Vorteil in der Formulierung von

E

0

als F

unk-tional der Elektronendi hte

ρ(x, y, z)

liegt dabei darin, dass die Elektronendi hte nur von

drei Raumkoordinaten abhängt, wohingegen die elektronis he Wellenfunktion eine Funktion

der Koordinaten aller

N

e

Elektronen des Systems ist. Die gegenwärtig angewendeten

dasvonKohn und Sham1965 eingeführtwurde,umdie S hwierigkeiten beiderdirekten

Be-stimmung deskinetis hen Anteilsvon

E[ρ]

zuüberwinden[32℄.DerKohn-Sham-Ansatz (KS) lässt si h analog zur Hartree-Fo k-Theorie formulieren, wobei nur der Austaus hteil

K

im

Zweielektronen-Teil derFo k-Matrix(sieheGlei hung 2.28)ersetztwerdenmuss, dur heinen

Austaus h-Korrelationsteil

V

XC

mitden Elementen

V

XC

_µν

_{= hµ|}

∂E

XC

∂ρ

|νi

(2.38)

undderAustaus h-Korrelationsenergie

E

XC

=

Z

f

XC

dxdydz.

(2.39)

Die Glei hung 2.39 muss man in der Regel dur h numeris he Verfahren integrieren (für ein

linear-skalierendes Verfahren siehe z.B. Referenz [33℄). Die exakte Form desFunktionals

f

XC

ist im Allgemeinen ni ht bekannt und die heute verwendeten Funktionale lassen si h gemäÿ

ihrer grundlegenden Abhängigkeiten von derElektronendi hteindrei Klassen einteilen[34℄:

ˆ L(S)DA-Funktionale

ˆ GGA-Funktionale

ˆ Meta-GGA-Funktionale

Die einfa hsten Näherungen innerhalb der vorges hlagenen heuristis hen Rangfolge [35,36℄, dieL(S)DA-Funktionale,basierenauf derlo al (spin) densityapproximation[37,38℄und hän-gennurvon derElektronendi hte selbstab.Die GGA-Funktionale (z.B.PBE [39℄)bauen auf dergeneralized gradient approximation auf,d.h. hängen zudem von derersten Ableitungder

Di hteab. Als Meta-GGAbezei hnet manFunktionale, diezusätzli h entweder von höheren

Ableitungen oder der kinetis hen Orbitalenergiedi hte

τ

[40℄ abhängen. Eine Verbesserung über die reinen Funktionale hinaus errei ht man typis herweise, indem man in

Hybridfunk-tionalen (z.B B3LYP [41,42℄, PBE0 [43℄, B97-2 [44℄) den Austaus h-Korrelationsteil

V

XC

auseinem Anteil

γ

exaktem Hartree-Fo k-Austaus h

K

und einem reinen DFT-Anteil

V

′

XC

zusammensetzt:

2.4 Bere hnung statis her Moleküleigens haften

Die vorausgehenden beiden Abs hnitte stellen Methoden vor, mit denen man die

elektroni-s heWellenfunktionbzw. EnergiedesGrundzustandsfüreinegegebene Molekülkonformation

und inAbwesenheit externer elektris herodermagnetis herFelder (näherungsweise)

bestim-men kann. Ausgehend von der elektronis hen Grundzustandsenergie lassensi h viele weitere

molekulare Eigens haften als Antwort desmolekularen Systems auf eine Störung bere hnen.

Die Störung kann dabei sowohl einen externen Ursprung haben (z.B. ein elektris hes oder

magnetis hes Feld) oderau h einen internen Ursprung (z.B. eine Kernverrü kung oder

kern-magnetis hesMoment)undkannprinzipiellzeitabhängig (dynamis heEigens haft)oder

zeit-unabhängig (statis he Eigens haft) sein. Als Voraussetzung der na hstehenden Kapitel wird

hier nur auf statis he Moleküleigens haften (z.B. molekulare Kräfte, NMR-Abs hirmungen)

eingegangen.BeideraufEnergieableitungenbasierendenStörungsmethodesetztmandie

elek-tronis heEnergie

E(ξ)

fürdenFalleiners hwa hen Störung

ξ

alsTaylor-Entwi klungumden

ungestörten Fall (0)an:

E(ξ) = E(ξ = 0)+

∂E

∂ξ









_ξ=0

ξ+

1

2

∂

2

E

∂ξ

2









_ξ=0

ξ

2

+

1

6

∂

3

E

∂ξ

3









_ξ=0

ξ

3

+· · · =

∞

X

k=0

1

k!

∂

k

E

∂ξ

k









_ξ=0

·ξ

k

(2.41) Allgemeinergilt dann fürdieStörungen dur heinexternes elektris hesFeld

Q

,magnetis hes

Feld

B

,einer Kernverrü kung

R

A

odereinem kernmagnetis hen Moment

µ

A

:

E(Q, B, R

A

, µ

A

) =

∞

X

k=0

∞

X

l=0

∞

X

m=0

∞

X

n=0

1

k! l! m! n!

∂

k

_∂

l

_∂

m

_∂

n

_E

∂Q

k

_∂B

l

_∂R

m

A

∂µ

n

A









Q=0,B=0,

R

A

=0,µ

A

=0

·Q

k

·B

l

·R

m

A

·µ

n

A

.

(2.42)

Indemmandiegrundlegendenphysikalis henAbhängigkeitenderWe hselwirkung

berü ksi h-tigt, lassensi hdiemolekularen Eigens haften alspartielleEnergieableitung

∂

k

∂

l

∂

m

∂

n

E

∂Q

k

_∂B

l

_∂R

m

A

∂µ

n

A









Q=0,B=0,

R

A

=0,µ

A

=0

inderTaylor-Entwi klung identizieren (siehez.B.Referenzen[1,45℄).Soliefertbeispielsweise dernegative Gradient derelektronis hen Energie na h der Kernkoordinate

R

A

die Kraft auf

den AtomkernA:

F

A

= −

∂

∂R

A

| {z }

∇

A

E

0

(2.43)

Diegemis hte zweite Ableitungna h äuÿerem Magnetfeld

B

und kernmagnetis hen Moment

µ

_A

ergibt denkernmagnetis hen Abs hirmtensor

σ

A

amAtomA:

σ

A

=

∂

2

E

0

∂B∂µ

A

=













∂

2

_E

₀

∂B

x

_∂µ

x

A

∂

2

_E

₀

∂B

x

_∂µ

y

A

∂

2

_E

₀

∂B

x

_∂µ

z

A

∂

2

_E

0

∂B

y

_∂µ

x

A

∂

2

_E

0

∂B

y

_∂µ

y

A

∂

2

_E

0

∂B

y

_∂µ

z

A

∂

2

_E

0

∂B

z

_∂µ

x

A

∂

2

_E

0

∂B

z

_∂µ

y

A

∂

2

_E

0

∂B

z

_∂µ

z

A













(2.44)

Dieisotrope Abs hirmung ergibt si hausdem Abs hirmtensor als

σ

A

=

1

3

· tr(σ

A

)

(2.45)

undNMR-Vers hiebungen gegenübereiner Referenz(Ref.) als

δ

A

=

ν

A

− ν

Ref.

ν

Ref.

=

(1 − σ

A

) − (1 − σ

Ref.

)

1 − σ

Ref.

≃

σ

Ref.

− σ

A

1

,

(2.46)

mitden Larmorfrequenzen

ν

A

und

ν

Ref.

.

DieTabelle 2.2gibteineÜbersi ht überweitereausgewählte molekulareEigens haftenbis

zur zweiten Ordnung.

Tab.2.2:Übersi htüber die Bere hnung ausgewählter molekularer Eigens haften als

k

-te Ableitung

na hdemelektris henFeld

Q

,

l

-teAbleitungna hdemmagnetis henFeld

B

,

m

-teAbleitung

na h einer Kernkoordinate

R

A

n

-te Ableitung na h einem kernmagnetis hem Moment

µA

derelektronis henEnergie

Molekulare Eigens haft

k

l

m

n

Elektris hesDipolmoment 1 0 0 0 Magnetis hes Dipolmoment 0 1 0 0 Hyperfeinkopplungskonstante 0 0 0 1 Molekulare Kräfte 0 0 1 0 Elektris he Polarisierbarkeit 2 0 0 0 Magnetisierbarkeit 0 2 0 0 Kernspin-Kopplung 0 0 0 2 Harmonis he S hwingungsfrequenzen 0 0 2 0 Kernmagnetis he Abs hirmung 0 1 0 1

2.4.1 Störungen erster Ordnung

Ausgehend von der Glei hung des RHF-Energieerwartungswerts 2.33 (oder dem analogen

Kohn-Sham-Ausdru k) ergibt si h für dieerste Ableitungna h einer Störung

ξ

∂

∂ξ

E

0

≡ E

ξ

0

= tr{P [2 · h + G(P )]}

ξ

+ V

ξ

NN

(2.47)

= tr{P

ξ

[2 · h + G(P )] + P [2 · h

ξ

+ G

ξ

(P ) + G(P

ξ

_{)]} + V}

_NN

ξ

= 2 · tr{P

ξ

[h + G(P )

|

{z

}

F

]} + 2 · tr{P h

ξ

} + tr{P G

ξ

(P )]} + V

NN

ξ

,

wobei

G

ξ

_{(P )}

deniertistals

G

ξ

_µν

(P ) =

X

κλ

P

κλ

n

2 (µν|κλ)

ξ

− (µλ|κν)

ξ

o

(2.48)

undausgenutzt wurde,dass

tr{P G(P

ξ

)} = tr{P

ξ

G

_{(P )}.}

(2.49)

Während si h in der Glei hung 2.47 die Ableitungen des Einelektronen-Teils

h

ξ

, der

Zwei-elektronen-Integrale

(µν|κλ)

ξ

und der Kern-Kern-We hselwirkungsenergie

V

ξ

NN

ges hlossen bere hnen lassen, kann man die gestörte Einteil hen-Di htematrix

P

ξ

in der Regel nur

ite-rativbestimmen(siehe Abs hnitt 2.4.2). Allerdings lassensi h dieEnergieableitungen erster

Ordnung so umformulieren, dass die Abhängigkeit von

P

ξ

entfällt. Dazu kann man von der

Zerlegung indie Subraumprojektionen,

P

ξ

= P SP

_|

_{z

ξ

SP

_}

P

ξ

oo

+ P SP

ξ

_{(1 − SP )}

|

{z

}

P

ξ

ov

+ (1 − P S)P

ξ

SP

|

{z

}

P

ξ

vo

+ (1 − P S)P

ξ

(1 − SP )

|

{z

}

P

ξ

vv

,

(2.50)

ausgehenund dieIdempotenzbedingung

P

= P SP ,

(2.51)

derenabgeleiteteForm

P

ξ

= P

ξ

SP

+ P S

ξ

P

+ P SP

ξ

(2.52)

unddieSCF-Konvergenzbedingung (siehe Glei hung 2.32)

berü ksi htigen.FürdieProjektionaufdenbesetzt-besetzten(o upied-o upied)Raumergibt si h

P

ξ

_oo

= P SP

ξ

SP

= P S

n

P

ξ

SP

+ P S

ξ

P

+ P SP

ξ

o

SP

(2.54)

= P SP

ξ

S P SP

_{| {z }}

P

+ P SP

_{| {z }}

P

S

ξ

P SP

_{| {z }}

P

+ P SP

_{| {z }}

P

SP

ξ

SP

= P SP

_|

_{z

ξ

SP

_}

P

ξ

oo

+P S

ξ

P

+ P SP

_|

_{z

ξ

SP

_}

P

ξ

oo

undfolgli h

P

ξ

_oo

_{= −P S}

ξ

P

.

(2.55)

Derunbesetzt-unbesetzte (v irtual-v irtual)Blo kvon

P

ξ

vers hwindet:

P

_vv

_{= (1 − P S)P}

ξ

_{(1 − SP ) = (P}

ξ

_{− P SP}

ξ

_{)(1 − SP )}

(2.56)

= P

ξ

_{− P}

ξ

SP

_{− P SP}

ξ

|

{z

}

P S

ξ

P

+ P SP

ξ

SP

|

{z

}

−P S

ξ

_P

= 0

Berü ksi htigt mandieGlei hungen 2.55und2.56,lässtsi hdervon

P

ξ

abhängige Termaus

der Glei hung 2.47nun s hreiben als

tr{P

ξ

F

_{} = tr{P}

ξ

_vo

F

+ P

ξ

_ov

F

_{− (P S}

ξ

P

_{)F }}

(2.57)

= tr{(1 − P S)P

ξ

SP F

+ P SP

ξ

_{(1 − SP )F − P S}

ξ

P F

_}

= tr{P

ξ

SP F

_{(1 − P S)}

|

{z

}

F

ov

+P

ξ

_{(1 − SP )F P S}

|

{z

}

F

vo

−P S

ξ

P F

_}

= −tr{P S

ξ

P F

_},

denn innerhalb derSpur lassensi h Matrixmultiplikationen zyklis h permutieren und esgilt

alsFolge desBrillouin-Theorems[46℄:

F

vo

= F

ov

= 0.

(2.58)

Insgesamt ergibt si h damit für dieEnergieableitung ersterOrdnung ausderGlei hung 2.47

E

₀

ξ

_{= −2 · tr{P S}

ξ

P F

_{} + 2 · tr{P h}

ξ

_{} + tr{P G}

ξ

_{(P )]} + V}

_NN

ξ

,

(2.59)

sodass Eigens haften ersterOrdnung (z.B.Dipolmomente, Kräfte)bere hnetwerdenkönnen

ohnedass diegestörte Einteil hen-Di htematrix

P

ξ

Während die Glei hung 2.59 allgemein für eine beliebige Störung gilt, hängt der Aufbau der Terme

S

ξ

,

h

ξ

,

G

ξ

_{(P )}

und

V

ξ

NN

vonder konkretenStörung

ξ

ab.FürKraftbere hnungen sinddieAbleitungen

S

ξ

,

h

ξ

,

G

ξ

_{(P )}

und

V

ξ

NN

unglei hnull,dennallevierGröÿen

S

,

h

,

G(P )

und

V

NN

hängenvondenKernkoordinatenab,fürdieKomponenteeinesKraftvektors

F

j

A

(mit

i = x, y, z

amAtomA)gilt also:

F

_A

i

_{= −E}

R

i

A

0

= 2 · tr{P S

R

i

A

P F

} − 2 · tr{P h

R

i

A

} − tr{P G

R

i

A

(P )]} − V

R

i

A

NN

(2.60)

Im Falle derAbleitung na h derKomponente eines kernmagnetis hen Moments

µ

j

A

(mit

j =

x, y, z

) isthingegenderTerm

V

ξ

NN

null.Ebensofallen dieTerme

S

ξ

und

G

ξ

_{(P )}

weg,denndie

Atomorbitale

χ

µ

hängen,beidertypis henWahlimRahmenvonGIAO-Methoden,ni ht von

denkernmagnetis henMomentenab;dieAbhängigkeitvon

µ

j

A

wirdalsonurim Einelektronen-Teil

h

(als Operatorableitung) bes hrieben:

E

µ

j

A

0

= 2 · tr{P h

µ

j

A

}

(2.61)

DabeiverstehtmanunterGIAOEi hursprungs-eins hlieÿendeAtomorbitale(gaugein luding

atomi orbitals, [4749℄), die standardmäÿig die Magnetfeldabhängigkeit einbeziehen (ni ht aber die Abhängigkeit von kernmagnetis hen Momenten), um so den Fehler dur h

Ei hur-sprungsvarianzfürabgebro hene Basisentwi klungen zureduzieren (füreineausführli he

Dis-kussionsiehe z.B.Referenz [50℄).

2.4.2 Störungen zweiter Ordnung, CPSCF-Verfahren

Um Störungen zweiter Ordnungen zu bere hnen, benötigt man neben der ungestörten

Ein-teil hen-Di htematrix

P

(ausdemSCF-Verfahren) diegestörteEinteil hen-Di htematrix

P

ξ

.

So folgt beispielsweise ausgehend von der Glei hung 2.61 für die Bere hnung eines Elements

des kernmagnetis hen Abs hirmungstensors

σ

ij

A

:

σ

_A

ij

=

∂

2

_E

0

∂B

i

_∂µ

j

A

≡ E

B

i

,µ

j

A

0

=

∂

∂B

i

h

2 · tr{P h

µ

j

A

}

i

= 2 · tr{P

B

i

h

µ

j

A

} + 2 · tr{P h

B

i

, µ

j

A

}

(2.62)

Diegestörte Einteil hen-Di htematrix

P

ξ

kannfür HF(und DFT-Funktionalen mitexaktem

HF-Austaus h) ni ht in einem S hritt bestimmt werden, sondern nur iterativ,

übli herwei-se mithilfedesCoupled-perturbed-self- onsistent-eld-(C PSCF) -Verfahr ens. Die dem

CPSCF-Verfahren zugrunde liegenden CPSCF-Glei hungen können inderDi hte-basierten

SCF-Konvergenzbedingung muss für eine beliebige Störung

ξ

bestehen bleiben, es muss also

gelten

∂

∂ξ

(F P S − SP F ) = 0

(2.63)

unddaraus folgt dur h Umstellen die Di hte-basierte Formulierung der CPSCF-Glei hungen

(für eine Molekülorbital-basierte Darstellung siehe z.B. Referenz [45℄ und darin enthaltene Referenzen):

F P

ξ

S

_{− SP}

ξ

F

|

{z

}

A

ξ

₁

+ G(P

ξ

)P S − SP G(P

ξ

)

|

{z

}

A

ξ

₂

= SP {

F

(ξ)

z

}|

{

h

ξ

+ G

ξ

(P )} − F

(ξ)

P S

+ S

ξ

P F

_{− F P S}

ξ

|

{z

}

b

ξ

(2.64)

Aufderre htenSeitederCPSCF-Glei hungstehendievon

P

ξ

unabhängigen,d.h.im

CPSCF-Verfahren konstanten, Terme (

b

ξ

-Teil). Die linke Seite hängt von

P

ξ

ab und muss in jeder

CPSCF-Iteration neu bere hnet werden. Der

A

ξ

1

-Teil ist dabei als Matrixprodukt der ak-tuellen gestörten Einteil hen-Di htematrix

P

ξ

mit der Fo kmatrix und Überlappmatrix

zu-gängli h. Für den

A

ξ

2

-Teil muss der Zweielektronen-Teil

G(P

ξ

₎

als Di htekontraktion der

Zweielektronen-Integrale gebildet werden (siehe Glei hung 2.48). Na h der Bere hnung des

A

ξ

₂

-Teilswird injederCPSCF-Iteration eineneuegestörteEinteil hen-Di htematrixgebildet. Dazu kannmandasGlei hungssystem

F P

ξ

S

_{− SP}

ξ

F

|

{z

}

A

ξ

₁

= b

ξ

_{− A}

ξ

₂

(2.65)

beispielsweise mittels eines konjugierten Gradientenverfahrens [6℄ lösen. Alternativ kann die neue gestörte Einteil hen-Di htematrix na h dem Di hte-basierteLapla e-transformierte

n-CPSCF-Verfahren(DL-CPSCF)au hdirektbestimmtwerden(sieheAbs hnitt 5.4.1und

Re-ferenz [7℄). Die erste CPSCF-Iteration gründet auf einer Anfangss hätzung, typis herweise beginnt dasCPSCF-Verfahren ausgehend nur vom

b

ξ

-Teil (der

A

ξ

2

-Teil auf derre hten Seite von derGlei hung 2.65 ist null). Die Konvergenz des CPSCF-Verfahren kann wie beim SCF

dur hdieDIIS-Methode[18℄bes hleunigtwerden.DieCPSCF-DIIS-Fehlermatrixindern-ten Iterationkannanalog zumSCF-Kriterium als

ε

(n)

_{≡ A}

ξ

₂

(n)

₋

A

ξ

₂

(n−1)

z

}|

{



b

ξ

_{− A}

ξ

₁

(n)



(2.66) bere hnet werden.

DieDi hte-basiertenFormulierungenermögli hen fürni htmetallis heSystemeeinlineares

Skalenverhalten der Re henzeit, wenn sie mit Sparse-Algebra-Algorithmen kombiniert

wer-den,d.h. mit Algorithmen[51℄, diedünne Besetzungen (sparsity) vonMatrizen (z.B.von

P

) ausnutzenkönnen [52℄.

2.5 Molekularme hanis he Methoden

Neben denquantenme hanis hen (QM)Methoden(sieheAbs hnitt 2.3) bildetdie

Molekular-me hanik (MM) dieGrundlagedervorgestellten QM/MM-Methoden und-Bere hnungen der

folgenden Kapitel. Daher sollen hier die Kernaspekte der molekularme hanis hen Methoden

vorgestelltwerden. Andersalsdiequanten hemis hen Methoden(z.B.Hartree-Fo koder

KS-DFT) basieren die molekularme hanis hen Methoden (oder au h Kraftfeldmethoden) ni ht

aufeiner Bes hreibung derElektronenstruktur,sondernverbleiben aufeineratomarenEbene:

Die elektronis he Energie wird ni ht als Funktional der elektronis hen Wellenfunktion oder

Elektronendi hte bestimmt, sondern als eine vorparametrisierte Funktion

E

MM

der

Atomko-ordinaten

R

.DieFunktion

E

MM

setztsi htypis herweise zusammenausEnergiebeiträgenfür

die Dehnung von Bindungen

E

1↔2

, der Beugung von Bindungswinkeln

E

1↔3

, der Rotation

von Bindungsebenen

E

1↔4

(Bindungswe hselwirkungen), sowie aus Van-der-Waals-Termen

E

vdW

undCoulomb-We hselwirkungs-Termen

E

qq

(Ni ht-Bindungswe hselwirkungen):

E

MM

=

N

_X

1↔2

n

E

_n

1↔2

|

{z

}

E

1↔2

+

N

_X

1↔3

n

E

_n

1↔3

|

{z

}

E

1↔3

+

N

_X

1↔4

n

E

1↔4

_n

|

{z

}

E

1↔4

+

N

_X

VdW

n

E

_n

vdW

|

{z

}

E

VdW

+

N

qq

X

n

E

_n

qq

| {z }

E

qq

(2.67)

Die S hreibweise 

1 ↔ 2

, 

1 ↔ 3

 bzw. 

1 ↔ 4

 gibt dabei an, ob si h die Terme über

zwei, drei bzw. vier miteinandergebundene Atome erstre ken. Die Formen der Energieterme

variierenfürvers hiedeneKraftfelder.Typis herweisewerdendieTerme

E

1↔2

n

und

E

1↔3

n

dur h einfa he Potenzreihen-Entwi klung en umeinen festgelegtenGlei hgewi hts-Bindungsabstand

R

equ

n

bzw. -winkel

Θ

equ

n

angesetzt, für das Beispiel des Amber-Kraftfeldes [53℄ als einfa he quadratis he Funktionen

E

_n

1↔2

= k

n

(R

n

− R

equ

n

)

2

,

(2.68)

E

_n

1↔3

= k

_n

′

(Θ

n

− Θ

equ

n

)

2

,

mit dem Bindungsabstand

R

n

, dem Bindungswinkel

Θ

n

und den Kraftfeldparametern

k

n

,

periodis hen Funktionen bes hrieben,

E

_n

1↔4

= k

′′

_n

(1 + cos[ω

n

]) + k

′′′

n

(1 − cos[2ω

n

]) + · · ·

(2.69)

mit den Diederwinkeln

ω

n

, teilweise zusätzli h mit Korrekturen für Änderungen von

Auÿer-Ebenen-Winkeln. Van-der-Waals-We hselwirkungen zwis hen Atompaaren im Abstand

R

n

sindhäug alsLennard-Jones-Potential [54℄wiedergegeben

E

VdW

_n

=

c

n

R

12

n

−

c

′

_n

R

6

n

(2.70)

unddielangrei hweitigen Coulomb-We hselwirkungen im einfa hsten Fall als

E

_n

qq

=

q

n,1

· q

n,2

R

n

(2.71)

mit den partiellen Punktladungsparametern

q

n,1

und

q

n,2

. Die Kraftfeldparameter werden

dur hFittenanexperimentelleDatenodermithilfevonab-initioRe hnungenAtomtyp-basiert

oder für Einzelfragmente (z.B. Aminosäuren, Nukleotide et .) bestimmt (für eine Übersi ht

undEinteilungderKraftfelder siehe z.B.Referenz[1℄).

DieKräftefür Molekulardynamik-Re hnungen (Abs hnitt 2.2) oder

Geometrieoptimierun-genbere hnet mananalogzu derGlei hung 2.12dann näherungsweise als

F

A

= −∇

A

E

MM

.

(2.72)

Nebenderprinzipiellen Fragena hderVerfügbarkeitundQualität derKraftfeldparameterfür

ein bestimmtes System lassensi h einigebesondere Problematiken einer einfa hen

Kraftfeld-bes hreibung festma hen:

ˆ Die Bindungstypen müssen vor derRe hnung festgelegt werden und können si h ni ht

ändern

ˆ Kopplungen zwis hen den Energietermen sindverna hlässigt (z.B. insbesondere

Polari-sationseekte beidenelektrostatis hen We hselwirkungen)

ˆ Glei hgewi htsferneStrukturen werdenin derRegels hle ht bes hrieben

ˆ Chemis he Reaktionen (z.B. s hon einfa he Protonenübertragungen) können im

2.6 QM/MM-Methoden

DerAbs hnitt2.3führtquanten hemis heMethodenein,mit derenHilfe molekulareSysteme

auf hohem theoretis hen Niveau bes hrieben werden können. Neben den quanten hemis hen

Methoden ermögli hen molekularme hanis he Methoden eine eziente Bes hreibung auf

ei-nemniedrigeren theoretis hen Niveau,wobeiallerdingsteilweise dieGenauigkeiten ni ht

aus-rei hen (sieheAbs hnitt 2.5). Umeine innere Molekülregion (z.B. ein gelöster Sto oderein

aktives Enzymzentrum) mögli hst genau zu bes hreiben und glei hzeitig groÿe Umgebungen

mögli hst ezient zu berü ksi htigten, kombinieren QM/MM-Ansätze diequanten hemis he

mitdermolekularme hanis hen Bes hreibung.

a)

i o

b) )

Abb.2.3: a) S hema einer inneren (i nner) und äuÿeren (o uter) Molekülregion als Grundlage des

QM/MM-Ansatzes;b)BeispieleinesgelöstenMoleküls[55℄ (innereRegion)inWasser (äu-ÿereRegion); )BeispieleinesDNA-Enzym-Komplexes[56℄;dasaktiveZentruminderNähe desblaugezeigtenResiduums bildetdieinnere Region,die DNA-Proteinumgebung die

äu-ÿereMolekülregion.

Prinzipiell lässt si h der QM/MM-Ansatz sowohlin einer additiven als au h subtraktiven

Weise formulieren [57℄. Die Bere hnungen und Untersu hungen der na hfolgenden Kapitel basieren auf der additiven Formulierung. Ausgehend vomquantenme hanis hen Standpunkt

setzt man für die additive Formulierung den Hamilton-Operator des Gesamtsystems formal

zusammen als

b

H = b

H

i

+ b

H

o

+ b

H

i,o

.

(2.73) DieOperatoren

H

b

i

und

H

b

o

bes hreibenjeweilsdieHamilton-Operatorenderisolierteninneren

bzw.äuÿeren Molekülregion und

H

b

i,o