Dft

Chapitre 1

Formalisme de la DFT 1. Introdu tion

Dans la bran he de la matière ondensée, les modèles servants à la des ription des problèmes traitants des systèmes atomiques, molé ulaires et solides sont développés de manière à permettre des al uls ave un nombre réduit (limité) d'atomes non-équivalents tout en essayant d'introduire le plus grand nombre possible d'intera tions entre les diérents onstituants (i.e. les éle trons et les noyaux) des systèmes à étudier. C'est e dernier ompromis qui n'est pas toujours aisé à satisfaire du fait de sa omplète dépendan e en la nature des systèmes. Les méthodes mises en ÷uvre pour la réalisation de e type de

al uls sont en général réparties en trois grandes lasses :  La première lasse onstituée par les méthodes de premier-prin ipe {ab-initio : LMTO, LAPW }1 développées par Andersen en 1975 et les méthodes du pseudopotentiel (PP)2 fon tionnant sur la base de l'idée mise à jour par Phillips & Kleinman en 1959 [1℄,[2℄,[3℄,[4℄,[5℄,[6℄,[7℄. Le point ommun partagé par toutes es méthodes est la possibilité d'un démarrage de leurs al uls du zéro sans pour autant faire appel à des paramètres d'ajustements expérimentaux (largeur de bande, moment magnétique, paramètre d'é hange ...) à l'opposé des méthodes semi-empiriques.  La se ond lasse est la méthode semi-empirique des liaisons fortes (TB )3 utilisant des orbitales onstruites à partir d'une ombinaison linéaire d'orbitales atomiques (LCAO )4 et permettant une meilleure des ription des éle trons fortement liés, du type des éle trons de valen e d des séries de transition. Des paramètres ajustables pour une meilleure adaptation des résultats al ulés ave les données de l'expérimentalmettant en ÷uvre des orbitales onstruites sur la base de la LCAO est l'une des méthodes semi-empiriques usuelles pour le traitement du magnétisme des systèmes métalliques à base des éléments des séries de transition [8℄,[9℄,[10℄,[11℄. 1 2 3 4

LMTO : Linear Mun-tin Orbitals, LAPW : Augmented Plane Waves Methods. PP : Pseudopotential Methods. TB : Tight-Biding Method LCAO : Linear Combination of Atomi Orbitals (ele troni eigen-fun tions of Hydrogen atom).

2

Chapitre 1. Formalisme de la DFT  La troisième lasse onstituée par les méthodes APW & KKR 5 développées respe tivement par Slater en 1937 et Korringa-Kohn-Rostoker en 1954 [12℄,[13℄,[14℄,[15℄.

Les al uls de premier-prin ipe réalisés dans le ontexte du présent travail ont pour prin ipal obje tif une des ription des propriétés stu turales, des stru tures éle troniques et du magnétisme des systèmes métalliques en multi ou hes et en alliages binaires de Cobalt et de Ni kel sur la base d'une appli ation de l'une des méthodes ab-initio de la première

lasse, à savoir la méthode des orbitales MTO's linéairisées développée par Andersen et al. et onstituant une base minimale susante pour e type d'études [1℄,[16℄. En fait, e

hoix de la base des MTO's est essentiellement soutenu par les nombreux avantages la

ara térisant, ette dernière onstituant à la fois :  Une base monoéle tronique entrée sur ha un des sites atomiques du réseau ristallin.  Une base se servant d'orbitales de ourtes portées (distan es).  Une base appropriée à des appli ations sur des massifs ainsi que sur des systèmes à faibles dimensions (surfa es et interfa es).  Une base susamment petite né essitant un nombre réduit d'orbitales par atome, ave neuf orbitales pour les as sp-d et seize orbitales pour les as spd-f.  Une base permettant la dépendan e en énergie né essaire pour la onstru tion des fon tions de Green, donnants lieu aux densités des états quantiques (partielles PDOS et totales TDOS ).  Une base omplète mettant en ÷uvre des potentiels mun-tins de symétries sphériques.  Une base appliquant des développements à un- entre de leurs orbitales en fon tion des ondes partielles radiales (PRW )6 , des harmoniques sphériques et des onstantes de stru tures [1℄,[17℄.

2. Théorie de la Fon tionnelle de la densité 2.1. Positionnement du problème

La résolution des problèmes quantiques onstitués par les systèmes physiques de diérentes natures (atomiques, molé ulaires ou solides) et traités dans la bran he de la matière

ondensée est réalisée sur la base d'une appli ation des on epts de la mé anique quantique et d'un nombre d'approximations simpli atri es mais rigoureuses (i.e. basées sur des réalités physiques) de manière à permettre une des ription mi ros opique des propriétés et des phénomènes de l'état fondamental des systèmes en question, à travers leurs mouvements éle troniques. En fait, le traitement exa t des systèmes quantiques réels est le plus souvent

onfronté au problème du nombre important de variables (éle troniques et nu léaires) les

ara térisant et rendant leur résolution exa te une tâ he quasiment impossible à réaliser. L'une des te hniques adoptées pour surmonter e type de problèmes est le re ours à des approximations des systèmes réels omplexes par des systèmes voisins de moindre omplexité, des orre tions sont par la suite apportées aux solutions approximatives résultantes de manière à mieux les on order ave les données expérimentales. La di ulté à soulever est liée au bon hoix des approximations à mettre en ÷uvre ainsi qu'à elui de la base des fon tions d'essai à appliquer, des hoix entièrement dépendants de la nature des problèmes à résoudre. 5 6

APW : Augmented Plane Waves, KKR : Korringa-Kohn-Rostoker Method. PRW : Partial Radial Wave

3 La théorie de la fon tionnelle de la densité (DFT )7 est utilisée omme un outil mathématique pour la résolution des problèmes à plusieurs orps, du type de eux ren ontrés dans les études des systèmes polyéle troniques orrélés en général et des solides ristallins en parti ulier [18℄,[19℄,[20℄. A l'opposé de la théorie de Hartree-Fo k (HFA)8 dé rivant des éle trons individuels en intera tion ave le reste des éle trons et des noyaux du milieu, la théorie de la DFT est basée sur une des ription du système en entier de manière à onstituer une meilleure approximation pour la résolution de e type de problèmes polyéle troniques [21℄,[22℄. Le formalisme de la DFT est une théorie développée sur la base des deux théorèmes de Hohenberg-Kohn et allant au-delà de la HFA, à travers une prise en ompte des eets de orrélation dans ses études des propriétés physiques de l'état fondamental des systèmes polyéle troniques orrélés. Les orre tions ainsi introduites en termes de ontributions d'é hange- orrélation (XC )9 ont révélé une meilleure pré ision de

al uls des énergies des systèmes polyéle troniques. En général, le traitement quantique d'un solide ristallin onstitué d'un nombre N de noyaux relativement lourds (de masses M et de harges +Ze ) et d'un nombre n d'éle trons (de masses m et de harges -e ) gravitant autour de es derniers est basée sur une résolution d'un problèmeo onséquent à n+N parti ules (i.e. à 3n+3N variables spatiales : n ~ 1, R ~ 2 ...R ~ N ) en intera tions mutuelles. L'hamiltonien total orrespondant ~r1 , ~r2 ...~rn , R à un tel système et dé rivant l'ensemble des intera tions s'y produisant est exprimé sous sa forme exa te suivante : ˆ = Tˆelec + Tˆnucl + U ˆelec−elec + U ˆelec−nucl + U ˆnucl−nucl H

(1.1)

Dans l'Eq. (1), le premier terme désignant l'opérateur énergie inétique des noyaux, le se ond elui des éle trons et les trois derniers dé rivants les intera tions respe tives (éle tron-éle tron), (noyau-noyau) et (éle tron-noyau), tels que : Telec =

X i

Uelec−elec =

−

X ~2 ∆2 ~2 ∆2i λ , Tnucl = − mi 2 Mλ 2 λ

e2 X −Zλ e2 X −1 , Uelec−nucl = 4πεo |~ri − ~rj | 4πεo ~ r − R ~ i λ i6=j λ, i Unucl−nucl =

e2 X +Zλ .Zµ ~ 4πεo ~ λ6=µ Rλ − Rµ

~ λ désignants respe tivement les ve teurs positions éle troniques et nu léaires servant ~ri et R à la lo alisation de ha un des éle trons i du système et de ha un de ses noyaux λ,

entrés sur ses sites atomiques. Les indi es i ={1...n } et λ ={1...N } sont ainsi appliqués de manière à distinguer les grandeurs éle troniques des grandeurs nu léaires. 7 8 9

DFT : Density Fun tional Theory HFA : Hartree-Fo k Approximation XC : Ex hange-Correlation Ee ts

4

Chapitre 1. Formalisme de la DFT

Compte tenu du nombre élevé 3 (n+N ) de degrés de liberté (variables spatiales) et d'intera tions mises en jeu dans e type problèmes (tridimensionnels et à plusieurs orps), leur traitement exa t s'est révélé une tâ he impossible à réaliser. L'une des solutions pré onisée est le re ours à des approximations approprées et simpli atri es, du type

elles de la densité lo ale et de la densité lo ale de spins {LDA & LSDA}10 appliquées dans le formalisme de la DFT. L'idée fondamentale visée à travers ette appli ation du formalisme de la DFT-LDA (LSDA) est le passage du al ul de fon tions d'onde en un

al ul de densités de harge sur la base du théorème de Hohenberg-Kohn développé en 1964 [23℄. Les deux niveaux d'approximations appliqués par la théorie de la DFT dans ses traitements des problèmes à plusieurs orps sont l'approximation de Born-Oppenheïmer (BOA) onstituant le premier niveau et l'approximation de la LDA (LSDA) omme se ond niveau.

2.2. Approximation de Born-Oppenheimer

L'approximation de Born-Oppenheimer (BOA)11 [24℄ est appliqué omme un premier niveau d'appro he des problèmes à plusieurs orps en mettant en éviden e le grand é art entre la masse des éle trons du système (plus légers, don de plus grande mobilité) et elle des noyaux (relativement plus lourds : M ≈ 1800×m ). Autrement dit, ette approximation est basée sur l'idée onsidérant les noyaux omme animés de mouvements susamment longs, relativement à eux des éle trons, de manière à les négliger sans grande erreur, en ~λ = R ~ oλ={1...N} sur elles des noeuds des réseaux les geant les noeuds de leurs sites R

ristallins. Cette approximation est ainsi appliquée de manière à ramener le problème à plusieurs orps (i.e. d'éle trons et de noyaux) de départ en un problème polyéle troniques (à éle trons seuls) à travers une disso iation de l'eet des noyaux, exprimé sous la forme d'une intera tion externe, de elui du nuage éle tronique. En onséquen e, ave un se ond terme nul et un dernier onstant, l'Eq.(1) est reé rite sous une forme plus réduite ave

omme seuls termes eux de l'énergie inétique éle tronique Tˆelec et de l'intera tion (éle tron-éle tron) Uelec−elec et de l'intera tion externe (éle tron-noyau) Uˆext : ˆ = Tˆelec + U ˆelec−elec + U ˆext H

(1.2)

L'expression résultante est elle d'un hamiltonien servant au al ul des systèmes à plusieurs éle trons (polyéle troniques) en intera tions mutuelles et en dépla ement dans ˆelec−elec onstituant un terme le potentiel externe nu léaire Uˆext . Le terme Tˆelec + U universel du fait de son entière indépendan e de la nature et de la stru ture des sysˆext . tèmes traités, des informations entièrement ontenues dans le terme externe des noyauxU En adoptant le système des unités atomiques ( f. ŸAppend.), adoptant le rayon de Bohr

omme unité de base des longueurs et le Rydberg (ou le Hartree ) omme unités des énergies, les termes de l'énergie inétique, de l'intera tion (éle tron-éle tron) et du potentiel externe des noyaux omposants l'Eq. (3) sont respe tivement exprimés sous la forme : Tˆelec = −

X ∆2

i

i

2

ˆelec−elec = , U

X i6=j

X 1 ˆext ≡ u ˆext (~ri ) , U |~ri − ~rj | i

LDA & LSDA} : Lo al and Lo al Spin Density Approximation. BOA : Born-Oppenheimer Approximation

10 { 11

(1.3)

5

ˆext représentant le potentiel externe engendré par les N noyaux du système Le potentiel U et ae tant ha un de ses éle trons. L'Eq. (3) est le plus souvent développée sous sa forme la plus usuelle :  X X  ∆i ˆ +U ˆelec−elec ˆelec−elec = h +u ˆext (~ri ) + U − H= (1.4) i 2 i i | {z }

(a)

u ˆext représentant le potentiel externe agissant sur ha un des éle trons du système traité {i=1...n } : X −Zλ u ˆext (~ri ) = ~ r − R ~ i λ λ

(b)

2.3. Théorie de la fon tionnelle de la densité

L'appli ation de l'approximation de la BOA a permis ertes de réduire la omplexité du problème de départ en le ramenant à un problème d'éle trons seuls mais sans pour autant le résoudre, l'Eq. (3) résultante restant toujours di ile à solutionner et né essitant des approximations supplémentaires. En général, le al ul des systèmes polyéle troniques résultants de l'approximation de a la BOA est réalisé sur la base d'une appli ation de l'approximation de la HFA ou du formalisme de la DFT. La première méthode orant une bonne des ription des as atomiques et molé ulaires s'est révélée inappropriée à l'étude des solides, ette dernière est ainsi altérée par la négligen e des eets de orrélations éle troniques mises en éviden e par la DFT. A l'opposé de la méthode de la HFA dé rivant l'intera tion de ha un des éle trons ave le hamp moyen du reste des éle trons et des noyaux du milieu, le formalisme de la DFT est basé sur une des ription du système tout entier12 . La théorie de la DFT, appropriée au traitement des solides, est développée par Hohenbergen 1964 sur la base de leurs théorèmes reposant sur l'idée d'une des ription du potentiel externe des noyaux u ˆext (~r) ≡ u ˆext [ρ] (et à travers lui l'énergie totale 13 et est déterminée à partir des solutions des équations Kohn-Sham propres à ha un des éle trons du système [23℄,[25℄,[26℄,[27℄. Kohn

2.3.1. Théorèmes de Hohenberg-Kohn En général, le traitement des problèmes des systémes quantiques polyéle troniques est basé une re her he des solutions (états) propres |Ψi et des énergies E (valeurs propres) orres-

MFA

12 : Mean Field Approximation 13 L'idée de base d'une é riture de l'énergie

totale du système omme une fon tionnelle de la densité de

harge est introduite pour la première fois par Thomas & Fermi en 1927, développée théoriquement par Hohenberg-Kohn en 1964 pour la formulation de leur théorie de la DFT et reprise par la suite par KohnSham en 1965 pour le développement d'une appro he de al ul des stru tures éle troniques omme une fon tionnelle unique de la densité de harge ρ(~r), une grandeur adoptée omme variable de base de la DFT

6

Chapitre 1. Formalisme de la DFT

pondantes à travers une résolution de l'équation stationnaire de du temps (TISE ) 14 .

S hrödinger

indépendante

ˆ |Ψi = E. |Ψi H h~r1 ...~rn |Ψi ≡ Ψ (~r1 ...~rn ) représentant la fon tion d'onde totale de système (à dé rite dans une représentation {|~r1 ...~rn i}.

(1.5)

n

éle trons)

Cette fon tion est développée dans l'approximation de la HFA, sous la forme d'un déterminent de Slater des orbitales atomiques AO's (solutions propres monoéle troniques ↑(↓) ↑(↓) pour l'atome d'Hydrogène ) et des fon tions de spins : ψi (~rj ) = φi (~rj ) × σj [21℄. ↑(↓) ↑(↓) ↑(↓) ψ 1 (~r1 ) ψ2 (~r1 ) . . . ψn (~r1 ) ↑(↓) ↑(↓) ↑(↓) 1 ψ1 (~r2 ) ψ2 (~r2 ) . . . ψn (~r2 ) HF A (~r1 ...~rn ) = √ Ψ (1.6) .. .. .. .. n! . . . . ↑(↓) ↑(↓) ψ (~rn ) ψ2 (~rn ) · · · ψn↑ (~rn ) 1 Pour e type de problèmes, le al ul de ette densité de harge a une importan e apitale du fait de son utilisation pour la détermination du nombre d'éle trons et du potentiel externe, deux grandeurs fondamentales né essaires pour la onnaissan e des propriétés mesurables de tout système éle tronique. Dans la base des fon tions d'onde |Ψi, la densité de harge ρ (~r) est développée sous la forme [21℄ : X Z ρ (~r) = hΨ | Ψi = (1.7) Ψ∗ (~r, ~r2 ...~rn ) .Ψ (~r, ~r2 ...~rn ) d~r2 ... d~rn σ1 σ2 ...σn

les o

upations éle troniques du système étant exprimées par : Z n = ρ (~r) d~r Le premier théorème développé par Hohenberg-Kohn est utilisé pour l'établissement d'une orrespondan e entre la densité de harge ρ d'un système polyéle tronique à l'état fondamental et le potentiel externe de ses noyaux vext [ρ] de telle sorte à pouvoir dé rire la ˆ omme une fon tionnelle unique de la densité éle tronique valeur de toute observable X exa te, propre à et état fondamental : ˆ (~r) |Ψi = X [ρ] hΨ| X ave

ρ ≡ ρ (~r)

14

TISE : Time Independent S hrödinger Equation

(1.8)

7 Cette possible utilisation de la densité de harge à la pla e de la fon tion d'onde est justiée par le fait que omme haque système polyéle tronique ne pouvant présenter qu'un seul et unique potentiel externe ( elui utilisé pour la onstru tion de son propre hamiltonien, Eq. (3)), une résolution de son équation de S hrödinger Eq. (5) ne onduisant qu'à une seule et unique fon tion d'onde et don qu'à une seule et unique densité éle tronique, ette dernière grandeur ontenant autant d'informations sur le système traité que sa fon tion d'onde. En fait, à haque potentiel externe va orrespondre une densité de harge unique et vi e-versa. Le se ond théorème de Hohenberg-Kohn est issu d'une appli ation du premier théorème ˆ ≡H ˆ ), 'est un théorème mettant en éviden e pour le as d'une observable hamiltonien (X une é riture de l'énergie totale de l'état fondamental du système polyéle tronique omme une fon tionnelle de sa densité de harge, E [ˆ uext [ρ]] ≡ E [ρ]. Une appli ation de es deux théorèmes de Hohenberg-Kohn aux systèmes polyéle troniques à l'état fondamental a donné lieu à des énergies totales des systèmes manifestant des minimums oïn idants ave les densités de harge propre à leur état fondamental. Autrement dit, parmi les nombreuses densités de harge ρ possibles, elle orrespondant à l'état fondamental est déterminée à partir du prin ipe variationnel de Rayleigh-Ritz sur la base d'une minimisation de l'énergie du système E [ρ], pour un potentiel externe uext [ρ] bien déni : . ˆ ˆelelc−elec |Ψi + hΨ| U ˆext |Ψi E [ρ] = hΨ H Ψi hΨ | Ψi = hΨ| Tˆelec + U | {z } Z = F [ρ] + u ˆext [ρ].ρ (~r) d~r | {z } =

F [ρ]

+

(1.9)

ˆext [ρ] U

(a) La fon tionnelle F [ρ] ne omportant au une information ni sur les noyaux ni sur leurs positions est ommune à l'ensemble des systèmes polyéle troniques, 'est une grandeur universelle entièrement indépendante du potentiel des noyaux et est exprimée sous la forme :

(b)

ˆelec−elec Ψi F [ρ] = hΨ Tˆelec + U

Dans le formalisme de la DFT, ette fon tionnelle universelle F [ρ] est déterminée à partir de la densité de harge de l'état fondamental al ulée à son tour par la méthode de Kohn-Sham à travers une résolution de leurs équations monoéle troniques. En n de

ompte, la onnaissan e de la densité de harge ρ du système à l'état fondamental et de la fon tionnelle orrespondante F [ρ] est susante pour une détermination des diérentes propriétés physiques (o

upations, moments magnétiques, énergies, DOS's, ompréssibilité, for es atomiques, ...) des milieux matériels traités.

8

Chapitre 1. Formalisme de la DFT

2.3.2. Méthode de Kohn-Sham L'opération mise en ÷uvre pour l'obtention d'une fon tionnelle universelle G [ρ] entièrement indépendante du potentiel des noyaux est basée sur une soustra tion de la partie éle trostatique (de Hartree ) de l'expression de la première fon tionnelle F [ρ]. La fon tionnelle G [ρ] ainsi obtenue est essentiellement onstituée d'une partie inétique et d'un é art d'énergie entre l'intera tion exa te Uelec−elec [ρ] et l'intera tion oulombienne UH [ρ] : Z G [ρ] = F [ρ] − vH (~r) .ρ (~r) d~r = F [ρ] − UH [ρ] (1.10) | {z } La forme de la fon tionnelle G [ρ] utilisée par Kohn-Sham est établie pour un système hypothétique d'éle trons sans-intera tion (S ) manifestant une même densité de harge ρ (~r) que elle du système réel (ave intera tion). En fait, e système hypothétique (S ) est introduit omme un système de référen e pour la onstru tion de des potentiels ee tifs ˆKS [ρ]. Kohn-Sham U Ave la prise en ompte de e on ept du système hypothétique (S ) de fon tionnelle de Hohenberg-Kohn est réé rite sous la forme : F [ρ] = Telec [ρ] + Uelec−elec [ρ]

Kohn-Sham,

la

(1.11)

= Telec [ρ] + Uelec−elec [ρ] + {TS [ρ] − TS [ρ]} = TS [ρ] + Uelec−elec [ρ] + Telec [ρ] − TS [ρ] | {z } = TS [ρ] + Uelec−elec [ρ] + EC [ρ] + {UH [ρ] − UH [ρ]} = TS [ρ] + UH [ρ] + EC [ρ] + Uelec−elec [ρ] − UH [ρ] | {z } = TS [ρ] + UH [ρ] + EC [ρ] + EX [ρ] | {z } = TS [ρ] + UH [ρ] + EXC [ρ]

donnant lieu ainsi à la fon tionnelle G [ρ] suivante : G [ρ] = F [ρ] − UH [ρ] = TS [ρ] + EXC [ρ] EXC [ρ] dé rivant le terme d'é hange- orrélation (XC ).

(1.12)

9 La fon tionnelle de l'énergie orrespondant au système d'éle trons sans-intera tion (S ) et en dépla ement dans le potentiel de Kohn-Sham vKS [ρ], est égale à elle du système réel (ave intera tion) et est exprimée omme : Z E [ρ] = TS [ρ] + vKS (~r) .ρ (~r) d~r (1.13) | {z } ˆKS [ρ] = TS [ρ] + U

Compte tenu des termes (éle trostatique, externe et de l'XC ) né essaires à la onstru tion des potentiels ee tifs de Kohn-Sham vKS [ρ], l'énergie totale de l'Eq. (13) est développée

omme elle d'un gaz d'éle trons sans-intera tion (S ) et soumis à l'a tion du potentiel externe (des noyaux) et de l'XC : Z Z E [ρ] = TS [ρ] + vH (~r) .ρ (~r) d~r + vext (~r) .ρ (~r) d~r +EXC [ρ] (1.14) | {z } | {z } = TS [ρ] + UH [ρ] + Uext [ρ] + EXC [ρ]

En se basant, d'une part, sur l'idée prin ipale du théorème de Kohn-Sham onsidérant la densité de harge ρ (~r) de leur système hypothétique sans-intera tion (S ) omme la densité exa te du système réel à l'état fondamental et, d'autre part, sur l'Eq. (14) révélant l'importan e d'une onnaissan e de la densité de harge ρ (~r) pour la détermination du reste des grandeurs énergétiques (potentiels et énergies) du système, le al ul de ette densité est réalisé à partir des solutions propres ψi (~r) des équations de Kohn-Sham [27℄. Ces dernières, développées et publiées en 1965, sont l'analogue des équations de S hrödinger monoéle troniques et sont utilisées pour le al ul de la densité de harge des systèmes polyéle troniques à l'état fondamental, pour ha un des états éle troniques {|jki} propre à ha un des éle trons i = {1, 2...n} [28℄ :   ∆ Hjk ψjk (~r) = − + vKS (~r) ψjk (~r) (1.15) 2 = Ejk .ψjk (~r)

(a) ave

ρ (~r) =

occ X j,k

(b)

2

|ψjk (~r)| =

occ occ X core 2 X val 2 φjk (~r) + ψjk (~r) j,k

j,k

Les formes non-relativiste et relativiste des hamiltoniens de Kohn-Sham Hi entrants dans le système d'équations pré édent (Eq. (15.a)) sont respe tivement données par : Hi = −

∆i + vKS (~r) 2

(1.16)

10

Chapitre 1. Formalisme de la DFT ∆ 1 Hi = − i + vˆKS (~r) + 2 2 c

(a)



∂v (~r) ∂ (E − vKS (~r)) + ∂~r ∂r 2



+ λ~ℓ.~s

En parti ulier, la forme s alaire-relativiste de l'hamiltonien monoéle tronique radial appliquée par la méthode de la LMTO est elle développée par Harmon-Koelling (1977 ) [29℄ : ˜ = H

1 ∂2 ∂v (r) ∂ 1 ℓ (ℓ + 1) 1 − + {v (r) − E} r + 2m ˜ (r) ∂r2 2m ˜ (r) .c2 ∂r ∂r r r2

(a) ave

m ˜ (r) = m +

1 c2 {E

− v (r)}SI = 1 +

1 c2 {E

− v (r)}a.u. (b)

c désignant la vitesse de la lumière dans le vide, ~ℓ le moment orbital et ~s éle tronique.

le spin

La résolution du système à n équations de Kohn-Sham (KS ) résultant est basé sur une détermination des n solutions propres individuelles ψi (~r) et des valeurs propres orrespondantes Ei (énergies), des grandeurs al ulées pour ha un des états éle troniques {|jki} (d'indi e de bande j et de ve teur d'onde éle tronique ~k ). Les ompli ations introduites par le terme d'intera tion (éle tron-éle tron) lors résolution de l'Eq. (5) sont omplétement surmontées ave la résolution de e système d'équations de KS, Eq. (15.a). La densité de

harge résultant de l'Eq. (15.b) pour le as du système hypothétique (S ) est assimilée à

elle du système polyéle tronique réel (ave intera tion) et est utilisée omme variable de base pour le al ul du potentiel ee tif de Kohn-Sham vKS [ρ] et du reste des grandeurs physiques intéressantes {o

upations éle troniques, spe tres des énergies, DOS's, moments magnétiques ...}.

2.3.3. Potentiel ee tif de Kohn-Sham L'expression du potentiel ee tif de Kohn-Sham vKS [ρ] est issue d'une minimisation, d'une part, de l'énergie totale de l'Eq. (14) donnant lieu à la densité de harge du système polyéle tronique réel à l'état fondamental : dTS [ρ] + vH (~r) + vext (~r) + vXC [ρ] = 0 dρ

(1.17)

ave

vXC [ρ] =

et, d'autre part, de elle du système

dEXC [ρ] dρ

hypothétique

(S ) :

dTS [ρ] + vKS (~r) = 0 dρ

(1.18)

la fon tionnelle du potentiel ee tif de Kohn-Sham déduite d'une identi ation de es deux dernières Eq's. (17) et (18) est développée sous la forme : vKS (~r) ≡ vKS [ρ] = vH (~r) + vext [ρ] + vXC [ρ]

(1.19)

11

Pour le as spin-polarisé, ave l'introdu tion du spin éle tronique, σ =↑ (↓) l'Eq. (19) de

e potentiel est réé rite : ↑(↓)

↑(↓)

vKS [ρ] = vH [ρ] + vext [ρ] + vXC [ρ]

(1.20)

les derniers termes de l'XC sont al ulés sur la base d'une appli ation des approximations de la LDA et de la LSDA, adaptées respe tivement pour les deux as non-polarisé et spin-polarisé sont utilisées omme premier niveau de al ul. Pour aller au-delà de es approximations de la densité lo ale, l'une des nombreuses versions de GGA's est mise en ÷uvre omme se ond niveau d'approximation du système traité. Le présent travail est basé sur une appli ation de la paramétrisation de von Barth-Hedin 15 pour l'évaluation des ↑(↓)

ontributions de l'XC au potentiel de Kohn-Sham vXC [ρ]. La paramétrisation de von Barth-Hedin (vBH72 ) est une généralisation au as spin-polarisé de elle de Hedin-Lundqvist (HL71 ) et est bien adaptée aux al uls des séries de métaux de transition non-magnétiques et magnétiques [30℄. Une fois évalués, les trois potentiels {d'XC vXC [ρ], éle trostatique vH (~r) 16 et externe vext [ρ] } sont sommés de manière à donner lieu au potentiels ee tifs (par spin) ↑(↓) vKS (~r) né essaires à la onstru tion des hamiltoniens de Kohn-Sham (par spin) rentrants dans les Eqs. (15.a) et dont la résolution est réalisée pour haque spin éle tronique σ =↑ (↓) :   ∆ ↑(↓) ↑(↓) ↑(↓) ↑(↓) Hjk ψjk (~r) = − + vKS (~r) ψjk (~r) (1.21) 2 ↑(↓)

↑(↓)

↑(↓)

= Ejk .ψjk (~r) ave

ρ↑(↓) (~r) =

occ X ↑(↓) 2 ψjk (~r) jk

↑(↓)

les densités de spins ρi (~r) issues d'une résolution de ha une de es Eqs. (21) sur la ↑(↓) base d'une re her he de leurs fon tions propres de spins ψi (~r) et de leurs énergies ↑(↓) sont utilisés pour le al ul de la densité de harge : individuelles Ei ρ (~r) = ρ↑ (~r) + ρ↓ (~r)

15 Suivant

la nature du milieu matérial à traiter et de la propriété à re her her, une multitude de paramétrisations est mise à jour pour une évaluation des ontributions d'é hange et de orrélation au potentiel de Kohn-Sham : HL71 [30℄, vBH72 [31℄, CA80 [32℄, VWN80 [33℄, PZ81 [34℄, LM83-86 [35℄, PW86-91 [36℄, [37℄, B88 [38℄, LYP [39℄, PBE96 [40℄, EV93 [41℄ 16 Le

potentiel atomique de Hartree est déterminé à partir d'une résoluation de l'équation de Poisson : −∆vH (~ r) = 4π ρ (~ r)

12

Chapitre 1. Formalisme de la DFT

Fig. 1.1  Pro essus itératif standard appliqué pour la résolution des équations monoéle troniques (KS ) de Kohn-Sham.

L'interdépendan e potentiel-densité (interprétée en termes, d'une part, de la dépendan e ↑(↓) des potentiels éle trostatique et de l'XC de la densité de harge ρi (~r) et, d'autre part, de la dépendan e de ette densité en es potentiels) est une ara téristique d'un problème self- onsistant (auto- ohérent), résolu par une pro édure itérative. Cette pro édure itérative (par spin) est initiée ave un potentiel externevext de forme bien dénie et une densité de σ(1) départ ρin servant de grandeur d'entrée à la première itération. Ces deux grandeurs sont σ(1) (1) σ(1) utilisées pour la onstru tion des potentiel ee tifs vKS (à partir de vext , vH et vXC ) et σ(1) σ(1) des hamiltoniens orrespondants Hi . La densité de sortie ρout est ainsi al ulée à partir

13 σ(1)

des solutions propres par spin ψi issues d'une résolution de l'Eq. (21). De la même σ(2) manière, la deuxième itération est initiée ave une densité ρin (résultant du mélange des σ(1) σ(1) densités d'entrée ρin et de sortie ρout de l'itération pré édente) servant au al ul du σ(2) σ(2) nouveau potentiel ee tif vKS et à la détermination de l'hamiltonien Hi . σ(2)

Les solutions propres ψi issues d'une résolution de l'Eq. (21) sont ainsi utilisées pour σ(2) la détermination à la densité de sortie ρout de ette se onde itération...et ainsi de suite. Pareille pro édure est répétée jusqu'à la satisfa tion de la ondition de onvergen e, une

ondition introduite à travers un test de omparaison des valeurs des densités (de sortie) des itérations su

essives {n & n + 1} et mettant n au pro essus itératif en ours, Fig. 1 : σ(n+1) σ(n) − ρout 6 10−6 ρout

La densité de harge résultant des densités de spins du pro essus est déterminée à partir de la somme : ρ (~r) =

X

ρσ (~r) = ρ↑ (~r) + ρ↓ (~r)

(1.22)

σ=↑,↓

l'o

upation éle tronique totale forme :

n

n=

est exprimée en fon tion de es densités de spins sous la

X

σ=↑,↓

nσ =

X Z

ρσ (~r) d~r

(1.23)

σ=↑,↓

(a) nσ représentant l'o

upation éle tronique par spin σ =↑ (↓) : Z n↑(↓) = ρ↑(↓) (~r) d~r (b)

2.4. Convergen e du pro essus self- onsistant

Les tests de onvergen e permettant de mettre n aux pro essus auto- ohérents sont souvent liés à la satisfa tion de la ondition imposée sur les densités de sortie de deux itérations su

essives : (n+1) (n) ρout − ρout 6 10−6

(1.24)

Les s hémas de mélange (mixage) usuels mis en ÷uvre pour l'a

élération de la onvergen e de es pro essus itératifs sont de deux types :

14

Chapitre 1. Formalisme de la DFT  le s héma d'Anderson se servant d'un mixage linéaire des grandeurs de al uls {densités de harge, o

upations, potentiels, énergies, moments magnétiques ...} (n) (n) d'entrée Xin et de sortie Xout d'une itération n (a hevée) pour la génération de la grandeur d'entrée de l'itération n+1 à venir : (n+1)

Xin

(a)

(n)

(n)

= α.Xout + (1 − α) .Xin

(1.25)

En termes de densités de harges (X ≡ ρ), e s hème fait l'usage du mixage suivant : (n+1)

ρin

(b)

(n)

(n)

= α.ρout + (1 − α) .ρin

une telle pro édure mettant en ÷uvre de faibles valeurs de paramètres de mixage (α < 0.5 ) et permettant une bonne onvergen e des systèmes massifs s'est révélé moins able pour le as des systèmes omplexes (du type des lms et des multi ou hes), sa onvergen e devenant problématique et parfois même in ertaine.  le s héma modié de Broyden basé sur la méthode développée par Broyden en 1965, adaptée par Bendt-Zunger à la onvergen e des solides, détaillée par Srivastava et modiée (améliorée) par Vanderbilt & Louie de manière à permettre une meilleure

onvergen e et une rapidité des al uls self- onsistants omplexes [42℄,[43℄,[44℄,[45℄. La méthode de Broyden est un algorithme de al ul des solutions (zéros) de fon tions à NDIM dimensions à partir d'un ve teur de départ X (n=0) C'est une méthode appliquée aussi bien au as des optimisations géométriques né essitant une détermination du zéro d'un ve teur de for e qu'au as des onvergen es di iles à travers un al ul de la diéren e des fon tions d'entrée et de sortie de deux itérations su

essives (n & n + 1), pour haque zéro. Autrement dit, pour e dernier as, le but re her hé est la minimisation des fon tions (ve teurs de for e) F . dé rites par la diéren e entre les grandeurs d'entrée Xin et de sortie Xout :

(a)

F ≡ F (Xin ) = Xout − Xin

(1.26)

Contrairement à la méthode linéaire qui ne faisant qu'au résultat de l'itération pré édente, l'algorithme de Broyden a et avantage de tenir ompte de l'ensemble des données issues des itérations antérieures (à elle en ours), pour la onstru tion itérative de ses Ja obiens (i.e. du gradient de leurs fon tions F ). Le problème lié à une re her he de la grandeur Xin permettant d'annuler es fon tions F est résolue à travers une appli ation de la méthode Newton-Raphson donnant lieu à une solution X (un zéro, annulant la fon tion F ≡ F (Xin ) = 0) développée sous la forme : (n+1)

Xin

i−1 h (n) F (n) = Xin − J (n)

(b) (n) (n) Xin et Xout représentant respe tivement les grandeurs d'entrée et de sortie d'une itération   (n+1) n (en ours), Xin la grandeur d'entrée de l'itération n + 1 (à venir), J (n) la matri e

15 . (n) (n) (n) ∂Xin, j ) et F (n) la fon tion à du Ja obien itératif (ave des éléments Jij = −∂Fi annuler dénie par l'Eq. (26.a) omme une diéren e des grandeurs d'entrée et de sortie, (n) (n) F (n) = Xout − Xin .   La méthode ainsi dé rite né essitant une première estimation de Ja obien de départ J (1) ,

e dernier est dé rit sous la forme du mixage linéaire d'Anderson suivant : i−1 i h i h h 1 J (1) = γ (1) . β (1) = − [1] α ave

h tel que

(

(2)

i β (1) = [1]

(1)

h i γ (1) = −α [1]

&

(1)

Xin = α.Xout + (1 − α) .Xin (1) (1) F (1) = Xout − Xin

h i F (1) F (1) → J (1) = = (2) (1) (1) ∆X Xin − Xin

( ) A partir de la deuxième itération, le s héma de Broyden modié est utilisé pour la
détermination de la solution (zéro) .X (n) annulant la fon tion F (n) ≡ F X (n) = 0 dé rite par la loi linéaire suivante : o i n h (n+1) − X (n) F (n+1) = F (n) − J (n) . Xin (d)

Pour le as des itérations d'ordre (n > 2), la forme générale de la solution (zéro) de l'Eq. (26.d) est donnée par : (n+1)

Xin

i−1 h (n) .F (n) = Xin + J (n)

(e) la matri e du Ja obien itératif orrespondant à la neme itération est exprimée  (n)   sous  la forme d'un produit d'une matri e gamma γ par l'inverse d'une matri e bêta β (n) : i−1 i h i h h J (n) = γ (n) × β (n)

(f)

La version modiée (améliorée) de la méthode de Broyden établie par Vanderbilt-Louie est basée sur un développement de la matri e du Ja obien itératif en faisant appel à deux fa teurs de poids wℓ et w′ (ae tés à ha une des itérations antérieures ℓ  (1)respe tivement  < n ainsi qu'au Ja obien initial J , n=1 [4℄ : i h J (n) =

( |

ℓ=1.n

)

{z

}

i h X (ℓ) wℓ2 .∆F (ℓ) .∆Xin w′2 . J (1) − [

γ (n)

]

(

× w′2 . [1] + |

X

2 (ℓ) wℓ2 . ∆Xin

ℓ=1.n

{z −1 [β (n) ]

)−1 }

16

Chapitre 1. Formalisme de la DFT

(g) (n)

∆Xin et ∆F (n) désignant respe tivement les diéren es des grandeurs d'entrée Xin et des for es orrespondantes, entre deux itérations su

essives n et n+1 : ( (n) (n+1) (n) − Xin ∆Xin = Xin ∆F (n) = F (n+1) − F (n)

(h)

une fois le Ja obien déterminé, le systéme linéaire d'équations issu de (26.e) est résolu (n+1) (n+2) (n+1) − Xin = Xin de manière à déterminer la diéren e ∆Xin et à travers elle la (n+2) : grandeur Xin (n+1)

[J (n+1) ].∆Xin

= F (n+1)

La onvergen e de e s héma de Broyden n'est atteinte qu'ave l'annulation de la diéren e (n+1) de densités ∆Xin to 0 : (n+2) (n+1) (n+1) − Xin 6 10−6 = Xin ∆Xin (j)

Ajouté au point fort de la méthode   de Broyden standard ( onventionnelle) mettant en oeuuvre des Ja obiens itératifs J (n) améliorés au fur et à mesure des itérations, les auteurs du s héma modié ont révélé qu'une utilisation de fa teurs de poids xes (indépendants de l'itération en question {w' = 0.01 & wℓ = 1 } orait à leur méthode une plus grande exibilité [4℄,[46℄. Les paramètres de mixage α utilisées pour faire onverger la plupart des

al uls sont en général de faibles valeurs 6 0.5. En fait, les valeurs de es paramètres sont

ara térisées par une dé roissan e rapide ave l'augmentation des dimensions (taille) de la

ellule al ulée et plus parti ulièrement dans le as des systèmes métalliques. En d'autres mots, les valeurs des paramètres de mixage sont d'autant faibles que la omplexité des systèmes traités est importante, Tableau 1. :

1.1  Valeurs des paramètres de mixage appliquées pour la onvergen e des diérents systèmes métalliques al ulés par la méthode de la TBLMTO-ASA

Tab.

17 2.5. Approximation de la densité lo ale

L'approximation de la densité lo ale est l'une des trois générations d'approximations mises à jour pour le traitement des systèmes polyéle troniques [47℄ :  l'approximation de la densité lo ale (LDA)17 .  l'approximation généralisée du gradient (GGA)18 .  la méthode du potentiel optimisé (OPM )19 [41℄,[48℄,[49℄. La approximation de la LDA est utilisée omme un premier niveau d'appro he pour le

al ul des potentiels de l'XC et des ontributions énergétiques (énergies) orrespondantes, pour des systèmes manifestant de très faibles variations de leurs densités. L'idée de base d'une telle approximation est liée au traitement du système polyéle tronique inhomogène

omme étant lo alement homogène à travers sa division en volumes élémentaires de densités lo ales uniformes ( onstantes) , la ontribution de ha un de es volumes à l'énergie de l'XC est égale à elle d'un même volume o

upé par un gaz d'éle trons homogène [25℄,[26℄,[50℄. En pratique, des volumes atomiques sont utilisés omme volumes élémentaires de manière ~ du système périodique (ave des à faire orrespondre à ha un des sites atomiques R

ellules unitaires polyatomiques) une densité atomique lo ale ρR (~r) propre et un potentiel ee tif atomique lo al vR (~r) bien déni [51℄. Ces grandeurs lo ales sont identiques pour tous haque type d'atomes (équivalents) à l'intérieur des ellules du réseau ristallin. Les rayons ve teurs éle troniques ~rR permettant la lo alisation des éle trons à l'intérieur ~ étant exprimés en fon tion de leurs ve teurs positions, ~r à l'intérieur de ha un de sites R ~ , Fig. 2. du solide sous la forme : ~rR ≡ ~r − R

Fig.

1.2  Lo alisation et repérage des éle trons à l'intérieur d'une ellule unitaire de al ul.

Compte tenu des onsidérations de la LDA et en absen e de toute ex itation extérieure, le potentiel externe engendré par le réseau des noyaux et ae tant haque éle tron est donné par : vext (~r) =

X R

17 18 19

LDA : Lo al Density Approximation GGA : Generalized Gradient Approximation OPM : Optimized Potential Method

vext (~rR ) =

X −ZR R

rR

(1.27)

18

Chapitre 1. Formalisme de la DFT

vext (~r) représentant le potentiel externe rée par le réseau périodique des noyaux, vext (~rR ) ~ ae tant l'éle tron lo alisé en ~rR (sur un même site). le potentiel du noyau d'un site R

Ainsi, le potentiel ee tif de est é rit sous la forme : vKS (~r) =

Kohn-Sham

XZ R

=

X R

ΩR

orrespondant à ha un des éle trons du solide

′ X ZR ρR (~rR ) ′ + vXC [ρR ] d~rR + − ′ |~rR − ~rR | rR

(1.28)

R

{vH (~rR ) + vext (~rR )} + vXC [ρR ] | {z } vR (~ rR )

~ , l'Eq. (28) est développée sur le reste des Pour ha un des éle trons ~r = ~rR de l'atome R ′ ~ ~ atomes R 6= R du réseau : ( ) Z ′ ′ X Z ρR (~rR ) ZR ZR′ ρR′ (~rR ′) ′ ′ + vXC [ρR ] vKS (~rR ) = + rR − rR′ − ′ | d~ ′ | d~ rR − ~rR rR rR′ − ~rR rR′ ′ ΩR |~ ΩR′ |~ R′ 6=R {z } | {z } | vR (~ rR )

vR′ (~ rR )

(1.29)

Le premier terme représentant le potentiel lo al (atomique) agissant sur l'éle tron ~r = ~rR ~ , le se ond terme dé rivant l'eet du reste des onstituants des dans son propre atome R ′ ~ ~ atomes R 6= R sur et éle tron et le dernier terme elui de l'XC. En général, tenant ompte des Eqs (10) et (11), l'énergie de l'XC est interprétée en termes de l'é art d'énergie entre l'énergie inétique TS [ρ] des éle trons du système hypothétique (S ) et la fon tionnelle universelle G [ρ] du système réel : (1.30)

EXC [ρ] = TS [ρ] − G [ρ]

Ave l'introdu tion des eets de orrélations dé rits omme la diéren e entre l'énergie totale exa te du système EExact et elle de la HFA [21℄ : (1.31)

EC [ρ] = E Exact [ρ] − E HF A [ρ]

le terme de orrélation résultant est ainsi donné par : EC [ρ] = E Exact [ρ] − EHF

      [ρ] = {Telec [ρ] + Uelec−elec [ρ]} − TS [ρ] + UH [ρ] +EX [ρ]  | {z }  {z } |  E[ρ] EH [ρ] | {z } EHF [ρ]

= Telec [ρ] − TS [ρ]

(1.32)

19 ave

Uelec−elec [ρ] = UH [ρ] + EX [ρ]

Telec [ρ] désignant l'énergie inétique exa te du système polyéle tronique réel (ave intera tion), Uelec−elec [ρ] le terme d'intera tion (éle tron-éle tron), .TS [ρ] l'énergie inétique d'un gaz d'éle trons (sans-intera tion), UH [ρ] le terme éle trostatique ( oulombien) de Hartree et EX [ρ] le terme d'é hange. Ce dernier étant déni omme la partie présente dans l'énergie de Hartree-Fo k et absente dans l'énergie de Hartree : EX [ρ] = EHF [ρ] − EH [ρ]. Du point de vue du magnétisme, es orrélations entre les éle trons de spins antiparallèles (négligées par la HFA) sont favorables à une stabilisation de l'ordre paramagnétique (PM) aux dépends du ferromagnétisme (FM ) [9℄. LDA [ρ] et des potentiels orresponLes expressions générales des ontributions de l'XC EXC dants mises en ÷uvre par les approximations de la LDA et de la LSDA dans leurs al uls approximatifs mais rigoureux des systèmes polyéle troniques (manifestant des densités variant lentement δρ (~r) → 0) sont exprimées sous la forme : Z LDA EXC [ρ] = εLDA r ) d~r (1.33) XC [ρ] .ρ (~

ave

∂EXC = vXC [ρ] = ∂ρ

  ∂εXC εXC [ρ] + ρ ∂ρ

↑(↓) ave la prise en ompte des densités (~r) omposant la densité de harge ρ (~r),  ↑ ↓de  spins ρ LSD l'expression de l'énergie EXC ρ , ρ résultante est donnée par : Z  ↑ ↓  ↑  ↑ ↓ LSD ρ , ρ . ρ (~r) + ρ↓ (~r) d~r EXC ρ , ρ = εLDA (1.34) XC | {z } ρ(~ r)

LDA [ρ] est basé sur une utilisation du résultat du gaz éle tronique Ce al ul des énergies EXC homogène (de densité uniforme, ρ (~r) = ρo ) et des fon tionnelle de la densité d'énergie εLDA XC [ρ] déterminées par une pro édure de paramétrisation (ave une fon tion analytique appropriée) et développées sous la forme [52℄ : hom LDA εLDA [ρ] + εLDA [ρ] XC [ρ] ≡ εXC [ρo ]ρo =ρ(~ C r ) = εX

(1.35)

LDA LDA LDA EXC [ρ] = EX [ρ] + EC [ρ]

(1.36)

e développement de l'Eq. (35) est utilisé pour une dé omposition de l'énergie de l'XC en un terme d'é hange et un autre de orrélation20 :

(a) Pour le as spin-polarisé, la fon tionnelle de l'énergie de l'Eq. (36) est réé rite sous la forme :

20 X

 ↑ ↓  ↑ ↓  ↑ ↓ LSD LSD LSD ρ ,ρ ρ , ρ + EC ρ , ρ = EX EXC

: eX hange, C : Correlation

20

Chapitre 1. Formalisme de la DFT

(b) L'énergie d'é hange (X) est expli itée dans l'appro he de la LDA par : Z LDA EX [ρ] = εLDA [ρ] .ρ (~r) d~r X

(1.37)

εLDA [ρ] représentant la fon tionnelle de la densité d'énergie d'é hange (en unité du X Ryd/ele ) dé rite sous sa forme analytique exa te de Dira -Slater suivante [52℄ : εLDA [ρ] X

 1/3 3 3 1/3 ρ(~r) =− 4 π

pour le as spin polarisé, l'énergie d'é hange orrespondante est obtenue par une sommation des énergies d'é hange par spin [41℄ :  ↓  ↑  ↑ ↓ LDA LDA LSD ρ ρ + EX ρ , ρ = EX EX

(1.38)

Une expression analogue à elle de l'Eq. (37) est adoptée pour l'énergie de orrélation (C ) : Z LDA [ρ] = εLDA [ρ] .ρ (~r) d~r EC (1.39) C Néanmoins, la forme analytique appropriée de la fon tionnelle de orrélation εC [ρ] s'est révélée plus ompliquée à mettre au point que elle de l'é hangeεX [ρ]. Dans une simulation stoe hastique de l'équation de S hrödinger par la méthode de Monte-Carlo quantique dans la région des (QMC )21 , Ceperley-Alder ont al ulé ette fon tionnelle εQMC C,CA [ρ]  3 basses densités (i.e. pour les grands rayons rs ≈ 20a.u., tel que 4πrs 3 = 1/ρ), pour un gaz homogène d'éle trons et pour diérentes valeurs de densités de harge ρ [32℄. Ces méthodes de la Monte-Carlo sont le plus souvent appliquées pour des al uls numériques des gaz d'éle trons homogènes en intera tion, ne pouvant être traités par le modèle des éle trons libres (adapté aux gaz d'éle trons sans intera tion) [53℄. En rapportant les résultats obtenus par Ceperley-Alder dans la limite des basses densités et eux de Gell-Mann-Brue kner obtenus pré édemment (1957 ) ave l'approximation de la phase aléatoire (RPA)22 dans la limite des hautes densités(la limite des basses densités est adoptée pour la négligen e des eets de orrélation) [54℄. Pour le as des faibles rayons rs < RP A 5a.u., EC,GB (rs ) ≈ 0.0622 ln (rs ) − 0.142Ryd, Perdew-Zunger ont développé une fon tion d'interpolation faisant le ra

ord entre es deux limites. Cette dernière étant ajustée aux 1/2 résultats de Ceperley-Alder, ompte tenu de l'appro he de Padé (en rs ) , elle a donné lieu à une paramétrisation de εSIC C,P Z [ρ] allant au-delà de la RPA [9℄,[34℄,[47℄. En parti ulier, la paramétrisation de vBH72 utilisée dans le présent travail pour la détermination des termes (densités d'énergie et potentiel de l'XC est expli itée sous la forme : LDA↑(↓)

εXC 21 22

1/3 i  h ↑(↓) ↑(↓) A + B [ρ] = A [ρ] 2x ρ, x ≡ εRP XC,V BH

QMC : Quantum Monte-Carlo RPA : Random Phase Approximation

(1.40)

21 (

ave


1/3  P A [ρ] = −2 3ρ + γ. εF C,HL − εC,HL π F P B [ρ] = µP C − γ. εC,HL − εC,HL

l'expression d'une telle énergie de l'XC est développée sous sa forme générale : A εRP XC,V BH [ρ, ξ] = εX [ρ, ξ] + εC [ρ, ξ]   εX [ρ, ξ] = ε0X [ρ] + ε1X [ρ] − ε0X [ρ] .f (ξ) εC [ρ, ξ] = ε0C [ρ] + ε1C [ρ] − ε0C [ρ] .f (ξ)

ave

(1.41)

le terme des densités d'énergie d'é hange (i ≡ X ) utilisé dans e as est elui donné par la formule de Dira -Slater (l'é hange de la LSDA) :
1/3 1/3
- pour le as non-polarisé : ξ = 0 → ρ↑ = ρ↓ : ε0X [ρ] = − 43 π3 ρ - pour le as

FM


polarisé : ξ = 1 → ρ↑ = ρ, ρ↓ = 0 : ε1X [ρ] = 21/3 ε0X [ρ]

Pour les densités d'énergie de orrélation (i ≡ adoptées :

ave

C ),

les expressions de

ε0C [ρ] ≡ εP C,HL = −cP .F



rs rP



ε1C



rs rF



[ρ] ≡

εF C,HL

F (z) = 1 + z

3




= −cF .F

Lundqvist-Hedin

sont

  z 1 1 + − z2 − .ln 1 + z 2 3

 représentant le rayon de la sphère ontenant en moyenne un éle tron (4πrs3 3 = ρ−1 )
et f x↑(↓) ≡ f (ζ) la fon tion d'interpolation (des densités d'énergies d'é hange et de orrélation) mise à jour par von Barth-Hedin, pour les as non-polarisé (ξ = 0 ) et spin-polarisé (ξ = 1 ), et développée sous la forme [47℄ : rs

4/3

4/3

+(1−ξ) −2 f (ξ) = (1+ξ) 24/3 −2 4/3 4/3 −2−1/3 ( 1+ξ ) +( 1−ξ 2 ) = 2 1−2−1/3
(x↑(↓) )4/3 +(1−x↑(↓) )4/3 −2−1/3 ≡ f x↑(↓) = 1−2−1/3

ξ est le paramètre du spin-polarisé : ξ=

ρ↑ − ρ↓ ρ↑ − ρ↓ = ρ↑ + ρ↓ ρ

&

x↑(↓) =

(1 ± ξ) ρ↑(↓) = ρ 2

22

Tab.

Chapitre 1. Formalisme de la DFT

1.2  Valeurs des paramètres valeurs des paramètres de

von Barth-Hedin

En termes des potentiels de l'XC, eux adoptés par la paramétrisation de sont dé omposés en une partie d'é hange et une autre de orrélation :

von Barth-Hedin

i h i h i h ↑(↓) RP A RP A = vX ρ, x↑(↓) + vC ρ, x↑(↓) vXC,V BH [ρ, ξ] ≡ vXC,V BH ρ, x

ave un premier terme suivant la formule de

Dira -Slater

(1.42)

:

1/3 i  6ρ h x↑(↓) vX ρ, x↑(↓) = π

(a)

(vBH72 )

(1.43)

et un se ond terme paramétrisé sous la forme :      1/3 i  h

P 4 F ↑(↓) P F P P .f x ε − ε +µ + µ − µ − 2x↑(↓) − 1 .γ. εF − ε vC ρ, x↑(↓) = C,HL C C C C,HL C,HL 3 C,HL

(b)

ave

µP C

  rP = −cP .ln 1 + rs

F & µC

  rF = −cF .ln 1 + rs

Les valeurs des diérents paramètres mis en ouvre dans la paramétrisation de rapportées sur la Tableau 2.

vBH72

sont

23 2.6. Approximation généralisée du gradient

L'approximation pré édente de la LDA mettant en ÷uvre des densités homogènes ave une variation spatiale très lente s'est révélé appropriée aux al uls des atomes libres et des molé ules servant la bran he de la himie quantique. L'appli ation de la LDA à es as atomiques et molé ulaires a révélé une sous-estimation de 0.5% des énergies de leurs états fondamentaux et de 5%. des distan es à l'équilibre des molé ules traitées. C'est surtout une méthode exa te pour le traitement des gaz d'éle trons uniformes, ave des densités d'énergie εhom XC [ρ] déterminées par laméthode de la QMC. Pour le as des solides, une bonne des ription des solides métalliques est oerte par l`approximation de la LDA, ave sous-estimation de leurs paramètres de réseau de 2%, de leurs énergies ohésives de 25% et de leurs bandes d'énergie de quelques pour ents. Néanmoins, les onditions d'homogénéité des densités de harge sont di iles à satisfaire par d'autres types de solides, expliquant ainsi les quelques limites (é he s) manifestées la LDA lors de son étude des systèmes en question. D'un oté, l'utilisation de ette approximation pour le al ul des solides semi ondu teurs a révélé une forte sous-estimation (40 à 50% ) des largeurs de bandes interdites (gaps ), 'est la raison du re ours à l'approximation généralisée du Gradient (GGA) permettant d'aller au-delà de la LDA en améliorant ses résultats, sans pour autant les orriger omplètement. En fait, e problème du gap est remédié ave l'utilisation de la méthode GW de Hedin développée sur la base sur la théorie des perturbations (au premier ordre) [55℄. De l'autre oté, la di ulté ren ontrée par la LDA lors de sa des ription des systèmes fortement orrélés { as des Oxydes : NiO, CoO, FeO...} manifestant de fortes orrélations (éle tron-éle tron) intrasites entre les bandes d des éléments des séries de transition (et f pour le as des terres rares) est surmontée par une orre tion en termes d'intera tion intrasite U de Hubbard apportée à

ette approximation et donnant lieu à la version LDA+U [56℄,[57℄,[58℄. D'un point de vue théorique, une prise en ompte de la non-homogénéité des densités de harge a né essité un re ours à des appro hes allant au-delà de la LDA. La première appro he développée dans le but d'améliorer la LDA fut l'approximation généralisée étendue (GEA)23 [59℄. Même si ette dernière a é houé dans la des ription de l'ensemble des solides, elle a onstitué néanmoins une bonne base de travail à de nombreux travaux ultérieurs du type du modèle du potentiel optimisé (OPM )24 [49℄,[61℄. En plus, les deux approximations de la densité moyenne (ADA)25 et de la densité pondérée (WDA)26 développées par la suite par Gunnarson et al. de manière à introduire les non-lo alités n'ont pas apporté une amélioration signi ative aux résultats de la LDA [20℄,[47℄. An de remédier aux limites de la LDA, la non-lo alité des densités de harge est introduite en termes de ontributions additionnelles à la densité d'énergie de l'XC (en gradients ~ ), des ontributions représentants l'eet des densités lo ales des atomes de la densité ∇ρ voisins. L'approximation mettant en ÷uvre e type de orre tions est dite Approximation généralisée du gradient (GGA) [60℄. Cette approximation de la GGA a pour base prin ipale l'introdu tion de la non-homogénéité ara térisant les densités de harge des systèmes polyéle troniques réels à travers une substitution de hla densité i d'énergie de l'XC de la ~ LDA GGA ρ, ∇ρ , dépendant non seulement LDA (LSDA) εXC [ρ] par une densité d'énergie εXC 23 24 25 26

GEA : Generalized Expanded Approximation. OPM : Optimized Potential Model. ADA : Average Density Approximation (1976 ). WDA : Weighted Density Approximation (1977 ).

24

Chapitre 1. Formalisme de la DFT

de la densité de harge mais aussi du gradient de ette densité. En général, les diérents versions de GGA introduisants semi-lo alement les non-homogénéités de la densité (à travers un développement de l'énergie de l'XC en série de puissan e de la densité et de ses gradients) ont amélioré les résultats de la LDA pour le as des des systèmes atomiques et molé ulaires, leur su

ès dans le as des solides est sensiblement dépendant de la nature des systèmes traités et des propriétés re her hées. En eet, pour les systèmes métalliques à base d'éléments de séries de transition, la GGA a oert une meilleure des ription de leurs propriétés à l'état fondamental que la LDA, ette derniére sous-estimant les paramètres de réseaux des séries {4d & 5d } ainsi que les énergies des éléments magnétiques 3d. Les énergies ohésives de es séries de transition al ulées par la GGA sont meilleures que

al ulées par la LDA, les deux approximations donnant aussi lieu à de bonnes stru tures de bandes (très similaires à elles de l'expérimental). Le traitement des semi- ondu teurs et des isolants massifs par la GGA a révélé une sous-estimation des valeurs de leurs gaps d'énergie. La diversité des versions de la GGA développées à nos jours27 et se diéren iants essentiellement par la forme de leurs termes en gradient est liée la diversité de hoix du gradient de leurs densités (et don des densités d'énergie de l'XC orrespondantes), un

hoix absent dans la LDA adoptant une fon tionnelle de l'XC dénie d'une seule et unique façon [34℄,[36℄,[40℄,[62℄,[63℄,[64℄. L'ensemble de es versions de la GGA a permis d'aller au-delà de la LDA (en remédiant à ertaines de ses limites) mais sans pour autant la

orriger omplètement. En fait, le rle joué séparément par ha un des termes additionnels mis en jeu dans les diérentes versions de la GGA dans le développement de l'énergie de l'XC est en ore mal déni (à la limite des basses densités, les termes en ore mal-dénis des non-lo alités de la GGA sont annulées de manière à donner lieu à la LDA), Perdew-Wang (1991 ) ont utilisé la paramétrisation de Rasolt-Geldart allant au-delà de la RPA) dans leurs

al uls de la limite des hautes densités. La version de l'approximation de la GGA adoptée dans le présent travail est elle LM-GGA développée par Langreth-Mehl-Hu [35℄,[37℄,[65℄. D'un point de vue théorique, les diérentes versions de la

ommune d'une é riture de l'énergie de l'XC sous la forme : GGA LDA EXC [ρ] = EXC [ρ] +

Z

GGA

sont basées sur l'idée

h i ~ . d~r gXC ρ, ∇ρ

(1.44)

le développement généralisé du gradient est interprété en termes d'un développement de

ette énergie en une somme innie des termes en gradient additionnels ( orre tifs) :

GGA LDA EXC [ρ] = EXC [ρ] +

Z

BXC [ρ].

2 ~ ∇ρ ρ4/3

d~r +

Z

CXC [ρ].

l'expression de l'énergie de l'XC mise en ÷uvre par la version

4 ~ ∇ρ ρ8/3

d~r + ...

LM-GGA

de

(1.45)

Langreth-Mehl

27 Tout au long de es dernières années, une grande variété de versions d'approximations du Gradient sont mises à jour dans le but d'améliorer autant que possible les résultats de al uls des diérents mileux matériels traités ( onduteurs, isolants, semi onduteurs, oxides ...) : PW-GGA, LM-GGA, PBE-GGA, Meta-GGA,

Hybrid-GGA, EXX ...

25 est expli itée sous la forme : 2  ∇ρ Z  ~ 7 −λ GGA LDA RP A 2e − . 4/3 d~r EXC,LM [ρ] = EX,exact [ρ] + EC [ρ] +a. 9 ρ {z } |

(1.46)

LDA [ρ] EXC

λ = b.

~ ∇ρ ρ7/6

, a=

π 4/3 16 (3π 2 )

= 2.143 × 10−3 , b = (9π)1/6 .f = 0.269, f = 0.15

(a) le terme de orrélation utilisé dans ette paramétrisation est elui de von ( al ulé par la RPA) : Z i h  ↑ ↓   RP A EC ρ , ρ = EC,V BH ρ↑ , ρ↓ = εC,V BH ρ, x↑(↓) .ρ (~r) d~r

Barth-Hedin

(b)

En termes de potentiels, eux utilisés par ette méthode LM-GGA sont déduits du potentiel RP A RP A ≡ vC,V de von Barth-Hedin vC BH :     ~ |2 |∇ρ  LDA RP A GGA −1/3 2 ∆ρ  − v [ρ] = v [ρ] + v [ρ] +2a.ρ −18f  X,exact C ρ ρ   XC,LM | {z } LDA vXC [ρ]     ~ ∇ ~ |∇ρ ~ |
~ 2
∇ρ.  11 7 2 |∇ρ| 2 −λ  −2e 1 − λ2 ∆ρ − − λ + λ + λ (λ − 3)  ~ | ρ 3 6 12 ρ2 2ρ.|∇ρ Cette version de Langreth-Mehl de la GGA (LM-GGA) a révélé son e a ité pour une investigation des propriéts physiques des systèmes atomiques, molé ulaires, solides et des surfa es à travers, d'un oté, une nette amélioration des résultats de la LDA et, de l'autre

oté, une possibilité de al ul ave une pré ision omparable sinon meilleure que elle de la

orre tion self-intera tion (SIC )28 dans le as des atomes, des molé ules, des solides et des surfa es métalliques [34℄,[66℄,[67℄. En fait, ette version LM-GGA adoptant la RPA omme méthode d'approximation de base est appropriée aussi bien aux systèmes uniformes qu'aux systèmes non-uniformes tout en manifestant une amélioration substantielle sur la RPA [35℄. La orre tion non-lo ale apportée par la LM-GGA ouplée ave les termes lo aux de la RPA de vBH72 ont permis une bonne approximation de al ul des termes des énergies de l'XC des systèmes métalliques [31℄.

28

SIC : Self Intera tion Corre tion

Bibliographie [1℄ O. K. Andersen, Phys. Rev. B12, 3060 (1975). 1, 2 [2℄ J. C. Phillips and L. Kleinman, Phys. Rev. 116, 287 (1959), L. Kleinman and J. C. Phillips, Phys. Rev. 116, 880 (1959). 1 [3℄ H. L. Skriver, in The LMTO Method (Springer-Verlag, Berlin, 1984), H. L. Skriver, in One Ele tron Theory of Metals : Cohesive and Stru tural Proeprties (Riso National Laboratory, Denmark, 1984). 1 [4℄ W. E. Pi kett, Comp. Phys. Rep. 9, 115 (1989). 1, 15, 16 [5℄ D. J. Singh, in Plane Waves, mi , Norwell, 1994). 1

Pseudopotentials and the LAPW Method (Klüwer

A ade-

[6℄ K. Horn and M. S heer, in Hand Book fo Surfa e S ien e Vol.II : Ele troni Stru ture, edited by S. Holloway and N. V. Ri hardson (Elsevier, Amsterdam, 2000). 1 [7℄ R. M. Martin, in Ele troni Stru ture bridge University Press, 2004). 1

: Basi Theory and Pra ti al Methods

[8℄ J. Friedel, in Transition Metals : Ele troni Stru ture of d-Band Metals, edited by J. Ziman (Cambridge University Press, 1969). 1

(Cam-

in The Physi s of

[9℄ F. Gautier, in Itinerant Magnetism (Université Louis Pasteur, Institut de Physique de Strasbourg, 1981), F. Gautier, in Métaux et Alliages de Transition (Université Louis Pasteur, Institut de Physique de Strasbourg, 1983). 1, 19, 20 [10℄ D. Stoeer, in

Stru ture Ele tronique, Croissan e et Propriétés Magnétiques des

Cou hes Métalliques Ultramin es et des Superréseaux,

Ph.D. thesis (Université Louis Pasteur, Strasbourg, 1992), D. Stoeer, G. Fabri ius, A. M. Llois and H. Dreyssé, Surf. S i. 269, 632 (1992). 1

[11℄ K. Louzazna, in Magnétisme des Films Ultramin es de Fer Epitaxiés dium(001), Thèse de Magistère (Université Ferhat Abbas, Sétif, 1998). 1

sur l'Iri-

[12℄ J. C. Slater, Phys. Rev. 51, 846 (1937). 2 [13℄ J. Korringa, Physi a 13, 392 (1947). 2 [14℄ W. Kohn and N. Rostoker, Phys. Rev. 94, 1111 (1954). 2 [15℄ S. Cottenier, in

Density Family Theory and the Family of (L)APW-Method : a step-

by-step Introdu tion

(Instituut voor Kern-en-Starlingsfysi a, K. U. Leuven, 2002) 2

[16℄ O. K. Andersen, O. Jepsen and O. Glötzel, in Highlights of Condensed Matter Theory, edited by F. Bassani, F. Fumi and M. P. Tosi (North-Holland, New York, 1985). 2 [17℄ W. R. L. Lambre ht and O. K. Andersen, Phys. Rev. B34, 2439 (1986). 2 [18℄ W. Kohn, Rev. Mod. Phys. 71, 1253 (1999), W. Kohn, in Ele troni Stru ture ter : Wave Fun tions and Density Fun tionals (Nobel Le ture, 1999). 3

of Mat-

28

Bibliographie

[19℄ R. Dreizler, in Relativisti Density Fun tional Theory, edited by S. Fiolhais, F. Nogueira and M. Marques (Springer-Verlag, Berlin, 2003). 3 [20℄ R. M. Dreizler and E. K. U. Gross, in Density Fun tional Theory : An the Quantum Many-Body Problem (Springer-Verlag, Berlin, 1990). 3, 23 [21℄ R. Parr and W. Yang, in Density Fun tional University Press, New York, 1989). 3, 6, 18 [22℄ J. Kohano, in in Ele troni Stru ture bridge Uiversity Press,......) 3

Approa h to

Theory of Atoms and Mole ules

(Oxford

Cal ulations for Solids and Mole ules

(Cam-

[23℄ P. Hohenberg and W. Kohn, Phys. Rev. 136, B864 (1964). 4, 5 [24℄ M. Born and J. R. Oppenheimer, Ann. der Phys. 84, 457 (1927). 4 [25℄ L. H. Thomas, Pro . Camb. Phi. So . 23, 542 (1927). 5, 17 [26℄ E. Fermi, Z. Phys. 48, 73 (1928). 5, 17 [27℄ W. Kohn and L. J. Sham, Phys. Rev. 140, A1133 (1965). 5, 9 [28℄ P. Giannozzi, in Computational

Approa hes to Novel Condensed Matter Systems : Ap-

pli ation to Classi al and Quantum Systems,

num Press, New York, 1995). 9

edited by D. Neilsonand M. P. Das (Ple-

[29℄ D. D. Koelling and B. N. Harmon, J. Phys. C : Sol. St. Phys. 10, 3107 (1977). 10 [30℄ L. Hedin and B. I. Lundqvist, J. Phys. C4, 2064 (1971). 11 [31℄ U. von-Barth and L. Hedin, J. Phys. C5, 1629 (1972). 11, 25 [32℄ D. M. Ceperley and B. J. Alder, Phys. Rev. Lett. 45, 566 (1980), D. M. Ceperley and M. H. Kalos, in Monte Carlo Methods in Statisti al Physi s, edited by K. Binder (Springer-Verlag, New York, 1979). 11, 20 [33℄ S. H. Vosko, L. Wilk and M. Nusair, Can. J. Phys. 58, 1200 (1980). 11 [34℄ J. P. Perdew and A. Zunger, Phys. Rev. B23, 5048 (1981). 11, 20, 24, 25 [35℄ D. C. Langreth and M. J. Mehl, Phys. Rev. B28, 1809 (1983). 11, 24, 25 [36℄ J. P. Perdew and Y. Wang, Phys. Rev. B33, 8800 (1986). 11, 24 [37℄ C. D. Hu and D. C. Langreth, Phys. Rev. B33, 943 (1986). 11, 24 [38℄ A. Be ke, Phys. Rev. A38, 3098 (1988). 11 [39℄ C. Lee, W. Yang and R. G. Parr, Phys. Rev. B37, 785 (1988). 11 [40℄ J. P. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996). 11, 24 [41℄ E. Engel and S. H. Vosko, Phys. Rev. B47, 2800 (1993), E. Engel and S. H. Vosko, Phys. Rev. B47, 13164 (1993), E. Engel and S. H. Vosko, Phys. Rev. B50, 10498 (1994). 11, 17, 20 [42℄ C. G. Broyden, Math. Comp. 19, 577 (1965). 14 [43℄ P. Bendt and A. Zunger, Phys. Rev. B26, 3114(1982). 14 [44℄ G. P. Srivastava, J. Phys. A : Math. Gen. 17, L317 (1984). 14 [45℄ D. Vanderbilt and S. G. Louie, Phys. Rev. B30, 6118 (1984). 14 [46℄ G. Krier, O. Jepsen, A. Burkhardt and O. K. Andersen, in vesrion-47 (Max-Plan k Institute, Stuttgart, 1995). 16 [47℄ E. K. U. Gross, in DFT I-II & TDDFT Physi s, Trieste, 2000). 17, 20, 21, 23

I-II

The TBLMTO-ASA ode

(International Centre for Theoreti al

[48℄ R. T. Sharp and G. K. Horton, Phys. Rev. 90, 317 (1953). 17

Bibliographie

29

[49℄ E. Engel and R. M. Dreizler, J. Comp. Chem. 20, 31 (1999). 17, 23 [50℄ J. Kohano, in

Ele troni Stru ture Cal ulations and First-Prin iples Mole ular Dy-

nami s Simulations

(International Centre for Theoreti al Physi s, Trieste, 1998). 17

[51℄ D. Knab, in Propriétés Ele troniques, Magnétiques et Optiques des Composés Ordonnés et des Matériaux Multi ou hes Métalliques Ph.D. thesis (Université Louis Pasteur, Strasbourg, 1991). 17 [52℄ P. Kratzer, in A Fast Guide to Density Fun tional for Theoreti al Physi s, Trieste, 1999). 19, 20

Cal ulations

(International Centre

[53℄ A. J. James, in Solving the Many Ele tron Problem Ph.D. thesis (University of London, 1995). 20

in Quantum Monte-Carlo Methods,

[54℄ M. Gell-Mann and K. A. Brue kner, Phys. Rev. 106, 364 (1957). 20 [55℄ L. Hedin, Phys. Rev. 139, A796 (1965). 23 [56℄ R. N. S hmid, E. Engel and R. M. Dreizler, in

Appli ation of Impli it Density Fun -

tional Theory to 3d Transition Metals Monoxides Dimension,

J. von Neumann (NIC series, 2002), Vol. 9, p. 213-223. 23

edited by H. Rollnik and

[57℄ J. Hubbard, Pro . Roy. So . A 276, 238 (1963). 23 [58℄ V. I. Anisimov, J. Zaanen and O. K. Andersen, Phys. Rev. B44, 943 (1991). 23 [59℄ F. Herman, J. P. van Dyke and I. B. Ortenburger, Phys. Rev. Lett. 22, 807 (1969). 23 [60℄ R. van Leeuwen, in

Kohn-Sham Potentials Density Fun tional Theory

(1994). 23

[61℄ E. Engel and S. H. Vosko, Phys. Rev. B42, 4940 (1990). 23 [62℄ S. Kurth, M. A. L. Marques and E. K. U. Gross, in Ele troni Fun tional Theory (Freie Universitat, Berlin, 2003) 24

Stru ture : Density

[63℄ J. P. Perdew, Y. Wang and E. Engel, Phys. Rev Lett. 66, 508 (1991), J. P. Perdew and Y. Wang, Phys. Rev. B45, 13244 (1992). 24 [64℄ A. Dal-Corso, Phys. Rev. B64, 235118 (2001). 24 [65℄ M. Rasdolt, G. Malmstrom and D. J. W. Geldart, Phys. Rev. B20, 3012 (1979). 24 [66℄ A. Savin, V. Wedig, H. Preuss and H. Stoll, Phys. Rev. Lett. 53, 2087 (1984). 25 [67℄ A. R. E. Mohamed and V. Sahni, Phys. Rev. B31, 4879 (1985). 25