Das Madelung-Problem

Anhang

7.1 Das effektive Einteilchenpotential

Bei gegebenem Gitter{T}sei ein Ortsvektor in derT-ten Zelle mitR+rbezeichnet (siehe auch Fig. 7.1) wobeiR zu einer Atomplatzposition zeige und r einen Punkt in der Kugel anR bezeichne. Dann kann das Potential anRals

VRσ(r) =V_Rσ^inter(r) +V_Rσ^intra(r) (7.1) geschrieben werden wobei eine Aufteilung in inter-atomare und intra-atomare Bei-tr¨age erfolgt

V_Rσ^intra(r) = −2Z_R r +

Z

(R)

2ρ_R(r⁰)

|r−r⁰|d³r⁰+V_xc,Rσ(r), (7.2) V_Rσ^inter(r) = X

R⁰

0Z

(R⁰)

2ρ_R⁰(r⁰)

|r+R−r⁰−R⁰|d³r⁰−X

R⁰

0 2Z_R⁰

|r+R−R⁰|. (7.3) Ein Strich an den Summenzeichen bedeutetR6=R⁰.

Der inter-atomare Beitrag (7.3) zum Einteilchenpotential (7.1) an R+r, ver-ursacht durch die Gesamtheit aller Ladungsdichten anR⁰+r⁰, wird als Madelung-Potential bezeichnet und soll mitVM,R(r) notiert werden.

Nun erfolgt eine Multipolentwicklung dieser Beitr¨age zum Einteilchenpotential Unter Benutzung von (7.23) und (7.21) ergibt sich f¨ur den ersten Term in (7.3):

X

R⁰

0Z

(R⁰)

2ρ_R⁰(r⁰)

|r+R−r⁰−R⁰|d³r⁰=

= 8πX

R⁰

0X

L

Z

(R⁰)

ρ_R⁰(r⁰)

|r⁰−R+R⁰|^l+1 r^l

2l+ 1YL(ˆr)Y_L^∗(r⁰−\R+R⁰)d³r⁰

= 8πX

R⁰

0X

LL⁰

4π r^l

2l+ 1Y_L(ˆr)Y_L^∗00(R\−R⁰)

|R−R⁰|^l00+1C_LL^{∗ 0}_L⁰⁰ (−1)^l(2l⁰⁰−1)!!

(2l−1)!!(2l⁰+ 1)!! ·

· Z

(R⁰)

r^0lY_L^∗0(ˆr⁰)ρ_R⁰(r⁰)d³r⁰

wobei die Summen auf die Terme mitl⁰⁰ =l⁰+lund m⁰⁰=m⁰−mbeschr¨ankt sind.

Definiert man jetzt ein MultipolmomentQRL welches sich aus der vollen Ladungs-dichte mittels

QRL=

r 4π 2l+ 1

Z

(R)

r^lY_L^∗(ˆr)ρR(r)d³r (7.4) ergibt, so kann dies als

8πX

R⁰

0X

LL⁰

4πr^l

2l+ 1Y_L(ˆr)Y_L^∗00(R\−R⁰)

|R−R⁰|^l00+1C_LL^∗ 0L⁰⁰

(−1)^l(2l⁰⁰−1)!!

(2l−1)!!(2l⁰+ 1)!!

r2l⁰+ 1 4π Q_R⁰_L⁰ geschrieben werden. Nutzt man jetzt (7.23), kann der zweite Term in (7.3) folgen-dermaßen entwickelt werden (Y_L(−ˆr) = (−1)^lY_L(ˆr)):

−X

R⁰

0 2ZR⁰

|r+R−R⁰| =−X

R⁰

0X

L

8πZR⁰

2l+ 1

(−r)^l

|R−R⁰|^l+1YL(ˆr)Y_L^∗(R\−R⁰) Damit ergibt sich schließlich f¨ur das Madelung-Potential:

V_M,R(r) = 8πX

R⁰

0X

LL⁰

4π 2l+ 1

(−r)^l

|R−R⁰|^l00+1Y_L(ˆr)Y_L^∗00(R\−R⁰)C_LL^∗ 0L⁰⁰·

· (2l⁰⁰−1)!!

(2l−1)!!(2l⁰+ 1)!!

r2l⁰+ 1 4π Q_R⁰_L⁰

−X

R⁰

0X

L

8πZ_R⁰ 2l+ 1

(−r)^l

|R−R⁰|^l+1YL(ˆr)Y_L^∗(R\−R⁰) Definiert man

q_RL=Q_RL−Z_Rδ_L,0, (7.5)

reduziert sich dies zu V_M,R(r) = 8πX

R⁰

0X

LL⁰

4π 2l+ 1

(−r)^l

|R−R⁰|^l00+1Y_L(ˆr)Y_L^∗00(R\−R⁰)C_LL^∗ 0L⁰⁰·

· (2l⁰⁰−1)!!

(2l−1)!!(2l⁰+ 1)!!

r2l⁰+ 1

4π qR⁰L⁰. (7.6)

In der ASA wird nur der sph¨arisch symmetrische Teil (l = 0) des Potentials benutzt, welchen man aus (7.6) durch sph¨arische Mittelung erh¨alt:

V¯_M,R(r) = 1 4π

Z

R

V_M,R(r)d²ˆr= Z

R

VM,R(r)Y00(ˆr)

√4π d²ˆr=

= X

R⁰L⁰

0M_Rs,R⁰_L⁰q_R⁰_L⁰ =const. (7.7)

Hierbei wurde folgende Definition der Madelung-Matrix benutzt M_Rs,R⁰_L⁰ = 2

r 4π

2l⁰+ 1Y_L^∗0(R\−R⁰)|R−R⁰|^−l0−1. (7.8) Der sph¨arisch symmetrische Teil des Madelung-Potentials an R wird also von Po-tentialen erzeugt, die sich ausL-Multipolen konstituieren, die sich an den Positionen R⁰ ungleichRbefinden. Dies ist eine Konstante f¨ur eine gegebene Kugel anR. Mit Gleichung (7.7) ist das Problem in einen strukturellen Teil, repr¨asentiert durch die Madelung-Matrix und einen elektronischen Anteil, charakterisiert durch Multipole, zerlegt worden.

Den Madelung-Anteil zur Gesamtenergie des Kristalls erh¨alt man durch Summie-ren aller Wechselwirkungen der Dichteverteilungen (Elektronendichte plus Kerne) und Madelung-Potentiale (7.6) an den OrtenR.

E^{M ad} = 1 2

X

R

Z

(R)

[ρR(r)−ZRδ(r)]VM,R(r)d³r

= 1

2 X

RR⁰

0X

L,L⁰

q_RL^∗ M_RL,R⁰_L⁰q_R⁰_L⁰ (7.9) wobei die volle Madelung-Matrix durch

M_RL,R⁰_L⁰ =X

L⁰⁰

8π(−1)^l(2l⁰⁰−1)!!C_LL^∗ 0L⁰⁰

(2l−1)!!(2l⁰−1)!!p

(2l+ 1)(2l⁰+ 1)

Y_L^∗00(R\−R⁰)

|R−R⁰|^l00+1. (7.10) gegeben ist. Hierbei sind die Summen auf die Terme mitl⁰⁰=l+l⁰ beschr¨ankt. F¨ur l = 0 f¨uhrt dies auf den Spezialfall (7.8). Gleichung (7.10) beschreibt die Wechsel-wirkung zweier (Einheits-)Multipole.

Repr¨asentiert man einen Atomplatz durch R=B+T(Bsei der Vektor der zu einem Basisatom in der Elementarzelle am Ursprung zeigt) so kann Gleichung (7.9) weiter unter Benutzung der Regeln

X

RR⁰

0{. . .} −→ X

BB⁰TT⁰

(1−δ_BB⁰δ_TT⁰){. . .}

transformiert werden. Das ergibt:

E^{M ad} = 1 2

X

BB⁰TT⁰

X

LL⁰

q_(B+T)L^∗ M_(B+T)L,(B⁰_+T⁰_)L⁰q_(B⁰_+T⁰_)L⁰(1−δ_BB⁰δ_TT⁰)

= 1 2

X

BB⁰TT⁰

X

LL⁰

q_BL^∗ M_BL,(B⁰_+T⁰_−T)⁰_L⁰q_B⁰_L⁰(1−δ_BB⁰δ_TT⁰)

= N

2 X

BB⁰

X

LL⁰

q^∗_BL

"

X

T

M_BL,(B⁰_+T)⁰_L⁰(1−δ_BB⁰δ_TT⁰)

# q_B⁰_L⁰

Definiert man jetzt

MfBL,B⁰L⁰ =X

T

M_BL,(B⁰_+T)⁰_L⁰(1−δ_BB⁰δ_TT⁰), (7.11) was die Fouriertransformierte vonMRL,R⁰L⁰ am Punk k= 0 = Γ ist, so erh¨alt man f¨ur die Madelungenergie pro Einheitszelle des Kristalls

E^{M ad} N = 1

2 X

BB⁰

X

LL⁰

q_BL^∗ Mf_BL,B⁰_L⁰q_B⁰_L⁰.

Da die auftauchenden Summen f¨ur niedrigel, l⁰-Werte (l+l⁰ <4) nur bedingt konver-gent sind, muß die Madelung-Matrix (7.11) unter Verwendung von Ewaldverfahren bestimmt werden, was im Abschnitt 7.2.1 ausgef¨uhrt wird.

In der numerischen Umsetzung ist die Berechnung von Ewaldsummen in Glei-chung (7.11) sehr zeitaufwendig. Außerdem skaliert die Gr¨oße der Matrix (7.11) wie (nal²)², wobei na die Anzahl der Atome der Einheitszelle und l die Maximalzahl der Drehimpulsentwicklung in den Kugeln ist. Beide Nachteile lassen sich umgehen, indem man die Ewaldsummen (7.12) vor Beginn der selbstkonsistenten Schleife be-rechnet und die Elemente der vollen Madelung-Matrix (7.11) in jeder Iteration der Rechnung neu aufstellt.

S_B,B⁰_,L⁰⁰ =X

T

Y_L^∗00(B−\B⁰−T)

|B−B⁰−T|^l00+1 (7.12) 7.2.1 Elektrostatische Felder - Monopolbeitr¨age

Beim Versuch, das elektrostatische Feld am Punktrder Einheitszelle zu berechnen, welches durch Monopole der Ladungq, die an Gitterpl¨atzen Tsitzen, erzeugt wird, ergibt sich ein Problem, da Summen der Form

X

T

q

|r−T| (7.13)

divergent sind (das kann man sich am leichtesten mit der Divergenz der harmoni-schen Reihe P

n1

n klar machen). Der Trick besteht darin, einen divergenten Term in (7.13) abzuspalten, welcher einen anderen Term ausl¨oscht, der auftaucht, wenn man Ladungsneutralit¨at fordert.

Summen der Art (7.13) lassen sich f¨ur translationsinvariante Systeme mit einer Methode auswerten, die auf Ewald [172] zur¨uckgeht und auf einer Fourieranalyse basiert (siehe auch Ham und Segall[173]).

Damit ergibt sich f¨ur das elektrische Feld an r, welches von Punktladungen an Gitterpl¨atzen Therr¨uhrt als:

ϕ₃(r) =X

T

erf c(η|r−T|)

|r−T| − π

Ωη² + 4π Ω

X

G6=0

e^−G

2 4η−iG·r

G² +A(G→0) (7.14) Dieses hat Singularit¨aten anr=T. Der regul¨are Teil anr=Tist durch

ϕ₃(r=T) =−2η

√π +X

T6=0

erf c(η|T|)

|T| − π Ωη² +4π

Ω X

G6=0

e^−G

2 4η

G² +A(G→0) (7.15) gegeben.

Der TermAsoll den Grenzwert der Summe im reziproken Raum bezeichnen, der f¨urGgegen Null auftaucht. Dieser Term r¨uhrt vom langreichweitigen Charakter des Coulombpotentials her (pr¨aziser formuliert ist es dieG= 0 Singularit¨at der Fourier-transformierten des Coulombpotentials). Man sieht leicht, daß dieser Term f¨urjedes festerwie 1/G²gegen unendlich geht. Berechnet man jetzt das elektrostatische Feld am Punktr, indem man alle Beitr¨age der Gitterpl¨atze ¨uberlagert

φ(r) =X

B

ϕ(r−B)q_Bs, (7.16)

ist sofort klar, daß sich diese Terme wegen der geltenden Ladungsneutralit¨atP

BqB = 0 ausl¨oschen. Folglich k¨onnen diese Beitr¨age in den Summen (7.14,7.15) weggelassen werden.

7.2.2 Multipolmomente in der LMTO-Method

Die Multipolmomente (7.4) werden durch die Elektronendichte, gegeben durch ρ(r) =

occ

X

kj

|ψ^kj(r)|² =X

kj

|ψ^kj(r)|²θ(E_F −E_j(k))

ausgedr¨uckt, wobeiψ_kj(r) die Kristallwellenfunktion ist. Von nun an verwenden wir die Notation von Skriver.[41] In diesem Abschnitt werden reelle Linearkombinatio-nen der Kugelfl¨achenfunktionen benutzt (siehe Abschnitt 7.3). Zuerst wird der Fall behandelt, daß ein Atom pro Elementarzelle vorhanden ist, welches am Ursprung sitzt (sp¨ater erfolgt eine Erweiterung). Die nach Drehimpulsen entwickelte Wellen-funktion ergibt sich damit zu:

ψ^kj(r) =X

L

ψ_L^kj(r) =X

L

[A^kj_Lφ_νL(r) +B_L^kjφ˙_νL(r)] (7.17) mitφ_νL(r) =Y_L(ˆr)ϕ_νl(r) und ˙φ_νL(r) =Y_L(ˆr) ˙ϕ_νl(r). Der Subskriptν bezeichnet die linearisierte Partialwelleφ_νL=φ_L(E,r)|_E=Eν. F¨ur das Betragsquadrat der Wellen-funktion erh¨alt man damit

|ψ^kj(r)|²= X

L⁰L⁰⁰

{ A^kj_L0^∗A^kj_L00φ^∗_νL⁰φνL⁰⁰+A^kj_L0^∗B_L^kj00φ^∗_νL⁰φ˙νL⁰⁰+ +B_L^kj0^∗A^kj_L00φ˙^∗_νL⁰φνL⁰⁰+B_L^kj0^∗B_L^kj00φ˙^∗_νL⁰φ˙νL⁰⁰},

womit sich die Multipolmomente (7.4) QL=

r 4π 2l+ 1

X

kj

X

L⁰L⁰⁰

θ(EF −Ej(k)) Z

sphere

r^lYL(ˆr){. . .}d³r ergeben. Nutzt man jetzt

• R

(R)d³r{. . .}=R

(R)drr²R

Ωd²ˆr{. . .} und

• die Definition der Gaunt-Koeffizienten (hier werden reelle Linearkombinatio-nen benutzt):C_LL⁰_L⁰⁰ =R

ΩY_L(ˆr)Y_L⁰(ˆr)Y_L⁰⁰(ˆr)d²ˆr, erh¨alt man f¨ur die Multipolmomente

QL = X

kj

X

L⁰L⁰⁰

CLL⁰L⁰⁰

Z

sphere

{A^kj_L0^∗A^kj_L00ϕ^∗_νl⁰ϕνl⁰⁰+A^kj_L0^∗B_L^kj00ϕ^∗_νl⁰ϕ˙νl⁰⁰+ +B_L^kj0^∗A^kj_L00ϕ˙^∗_νl⁰ϕ_νl⁰⁰+B_L^kj0^∗B_L^kj00ϕ˙^∗_νl⁰ϕ˙_νl⁰⁰}r^l+2θ(EF −Ej(k))dr Erweiterung des Formalismus

Hier soll der bis oben aufgef¨uhrte Formalismus weiterentwickelt werden, so daß meh-rere Atome pro Elementarzelle ber¨ucksichtigt werden und Symmetrieoperationen zur Reduktion von Integrationen auf den irreduziblen Teil der Brillouinzone einbezogen werden.

Der Ausdruck f¨ur die Wellenfunktion (7.17) wird mit einem weiteren Index ver-sehen, welcher die Basisatome numeriert:

ψ^kj(r) =X

LB

[A^kj_LBφ_νBL(r_B) +B_LB^kj φ˙_νBL(r_B)]Θ(s_B−r_B) (7.18) Dabei istr_B =r−Bein Vektor in der Kugel anBunds_B sei der Radius der Kugel.

Seig={p|v}eine Symmetrieoperation der Raumgruppe G des Kristalls, wobei p eine Punktgruppenoperation des Gitters {T} bezeichne. v sei ein Vektor, der nicht in der Menge {T} enthalten sei, jedoch ergebe ein ganzzahliges Vielfaches von v einen Vektor des Gitters (v = 0 f¨ur symmorphe Raumgruppen). Sei eine Brillouinzone (BZ) gegeben; dann existiert eine Untermenge von Punkten Ω_IBZ so, daß P

p{p|0}Ω_IBZ gleich der Punktmenge der ganzen BZ ist. Dies nutzend kann man Summen ¨uber die volle BZ mit folgender Regel umschreiben:

BZ

X

k

f(k) −→ 1 N^so

X

p IBZ

X

k

f(pk)

wobei N^so die Anzahl der Symmetrieoperationen der Punktgruppe ist. Dann kann dasL-te Moment an B folgendermaßen geschrieben werden

Q_BL= 1 N^so

X

p IBZ

X

kj

Z

(B)

dr_B Z

dˆr_Br^l+2_B Y_L(ˆr_B)|ψ^(pk)j(r)|².

Wendet man die Raumgruppenoperation {p|v} auf eine Bloch-Wellenfunktion an, so erh¨alt man (siehe z. B. [174])

{p|v}ψ^k(r) =ψ^pk(r) =ψ^k(p⁻¹r−v). (7.19)

In Gleichung (7.18) wurde die Kristall-Wellenfunktion in atomzentrierte Wellenfunk-tionen entwickelt, die am Kugelrand abgeschnitten werden. F¨ur eine gegebene Kugel wurde der r-Vektor als r = B+rB geschrieben. Damit und mit Gleichung (7.19) erh¨alt man

ψ^pk(B+r_B) =ψ^k(p⁻¹r_B−p⁻¹B−v),

was nichts anderes als der Wert der Blochwelle mit dem Index k in einer Kugel zentriert anB⁰ =p⁻¹B+v am rotierten Punkt rB ist. Folglich ergibt sich f¨ur das Integral ¨uber die Kugel

Z

(B)

dr_B Z

dˆr_Br_B^l+2Y_L(ˆr_B)X

L⁰L⁰⁰

[A^jk_L0BA^jk_L00^∗Bφ_ν(p⁻¹_B−v)L⁰(p⁻¹r_B)φ^∗_ν(p−1B−v)L⁰⁰(p⁻¹r_B)+. . .], wobei die Auslassungspunkte die anderen Terme des Produkts der Blochwellen an-deuten. Nutzt man jetzt folgende Relation f¨ur das Integrieren ¨uber den Winkelanteil

Z

Ω

dˆrY_L(ˆr)Y_L⁰(p⁻¹ˆr)Y_L⁰⁰(p⁻¹ˆr) = Z

Ω

dˆrY_L(pˆr)Y_L⁰(ˆr)Y_L⁰⁰(ˆr) =X

L⁰⁰⁰

D_LL^p 000C_L⁰⁰⁰_L⁰_L⁰⁰,

kann dasL-te Moment anB als Q_BL=X

p

X

L⁰

D^p_LL0Q^irr_(p−1B−v)L⁰

ausgedr¨uckt werden. Hierbei istQ^irr_BLdas Multipolmoment, welches aus der Integra-tion ¨uber die irreduzible BZ erhalten wurde. Die Wigner-Funktionen werden in 7.3 n¨aher erl¨autert.

Im Dokument Beschreibung der elektronischen Struktur von Übergangsmetalloxiden mittels selbstwechselwirkungskorrigierter Dichtefunktionaltheorie: Volumenkristalle, Oberflächen und Punktdefekte (Seite 96-102)