Mouvement des ions

Lorsque nous avons présenté la théorie de la fonctionnelle statique ou dépendant du temps, il avait toujours été supposé que le mouvement des noyaux atomiques pouvait être négligé. Dans les exemples qui nous intéressent, à savoir les petits agrégats de sodium en particulier, ce mouvement doit être pris en compte lorsqu'on s'intéresse à des phénomènes sur une échelle de temps supérieure à une centaine de femtosecondes, surtout, lorsque comme ici, on considère des processus déposant une grande quantité d'énergie dans le nuage électronique. Une partie de cette énergie pourra alors se transférer aux degrés de liberté ioniques (en assimilant électrons de cur et noyau en un ion) de façon à déclencher leur mouvement. Il faut donc écrire des équations du mouvement pour ces ions. Nous allons résumer pour ce faire les discussions de [31].

Considérons d'abord de façon générale le système sous la forme de $N$ électrons et de $N_A$ noyaux, de charge $Z_A$ et de masse $M_A$, $A=1\ldots K$. Soit ${\bf R}_{A\alpha}$ le vecteur position du $\alpha$-ième noyau de l'espèce $A$. L'équation de Schrödinger du système s'écrit alors

$$i\hbar \frac{\partial}{\partial t} \Psi = \left[ \widehat{H}_e({... ...n}({\bf r}_1,\ldots,{\bf r}_N,\left\{{\bf R}_{A\alpha}\right\}t) \right] \Psi$$

$\Psi$ étant la fonction d'onde exacte du système complet:

$$\Psi=\Psi({\bf r}_1,\ldots,{\bf r}_N, \left\{ {\bf R}_{A\alpha}\right\},t)$$

Le Hamiltonien électronique est toujours :

$$\widehat{H}_e=\sum_{j=1}^{N}\left(-\frac{\hbar^2}{2m_e}\nabla^2_{... ... + \frac{1}{2}\sum_{j,k=1,j\ne k}^{N}\frac{e^2}{\vert{\bf r}_j-{\bf r}_k\vert}$$

le Hamiltonien nucléaire

$$\widehat{H}_n=\sum_{A=1}^{K}\sum_{\alpha=1}^{N_\alpha} \left( -\f... ...ta=1}^{N_\beta} \frac{e^2Z_AZ_B}{\vert{\bf R}_{A\alpha}-{\bf R}_{B\beta}\vert}$$

les indices $A\alpha$ et $B\beta$ étant naturellement différents. L'interaction électron-noyau, elle, s'écrit simplement

$$\widehat{H}_{en}=-\sum_{j=1}^N\sum_{A=1}^{K}\sum_{\alpha=1}^{N_A} \frac{e^2Z_A}{\vert{\bf r}_{j}-{\bf R}_{A\alpha}\vert}$$

En analogie avec les théorèmes de la fonctionnelle de la densité dépendant du temps pour les électrons seuls, on peut établir à partir d'un état initial $\Psi(t_0)$ une relation biunivoque (du type du théorème de Hohenberg-Kohn) entre l'ensemble des potentiels extérieurs et l'ensemble des densités électroniques et nucléaires (pour chaque espèce), à un instant donné.

De là découle un ensemble d'équations de Kohn-Sham dépendant du temps pour les noyaux et les électrons.

$$i\hbar\frac{\partial}{\partial t}\phi_j({\bf r},t)= \left(-\frac{\hbar^2}{2m_e}\nabla^2_{\bf r}+v_s[n,\{n_B\}]({\bf r},t)\right) \phi_j({\bf r},t)$$

$$i\hbar\frac{\partial}{\partial t}\psi_{A\alpha}({\bf R},t)= \left... ...}\nabla^2_{\bf R}+V_s^A[n,\{n_B\}]({\bf R},t) \right) \psi_{A\alpha}({\bf R},t)$$ (43)

$n$ représentant la densité électronique, et $n_A$ les densités nucléaires, définies par  $n_A=\sum_{\alpha=1}^{N_A}\big\vert\psi_{A\alpha}({\bf R},t)\big\vert^2$. $v_s$ a la même expression que lorsque les noyaux étaient considérés fixes, et $V_s^A$, potentiel de Kohn-Sham pour les noyaux, peut être décomposé de façon similaire en

$$V_s^A[n,\{n_B\}]({\bf r},t) = V^A_{\mbox{\scriptsize ext}} ({\bf R},t) -e^2Z_A\int \!d{\bf r}'\,\frac{n({\bf r},t)} {\big\vert{\bf R}-{\bf r}\big\vert}$$

$$+ e^2Z_A\sum_{B=1}^{K}\int\!d{\bf R}'\, \frac{Z_Bn_B({\bf R}',t)}... ...t{\bf R}-{\bf R}'\big\vert} + V_{\mbox{\scriptsize xc}}^A[n,\{n_B\}]({\bf R},t)$$

Les potentiels d'échange-corrélation pour les noyaux sont eux aussi des dérivées fonctionnelles de l'action, comme dans le cas avec noyaux fixes. Dans les systèmes moléculaires que nous considérons, l'échange entre noyaux identiques est faible car aux séparations nucléaires usuelles le recouvrement des fonctions d'onde nucléaires est peu important. Cependant, il faudrait en principe appliquer une procédure du type SIC (Self-Interaction Correction) pour annuler les termes d'auto-interaction contenus dans les termes de type Hartree de l'équation précédente.

Même avec ces simplifications, la résolution de cet ensemble d'équations reste très lourde. On peut cependant [31] reprendre la discussion précédente en ne considérant que les électrons de valence et en remplaçant le potentiel coulombien créé par les noyaux par un pseudopotentiel pertinent. C'est ce que nous allons faire dans toute la suite, employant les pseudopotentiels que nous avons présentés ; de plus, nous allons encore plus simplifier le problème en traitant le mouvement des noyaux en mécanique classique. Nous continuerons cependant à utiliser le terme coulombien dans les équations au lieu de pseudopotentiels, dans un souci de simplicité.

En appliquant le théorème d'Ehrenfest à l'équation C.1 de Kohn-Sham pour les noyaux, la trajectoire classique

$${\bf R}_{A\alpha}^{\mbox{\scriptsize class}}= \langle \psi_{A\alp... ...R}}\big\vert\psi_{A\alpha}(t) \rangle = \int \!d{\bf R}\,n_{A\alpha}({\bf R},t)$$ (44)

du noyau $A\alpha$ vérifie l'équation du mouvement

$$M_A\frac{d^2}{dt^2} {\bf R}_{A\alpha}^{\mbox{\scriptsize class}} (t)={\bf F}_{A\alpha}(t)$$

avec la force donnée par

$${\bf F}_{A\alpha}(t)=- \langle \psi_{A\alpha}(t)\big\vert\nabla_{\bf R}V_s^A\big\vert\psi_{A\alpha}(t) \rangle$$

Remarquons au passage que puisqu'en principe les équations de Kohn-Sham dépendant du temps donnent les densités exactes, l'équation C.2 donne la trajectoire exacte classique s'il n'y a qu'un noyau de l'espèce $A$. Sinon, on a un ensemble de particules indiscernables et les trajectoires de noyaux individuels ne peuvent pas être définies classiquement ; seule la trajectoire du centre de masse a un sens. On utilisera toujours cependant l'équation C.2 afin d'imposer une définition des trajectoires classiques des noyaux.

Avec l'approximation SIC dont nous parlions plus haut, l'expression de la force ${\bf F}_{A\alpha}$ devient

$${\bf F}_{A\alpha}(t) = -\int \!d{\bf R}\,n_{A\alpha}({\bf R},t) \Bigg[ \nabla_{{\bf R}} V^A_{\mbox{\scriptsize ext}} ({\bf R},t)$$

$$-e^2Z_A\int \!d{\bf r}\, \nabla_{{\bf R}}\frac{n({\bf r},t)} {\bi... ...\nabla_{{\bf R}} \frac{Z_Bn_B({\bf R}',t)} {\big\vert{\bf R}-{\bf R}'\big\vert}$$

$$+ e^2Z_A\sum_{\beta=1,\beta \neq \alpha}^{N_A}\int\!d{\bf R}'\, \... ...}} \frac{Z_An_{A\beta}({\bf R}',t)} {\big\vert{\bf R}-{\bf R}'\big\vert} \Bigg]$$

où $n_{A\beta}({\bf R},t)$ vaut $\big\vert\psi_{A\beta}({\bf R},t)\big\vert^2$. Les densités nucléaires seront cependant souvent très piquées autour des trajectoires classiques. Dans ce cas, toutes les intégrales de la forme $\int\!d{\bf R}\,n_{A\alpha}({\bf R},t)G({\bf R})$ peuvent se simplifier par un développement de Taylor de $G({\bf R})$ autour de la densité :

$$\int\!d{\bf R}\,n_{A\alpha}({\bf R},t)G({\bf R})= G({\bf R}_{A\al... ...}(t)) + \hbox{\scr{O}}({\bf R}-{\bf R}_{A\alpha}^{\mbox{\scriptsize class}})^2$$

car le terme d'ordre 1 disparaît de par la définition C.2 de la trajectoire classique.

Négliger les termes d'ordre supérieur revient à écrire $n_{A\alpha}({\bf R},t)=\delta({\bf R}-{\bf R}_{A\alpha}^{\mbox{\scriptsize class}})$. Les équations du mouvement se réduisent alors à

$$M_A\frac{d^2}{dt^2} {\bf R}_{A\alpha}^{\mbox{\scriptsize class}} (t) = -\nabla_{{\bf R}_{A\alpha}^{\mbox{\scriptsize class}}} \Bigg[ V^A... ...r},t)} {\big\vert{\bf R}_{A\alpha}^{\mbox{\scriptsize class}}-{\bf r}\big\vert}$$

$$+ e^2\sum_{B=1}^{K}\sum_{\beta=1}^{N_\beta} \frac{Z_AZ_B} {\big\v... ...criptsize class}}- {\bf R}_{B\beta}^{\mbox{\scriptsize class}}\big\vert} \Bigg]$$ (45)

avec $(B\beta)\neq(A\alpha)$ dans la somme.

Ces équations sont donc à résoudre en même temps que les équations de Kohn-Sham dépendant du temps pour les électrons, auxquelles elles sont couplées. Malgré les nombreuses hypothèses et simplifications qui ont été faites, le traitement ici présenté ne fait pas appel à l'approximation de Born-Oppenheimer : on résout des équations dépendant du temps pour les électrons et non pas des équations de Kohn-Sham statiques. D'autre part, l'orthogonalité des fonctions d'onde électroniques au cours du temps n'a pas besoin d'être imposée comme dans le schéma de Car-Parinello.