Etude ab initio des glaces d’ammoniac fluoré et hydraté sous conditions thermodynamiques extrêmes

Théorie de la fonctionnelle de la densité

 Dans les années 1960, les travaux de Paul Hohenberg, Walter Kohn et Lu Sham ont permis de donner une approche théorique quant à la résolution de l’équation électronique aux valeurs propres de He. Les concepts fondamentaux de la théorie de la fonctionnelle de la densité ont été rédigés dans les articles originaux de P. Hohenberg, W. Kohn et L. Sham [7, 8] qui font notamment référence aux théorèmes de Hohenberg-Kohn. Cette théorie a eu un impact majeur dans le calcul des structures électroniques des solides et à ce jour elle continue à être améliorée. Le physicien W. Kohn, fondateur principal de la théorie de la fonctionnelle de la densité, a obtenu le prix Nobel de Chimie en 1998 avec le chimiste et numéricien J. A. Pople.

Les théorèmes de Hohenberg-Kohn

 L’idée de départ de la théorie de la fonctionnelle de la densité consiste à reformuler la question de la résolution de l’équation de Schrödinger dans l’approximation adiabatique. Au lieu de chercher directement à résoudre cette équation qui dépend de la fonction d’onde poly-électronique ψ e , la théorie 6 Chapitre 1: Méthodes de modélisation de la fonctionnelle de la densité propose une alternative en considérant comme quantité de base pour la description du système la densité électronique n(r). Ceci permet de réduire le nombre de degrés de liberté de 3Ne à 3 puisque la densité électronique est définie en un point donné de l’espace. La théorie de la fonctionnelle de la densité est basée sur deux théorèmes qui ont d’abord été présentés et démontrés par P. Hohenberg et W. Kohn en 1964 [7] puis reformulés par R. M. Martin [9] en 2004 : Théorème. Pour un système d’électrons en interaction, le potentiel externe Vext(r) est déterminé de façon unique, à une constante près, par la densité électronique de l’état fondamental notée n0(r). Corollaire. Si le Hamiltonien du système est connu, il s’ensuit que la fonction d’onde poly-électronique associée à l’état fondamental et aux états excités peut être déterminée. Toutes les autres propriétés du système sont déterminées par la densité électronique de l’état fondamental. Ce théorème permet d’exprimer un certain potentiel Vext(r) comme une fonctionnelle de la densité électronique n(r). Ceci étant dit, revenons à l’expression mathématique du Hamiltonien électronique He (1.4). Dans l’expression de ce dernier, le terme Ve,n décrit un potentiel externe créé par les noyaux qui agit sur chaque électron indépendamment des autres électrons : Ve,n = X Ne i=1 Vext(ri) (1.14) avec Vext(ri) = − X Nn I=1 ZIe 2 |ri − RI | (1.15) Puisque le potentiel Vext(r) (1.15) dépend d’une configuration nucléaire donnée, la fonction d’onde poly-électronique dépend elle aussi paramétriquement de cette même configuration des noyaux. D’après le théorème précédent, le potentiel Vext(r) s’écrit donc comme une fonctionnelle de la densité électronique qui est fonction de la position des électrons : Vext[n(r)] = Z Vext(r)n(r)dr (1.16) où l’intégration se fait sur un volume infinitésimal dr. D’après le corollaire, l’énergie électronique totale Ee 2 est également une fonctionnelle de la densité électronique : E e [n(r)] = Te[n(r)] + Ve,e[n(r)] + Vext[n(r)] (1.17) En regroupant les fonctionnelles électroniques au sein de la fonctionnelle F[n(r)] nous obtenons : E e [n(r)] = F[n(r)] + Vext[n(r)] (1.18) L’Éq. (1.18) ne permet pas de connaître la densité électronique du système étudié. Un deuxième théorème basé sur le principe variationnel et les propriétés des fonctionnelles donne l’énergie du niveau fondamental si la densité électronique correspondante est la densité électronique de l’état fondamental. 

 Le potentiel Vn,n qui traduit l’interaction coulombienne classique entre les noyaux, doit être pris en compte dans le calcul de l’énergie totale du système. Mais ce terme n’intervient pas explicitement dans l’expression de la fonctionnelle E e . 

Théorie de la fonctionnelle de la densité 

Théorème. Pour un potentiel externe Vext(r) donné, l’énergie de l’état fondamental du système est donnée par le minimum global de la fonctionnelle Ee [n(r)]. La densité électronique qui minimise cette fonctionnelle est la densité électronique de l’état fondamental n0(r) : δE[n(r)] δn(r)     n0(r) = 0 Corollaire. La connaissance de la fonctionnelle Ee [n(r)] suffit pour déterminer l’énergie et la densité électronique de l’état fondamental. Dérivons l’Éq. (1.18) par rapport à la densité électronique : δEe [n(r)] δn(r) = δF[n(r)] δn(r) + Vext(r) (1.19) Il faut ensuite utiliser la méthode des multiplicateurs de Lagrange, qui dit qu’avec une fonctionnelle A[n(r)] nous construisons une nouvelle fonctionnelle B[n(r)] telle que : B[n(r)] = A[n(r)] − µC[n(r)] (1.20) avec : • C[n(r)] la contrainte extérieure, ici C[n(r)] = R n(r)dr = Ne. • µ est le multiplicateur de Lagrange et il s’agit du potentiel chimique associé à la conservation du nombre d’électrons. • A[n(r)] est égale à la fonctionnelle Ee [n(r)]. • B[n(r)] est une fonctionnelle quelconque. Ainsi l’Éq. (1.20) s’écrit : B[n(r)] = E e [n(r)] − µ Z n(r)dr (1.21) La méthode des multiplicateurs de Lagrange nous dit que δB[n(r)] = 0. En utilisant la définition de la dérivée d’une fonctionnelle, nous obtenons que : δB[n(r)] = Z δB[n(r)] δn(r) δn(r)dr = 0 (1.22) En dérivant l’Éq. (1.20) par rapport à la densité électronique il vient que : δB[n(r)] δn(r) = δEe [n(r)] δn(r) − µ (1.23) En remplaçant dans l’Éq. (1.22) on obtient : δB[n(r)] = Z  δEe [n(r)] δn(r) − µ  δn(r)dr = 0 (1.24) Ceci implique que : δEe [n(r)] δn(r) = µ (1.25) En remplaçant dans l’Éq. (1.19) il en résulte que : Vext(r) + δF[n(r)] δn(r) = µ (1.26) 8 Chapitre 1: Méthodes de modélisation L’Éq. (1.26) correspond à l’équation fondamentale du formalisme propre à la théorie de la fonctionnelle de la densité. Ainsi, pour accéder à l’énergie totale du système il faut pouvoir être en mesure de déterminer F[n(r)]. La méthode des multiplicateurs de Lagrange est un outil mathématique permettant de prendre en compte la contrainte du problème qui est la conservation du nombre d’électrons. Cependant, l’application de cette méthode n’est pas sans signification physique puisque µ correspond au potentiel chimique électronique, or à température nulle µ n’est autre que l’énergie de Fermi. En toute rigueur, le spin doit être introduit que si le potentiel externe dépend su spin. Sinon, l’énergie reste une fonctionnelle de la densité totale y compris pour les systèmes à multiplicité de spin différente de zéro. Cependant, l’introduction du spin [10] au sein de la théorie de la fonctionnelle de la densité amène une flexibilité qui améliore la qualité des approximations qui seront présentées dans la section 1.2.3. 

 Équations de Kohn-Sham

 La fonctionnelle F[n(r)] dépend de l’énergie cinétique des électrons Te[n(r)] et de l’interaction coulombienne mutuelle des électrons du système, Ve,e[n(r)]. L’approche choisie par W. Kohn et L. Sham en 1965 pour déterminer F[n(r)] et obtenir l’énergie du niveau fondamental consiste à introduire l’énergie cinétique totale Ts[n(r)] d’un système d’électrons fictifs Ss sans interaction mais de même densité électronique que le système réel noté Se. La fonctionnelle F[n(r)] est construite à partir d’une somme de trois termes dont l’énergie cinétique Ts[n(r)], un terme d’interaction coulombienne nommé également terme de Hartree EH[n(r)] et une fonctionnelle notée Exc[n(r)] : F[n(r)] = Ts[n(r)] + EH[n(r)] + Exc[n(r)] (1.27) où la fonctionnelle Exc[n(r)] est la fonctionnelle d’échange et corrélation qui contient des termes de correction sur l’énergie cinétique et l’interaction coulombienne : Exc[n(r)] = Ex[n(r)] + Ec[n(r)] = (Ve,e[n(r)] − EH[n(r)]) + (Te[n(r)] − Ts[n(r)]) (1.28) La fonctionnelle Exc[n(r)] est un terme qui, premièrement, prend en compte les effets d’échange. La fonction d’onde poly-électronique doit être antisymétrique par l’échange de deux électrons. Une conséquence importante est que deux électrons de spins parallèles ne peuvent pas occuper la même position. Ainsi, l’antisymétrisation de la fonction d’onde a pour effet de produire une séparation spatiale des électrons de spins parallèles ce qui réduit l’énergie coulombienne du système. On donne parfois de ce phénomène une représentation imagée selon laquelle un électron est entouré d’un “trou de Coulomb” et d’un “trou de Fermi” qui annihile, plus ou moins, l’apparition dans son voisinage d’un autre électron même spin. Cette réduction de l’énergie coulombienne du système correspond à l’énergie d’échange qui s’exprime comme la différence de l’énergie d’interaction Ve,e propre au système réel moins l’énergie de Hartree EH. Deuxièmement, l’énergie coulombienne du système peut encore être modifiée en augmentant la distance de séparation des électrons présentant des spins anti-parallèles. Cependant, la diminution des interactions coulombiennes s’accompagne d’une augmentation de l’énergie cinétique du système réel. La différence d’énergie entre le système réel et l’énergie du système calculée dans l’approximation Hartree-Fock est appelée énergie de corrélation (Ec[n(r)] = Te[n(r)] − Ts[n(r)]). La fonctionnelle Exc[n(r)] n’est pas connue exactement, et le choix d’une fonction d’échange et corrélation approximée constitue l’un des principaux choix d’approximation en théorie de la fonctionnelle de la densité. Nous reviendrons sur les différentes approximations possibles pour Exc[n(r)] dans la section 1.2.3. Si on utilise de nouveau l’équation fondamentale du formalisme propre à la théorie de la fonctionnelle de la densité (1.26) et on y injecte la nouvelle expression de la fonctionnelle F[n(r)] (1.27), nous obtenons l’expression suivante : 1.2 Théorie de la fonctionnelle de la densité 9 Vext(r) + Ts[n(r)] δn(r) + EH[n(r)] δn(r) + Exc[n(r)] δn(r) = µ (1.29) On regroupe le potentiel externe, la dérivée fonctionnelle du potentiel de Hartree et la dérivée fonctionnelle du potentiel d’échange et corrélation dans un même terme noté Vs(r) : Vs(r) + δTs[n(r)] δn(r) = µ (1.30) où le potentiel Vs(r) représente le potentiel externe qui agit sur les électrons du système Ss. La fonctionnelle de l’énergie de Hartree a pour expression : EH[n(r)] = Z Z n(r)n(r 0 ) |r − r 0 | drdr 0 (1.31) La dérivée fonctionnelle de l’Éq. (1.31) donne : EH[n(r)] δn(r) = Z n(r 0 ) |r − r 0 | dr 0 = VH(r) (1.32) Enfin, le terme associé à la dérivée fonctionnelle de l’énergie d’échange et corrélation est noté Vxc(r) : Vxc(r) = δExc[n(r)] δn(r) (1.33) Nous obtenons donc pour le potentiel Vs : Vs(r) = Vext(r) + VH(r) + Vxc(r) (1.34) La densité électronique du système Se d’électrons en interaction est égale à la densité électronique du système Ss d’électrons sans interaction si le potentiel Vs satisfait l’égalité (1.34). Rappelons que le potentiel Vs est le potentiel externe au sein duquel les électrons du système Ss évoluent. Ainsi, la signification physique de l’Éq. (1.34) est la suivante : chaque électron du système réel Se évolue dans un potentiel effectif de champ moyen. L’énergie cinétique Ts s’écrit comme : Ts[n(r)] = X Ne i=1 hϕi |Ts[n(r)]|ϕii = − ~ 2 2me X Ne i=1 hϕi |∇2 ri |ϕii (1.35) où les ϕi sont les orbitales mono-électroniques du système virtuel. La nouvelle équation de Schrödinger à résoudre prend alors la forme suivante :  » − ~ 2∇2 ri 2me + Vs(r) # ϕi(r) = iϕi(r) (1.36) Et nous définirons un Hamiltonien de Kohn-Sham comme étant : HKS = − ~ 2∇2 ri 2me + Vs(r) (1.37) avec le potentiel effectif Vs donné par l’Éq. (1.34). La densité électronique a pour expression : n(r) = X Ne i=1 |ϕi(r)| 2 (1.38) Ainsi, la résolution de l’équation de Schrödinger pour Ne électrons en interaction se ramène à la résolution d’une équation de Schrödinger pour un électron évoluant dans le potentiel de champ moyen 10 Chapitre 1: Méthodes de modélisation Vs. Ces équations sont auto-cohérentes au sens où VH et Vxc dépendent de la densité électronique par l’intermédiaire des orbitales mono-électroniques ϕi (elles sont résolues de manière itérative). Une densité électronique initiale est générée pour obtenir le potentiel effectif Vs, une nouvelle densité électronique est alors déterminée. À partir de cette nouvelle densité électronique, un nouveau potentiel effectif est calculé. Cette opération est répétée jusqu’à obtenir une densité électronique similaire, à un seuil de tolérance près, à la densité électronique du pas précédent. Pendant chaque itération, les potentiels Vxc et VH sont calculés. Ce processus de résolution des Éq. (1.36) correspond à ce qu’on appele la boucle auto-cohérente. À la fin, on obtient l’énergie électronique ainsi que la fonction d’onde poly-électronique par l’intermédiaire des orbitales mono-électroniques. Néanmois, la résolution des équations de Kohn-Sham nécessite de connaître le potentiel d’échange et corrélation, et par conséquent l’énergie associée. À ce jour, il n’existe pas une expression analytique exacte de l’énergie d’échange et corrélation. Cette dernière doit donc être approximée. 

 Approximations sur la fonctionnelle d’échange et corrélation 

Dans l’article original de W. Kohn et L. Sham [8] il est écrit que, dans la limite du gaz d’électrons homogène, les effets d’échange et corrélation ont un caractère local. L’approximation de la densité locale (ou LDA pour Local Density Approximation en anglais) a historiquement été la première approximation utilisée pour l’énergie d’échange et corrélation. Approximation de la densité locale La fonctionnelle de l’énergie d’échange et corrélation peut être approximée simplement comme l’intégrale, sur tout l’espace, de l’énergie d’échange et corrélation par particule  hom xc obtenue pour un gaz d’électrons homogène [11] : E LDA xc [n(r)] = Z n(r) hom xc (n(r))dr (1.39) Dans la littérature, il existe différentes méthodes de construction pour  hom xc [12–15]. Ces méthodes se basent sur des expressions mathématiques pour l’énergie d’échange et sur des approximations pour l’énergie de corrélation. Les travaux de D. M. Ceperley et B. J. Alder [16] ont permis de donner des valeurs précises de l’énergie de corrélation via des simulations de Monte Carlo quantique. L’approximation de la densité locale reste une bonne approximation pour l’énergie d’échange et corrélation dans des solides où la densité électronique est similaire à celle d’un gaz d’électrons uniforme. En revanche, la LDA présente des limites lorsqu’on s’intéresse à des systèmes isolés tels que des molécules. Le plus grand défaut de la LDA est qu’elle tend à rendre les électrons localisés. Ainsi, les systèmes d’électrons fortement corrélés sont mal décrits. De plus, la fonctionnelle LDA est connue pour surestimer l’énergie de la liaison hydrogène [17] et donne des gaps trop faibles pour les semiconducteurs et les isolants. Approximation du gradient généralisé Contrairement à la LDA, l’approximation du gradient généralisé (ou GGA pour Generalized Gradient Approximation) permet de prendre en compte la variation spatiale de la densité électronique. Les systèmes réels sont le plus souvent très différents d’un gaz d’électrons uniforme et c’est pourquoi la GGA prend en compte la variation spatiale de la densité électronique en exprimant l’énergie d’échange et corrélation en fonction de son gradient : E GGA xc [n(r)] = Z n(r) GGA xc (n(r), ∇n(r))dr (1.40)

 Théorie de la fonctionnelle de la densité

 Plusieurs méthodes de paramétrisation existent pour  GGA xc . On distingue deux catégories de fonctionnelles basées sur l’approximation du gradient généralisé. Dans la première catégorie, la fonction  GGA xc est ajustée sur des données expérimentales pour un grand nombre de molécules. Il s’agit de la fonctionnelle semi-empirique BLYP qui est formée de la somme d’une fonctionnelle d’échange, celle de Becke [18], plus la fonctionnelle de corrélation de Lee, Yang et Parr [19]. Concernant la deuxième catégorie, la construction de la fonction  GGA xc repose sur des principes de base de la mécanique quantique et cette construction n’inclut pas de données expérimentales. La fonctionnelle PBE [20] est construite de cette manière. Pour les simulations ab initio effectuées dans cette thèse, la fonctionnelle utilisée est une fonctionnelle de type GGA, en l’occurence la fonctionnelle PBE. La fonctionnelle PBEsol [21], considérée comme une fonctionnelle dite hybride, sera utilisée pour certains calculs.

 Application de la théorie de la fonctionnelle de la densité au cas des solides cristallins

 Dans toute notre étude, la théorie de la fonctionnelle de la densité a été utilisée pour réaliser des simulations sur des cristaux moléculaires et ioniques. Dans l’hypothèse d’un cristal parfait, une maille conventionnelle est reproduite à l’infini dans les trois directions de l’espace afin de générer la structure cristalline. La simulation d’un nombre donné d’atomes composant le cristal se fait en considérant que ces derniers sont contenus dans un volume fini communément appelé boîte de simulation. Cependant, ce nombre d’atomes contenus dans la boîte de simulation reste très petit devant le nombre d’atomes présents dans un échantillon expérimental et ce pour des raisons d’ordre calculatoire et de symétrie. La conséquence directe de l’utilisation d’une boîte de simulation est qu’une proportion non négligeable d’atomes se trouve à proximité de la surface de la boîte donnant lieu à des effets de bords, c’est à dire la présence de surfaces libres. Des conditions aux limites périodiques, appelées également conditions de Born-Von Kármán [22], sont nécessaires pour remédier à ce problème. Elles sont préférables aux conditions aux limites rigides (réflexion parfaite des parois) puisque ces dernières introduises d’autres effets de bords. Des images de la boîte de simulation sont construites autour de la-dite boîte afin de paver tout l’espace autour de la maille conventionnelle. Le système infini et périodique ne représente en aucun cas le système expérimental mais permet tout de même de réaliser des simulations avec un nombre raisonnable d’atomes.