Master II, module 524 : la corrélation électronique

1.4  Langage diagrammatique

Pour travailler plus profondement avec les équations de la théorie de la perturbation, et pour notamment prouver une propriété importante, nous allons introduire un système graphique - un langage diagrammatique, d'après les travaux de Goldstone (1957).

Figure
#1Un diagramme représentant une contribution au quatrième order de perturbation Møller-Plesset. A gauche il n'y a que les points d'interaction biélectronique, au milieu le diagramme complet, et à droite le diagramme avec les lignes horizontales pour déterminer le dénominateur.

Une contribution à l'ordre n est représentée par n points dédoublés empilés (étape 1 sur la figure). Ensuite nous mettons dans chaque point une ligne descendante (trou ou orbitale occupée) et une ligne montante (particule ou orbitale virtuelle) pour relier tous les points du diagramme par des chemins fermés. Nous devrons établir toutes les possibilités de relier les points une par une, et trouver la contribution à l'énergie de corrélation de chaque diagramme. Pour cela il y a un jeu de règles à appliquer :

1. Nous sommons sur toutes les lignes internes du diagramme
   

   å i  å j  å k  å l  å a  å b  å c  å d

2. chaque liaison horizontale entre 2 points donne lieu à une intégrale biélectronique en notation des chimistes regroupant les orbitales qui se retouvent dans un point
   

   (ib|kc)(kl|cd)(ja|ld)(ia|jb)

3. nous mettons une ligne horizontale entre chaque paire de points. Chaque ligne donne un facteur au dénominateur par les énergies des orbitales å_{i Î occ} e_i - å_{a Î virt}e_a dont les lignes sont coupées par la ligne horizontale :
   

   (ei+ek-eb-ec)(ei+el-eb-ed)(ei+ej-ea-eb)

4. le signe de la contribution est donnée par le nombre de boucles fermées et le nombre de trous comme
   

   (-1)boucles + trous

5. si le diagramme possède d'un plan de symétrique vertical, nous devrons le multiplier par un facteur 1/2.
6. pour tenir compte des différents combinaisons de spins pour un système à couches fermées nous multiplions la contribution par un facteur 2^boucles

Appliqué au diagramme du quatrième ordre de la Figure 1.4 cela donne

(-1)4+2  22  å ijkl  å abcd  (ib|kc) (kl|cd) (ja|ld) (ia|jb) (ei+ek-eb-ec) (ei+el-eb-ed) (ei+ej-ea-eb)        .

Au deuxième ordre nous avons seulement deux diagrammes

Figure
#1Les deux diagrammes fondamentaux de MP2.

qui nous donnent tout de suite l'expression dérivée auparavant (équation 6) :

E(2)   =   1 2 é ê ë 2  å ij  å ab  (ia|ib)2 ei+ej-ea-eb - å ij  å ab  (ia|ib)(ib|ja) ei+ej-ea-eb ù ú û

A quoi bon ce langage diagrammatique - nous avons déjà établi les expressions analytiques ? Parfois il est utile, surtout si on s'intéresse aux développements méthodologiques. Tout de suite nous établissons qu'au deuxième order il n'y a que des déterminants diéxités (comme au troisième ordre parce que nous ne pouvons qu'aller de la réference |F₀ñ vers des diexcités |F_Dñ, ensuite de |F_Dñ®|F_Dñ pour refermer le chemin par |F_Dñ®|F₀ñ), et au quatrième ordre il y a des déterminants mono-, di-, tri- et quadriexcités à respecter : Figure

Nous voyons par la suite autres deux autres applications des diagrammes.

1.5  La compensation de termes au troisième ordre

Figure

Nous n'avons que deux énergies orbitalaires

e1   =  h11 + ( 2 J11-K11) = h11+ (11|11) e2   =  h22 + 2 J12 - K12 = h22 + 2 (11|22) - (12|12)

et donc qu'un seul dénominateur

D = E0(0) - Ek(0) = 2 (e1 - e2) < 0        .

L'énergie en ordre zéro est de E₀ = 2 N e₁ , et l'énergie Hartree-Fock de

EHF = 2 N h11 + N (2 J11 - K11) = N (2 h11 + J11)       ,

les deux bien proportionnelles au nombre de molécules N. La correction de l'énergie au premier ordre est

E0(1) = áF0|V|F0ñ = áF0|H-H0|F0ñ = EHF-E0(0) = -N J11        .

L'énergie au deuxième ordre nous fournit l'élément de matrice áF₀|V|F_kñ :

E0(2)  =  - å ijab  (2 (ia|jb)-(ib|ja)) (ia|jb) D = N K122 D   ® áF0| V|Fkñ = K12        .

Avec ces préliminaires nous pouvons exprimer tout de suite le deuxième terme de l'équation 7 (E⁽³⁾ = A⁽³⁾+B⁽³⁾) :

B(3)  =  -áF0|V|F0ñ  å k    áF0| V|Fkñ2 D2   =  -(-N J11) N  K122 D2 ~ N2

(8)

et nous nous apercevons que ce terme a une croissance quadratique avec le nombre de particules. Nous devons alors trouver dans A⁽³⁾ une contribution qui contrebalance cette proportionnalité inattendue.

Calculons alors le terme A⁽³⁾ :

A(3)  = å k  å m  áF0|V |FkñáFk| E0(0)-Ek(0) V |FmñáFm| E0(0)-Em(0) V |F0ñ = å k  å m  K122 D2 áFk|V|Fmñ

Par construction l'élément áF_k|V|F_mñ ne peut être différent de zéro qu'à condition que l'excitation k et l'excitation m soient sur la même molécule, donc

å k  å m  áFk|V|Fmñ  ®   N áF1[`1]2[`2] |V|F1[`1]2[`2] ñ

(9)

Nous savons que

á11 _ 1   1  ¼2i _ 2   i  ¼|H | 11 _ 1   1  ¼2i _ 2   i  ¼ñ = (N-1) ( _ 1     _ 1   |11) +( _ 2     _ 2   |22) + 2(N-1) h11 + 2 h22 á11 _ 1   1  ¼2i _ 2   i  ¼|H0 | 11 _ 1   1  ¼2i _ 2   i  ¼ñ = 2(N-1) e1 +2 e2

L'insertion des e et la différence H-H₀ nous donnent l'élément de matrice cherché :

áFk|V|Fmñ = -(N-1) J11 + J22 -4 J12 +2 K12

et alors

A(3)  =  - N2J11K212 D2 + N  J11+J22-4 J12+2 K12 D2 K212        .

(10)

La somme d'équation 8 et équation 10 et (E⁽³⁾) a donc une croissance correcte avec le nombre de particules par la compensation des termes quadratiques en N.

Dans un développement diagrammatique nous n'avons pas cette compensation : la partie áF₀|V|F₀ñ du terme B⁽³⁾ n'ayant pas d'équivalent diagrammatique, ce terme doit être nécessairement compensé quelque part. Dès le début nous n'incluons que des termes à croissance proportionnelle au nombre de particules, mais nous avons 12 diagrammes différents à considérer au lieu de deux termes analytiques.

Figure
#1Les 12 diagrammes à respecter au troisième ordre de perturbation. Pour chaque diagramme le signe et le poids est indiqué.

Dans notre système modèle les interactions biélectroniques ne sont différents de zéro qu'à condition que les 4 orbitales impliquées soient sur la même molécule. Par conséquent les connexions et les indices dans un diagrammes ne peuvent jamais sortir de cette molécule ; nous avons alors un diagramme par molécule pour n'importe quel ordre de perturbation et alors toujours une croissance linéaire avec le nombre de particules. Ceci constitue le théorème des agrégats liés (Linked Cluster Theorem) de Goldstone.

Figure
#1Le diagramme à gauche donne lieu à une contribution toujours proportionnelle à N, tandis que le diagramme non-lié à droite amène à une double somme et, par conséquent, à une croissance quadratique de la contribution. Par construction les diagrammes non-liés sont exclus de la série diagrammatique.

1.6  Sommation infinies

Voyons alors une façon de faire des sommations sans effort supplémentaire. Nous prenons les deux diagrammes fondamentaux (Figure 1.4) et nous introduisons des interactions entre les propagateurs sans changer d'indices du diagramme. La figure montre un exemple que nous allons regarder plus en détail : Figure

Par la sommation immédiate d'une série géométrique la sommation infinie peut s'écrire comme

å ijab  (ia|jb)2 Dijab æ ç è 1- (ii|aa) Dijab + (ii|aa) Dijab (ii|aa) Dijab - (ii|aa) Dijab (ii|aa) Dijab (ii|aa) Dijab + ¼- ¼ ö ÷ ø =   å ijab  (ia|jb)2 Dijab é ê ë ¥ å n = 0  æ ç è - (ii|aa) Dijab ö ÷ ø n   ù ú û   =   å ijab  (ia|jb)2 Dijab é ê ê ê ê ë 1 1+ (ii|aa) Dijab ù ú ú ú ú û =   å ijab  (ia|jb)2 Dijab+Jia

avec l'intégrale biélectronique J_ia = (ii|aa) et le dénominateur D_ij^ab = e_i+e_j-e_a-e_b. La sommation infinie se réduit alors en une intégrale supplémentaire dans le dénominateur. En incluant toutes les possibilités de dessiner des diagrammes sans changer les indices Figure

nous réalisons un shift général

Dijab ® Dijab - ~ J   ij  - ~ J   ab  + ~ J   ia  + ~ J   ib  + ~ J   ja  + ~ J   jb

où [J\tilde]_mn est la combinaison J_mn-K_mnd_smsn. Elle est dépendant des spins des deux orbitales impliquées. Les éléments d'échange K ne sont présents que si les orbitales sont occupées par des électrons de même spin.

Le nouveau dénominateur peut s'écrire comme

áFijab|H-EHF|Fijabñ

qui nous amène à la perturbation dite Epstein-Nesbet (1955) avec une autre décomposition de l'hamiltonien complet :

H0EN   =   å I   |FIñáFI| H|FIñáFI| VEN   =   å I ¹ J   |FIñáFI| H|FJñáFJ| E(0)0   =  áF0| æ è å I    |FIñáFI|H|FIñáFI| ö ø |F0ñ = áF0|F0ñáF0|H|F0ñáF0|F0ñ = EHF E0(1) =  0

et les numérateurs áF₀|V^EN|F_Iñ = áF₀|H|F_Iñ = áF₀|V^MP|F_Iñ de la série Møller-Plesset. Cette série de perturbation implique déjà au deuxième ordre des sommations de diagrammes de tous les ordres de la série Møller-Plesset. Néanmoins, ses propriétés sont très difficiles à évaluer, ce qui fait que cette décomposition de l'hamiltonien n'est que rarement utilisée aujourd'hui.

2  Interaction de configurations

Nous gardons la condition de normalisation áY|Yñ = 1 par un multiplicateur l de Lagrange, ce qui nous donne un fonctionnel L( { c_I};l) à minimiser :

E({cI})  =  ác0F0+ å I   cIFI | H|c0F0+ å J   cJFJ ñ L({cI};l)   =  E({cI}) -  l (áY|Yñ-1)

Les dérivées nous amènent avec áF_I|H|F_Jñ = áF_J|H|F_Iñ à un système d'équations linéaires

0  =   ¶L({cI};l) ¶cI   =  2  æ è cIáFI|H|FIñ+ å J ¹ I  cJ áFI|H|FJñ ö ø - 2 l æ è å J  cJ áYI|YJñ dIJ  ö ø =  2  æ è c0 áF0|H|FIñ+ å J  cJ áFI|H|FJñ- l cI ö ø 0  =   ¶L({cI};l) ¶c0   =  2  æ è c0áFI|H|FIñ+ å J  cJ áF0|H|FJñ- l c0 ö ø

Nous recherchons alors un vecteur propre (c₀, {c_I}) et une valeur propre l (la plus petite pour l'état fondamental) de la matrice de H dans la base de tous les déterminants excités

æ ç ç ç ç è H00 ¼ H0I ¼ : ··· H0I ¼ HII ¼ : ··· ö ÷ ÷ ÷ ÷ ø æ ç ç ç ç è c0 : cI : ö ÷ ÷ ÷ ÷ ø   =  l  æ ç ç ç ç è c0 : cI : ö ÷ ÷ ÷ ÷ ø

(11)

Nous nommons l l'énergie totale E et nous normons la fonction en divisant par c₀ (normation intermédiaire). La première ligne du système à vecteur propre devient alors une équation pour déterminer l'énergie de corrélation

áF0 | H | F0 ñ+ å I  cI áF0 | H |FI ñ  =  l  = Etot   =  EHF + ECorr        .

D'où l'expression générale de l'énergie de corrélation variationnelle

ECorr   =   å I   cI áFI|H|F0ñ

(12)

où I ne sont que les déterminants diexcités. L'argument de la restriction sur les diexcités est le même que pour la perturbation MP2 : les monoexcitations n'interagissent pas avec la référence grace au théorème de Brillouin, et les excitations au-delà des diexcitations ne peuvent pas non plus interagir avec la référence par un opérateur biélectronique.

Equation 12 est donc toujours utilisable, même si la fonction d'onde contient des déterminants triexcités, quadriexcités etc. Par contre, en perturbation, les déterminants quadriexcités contribuaient explicitement dans l'expression de l'énergie à partir du quatrième ordre - nous n'avions pas de fonction propre de l'hamiltonien. En perturbation (ou autre méthode où HY ¹ E Y) nous sommes en principe amenés à calculer explicitement pour Y = F₀+å_I c_IF_I

ECorr = Etot - EHF = áY|H|Yñ áY|Yñ - áF0|H|F0ñ = EHF+ å I  [ 2cI áFI|H|F0ñ+ cI2 áFI|H|FIñ] + å I < J  2 cIcJ  áFJ|H|FIñ 1+ å I  cI2 - EHF (13)

Pour dériver des expressions bien plus simples de l'énergie E⁽²⁾ (éq. 6) et E⁽³⁾ (éq. 7), en perturbation, nous avions dévéloppé les coefficients de la fonction d'onde et l'énergie en même temps - d'où le besoin de diagrammes impliquant déterminants mono-, di-, tri- etc. excités explicitement. Mais nous n'avions pas une fonction d'onde, puis une énergie pae équation 13, que l'ensemble des amplitudes et énergies après convergence satisfera rigoureusement l'équation de Schrödinger et fournira une fonction d'onde corrélée. Équation 12 nous donne une énergie de corrélation variationnelle, équation 13 une énergie associée à une fonction d'onde quelconque, et équations 6 et 7 des énergies de corrélation perturbatives.

2.1  Réduction de la taille de la matrice d'IC

Une restriction immédiate et intuitive est de réduire le dégré d'excitation en n'utilisant que des déterminants mono- et di-excités - nous avons vu que les coefficients des diexcitations déterminent l'énergie de corrélation. Grace au théorème de Brillouin (interaction nulle entre référence et des monoexcités) notre matrice H prend la forme suivants :

H = æ ç ç ç è á0|H|0ñ 0 á0|H|Dñ 0 áS|H|Sñ áS|H|Dñ á0|H|Dñ áS|H|Dñ áD| H|Dñ ö ÷ ÷ ÷ ø

(15)

Des 1.8×10⁹ déterminants de l'IC complet il nous ne resteront que 12636 » 1.8 ×10⁴ déterminants des mono- et diexcitations, c'est-à-dire le problème à resoudre a été réduit d'un facteur 100 000.

Au-delà de cette réduction importante on peut se rendre compte qu'en fait la plupart des éléments de la matrice sont strictement zéro. Nous verrons pourquoi.

2.1.1  Par symétrie de spin

nous renseigne sur combien de fonctions propres de spin S nous pouvons construire à partir d'un nombre n d'électrons non-appariés donné. Le nombre N de déterminants excités à utiliser pour construire ces fonctions propres de spin est

N = n/2 å na = 0  æ ç è (n/2)! na!(n/2-na)! ö ÷ ø 2          .

Ceci ressemble légèrement au triangle de Pascal. Un exemple : quatre couches ouvertes, alors 1+4+1 = 6 simples déterminants disponibles pour former six fonctions propres de spin linéairement indépendantes, dont 2 singlets, trois triplets et un quintet. Avec huit électrons (déterminants au moins quadriexcités) nous aurons 1+16+36+16+1 = 70 différents déterminants pour 14 singlets indépendants ; la réduction de la matrice à diagonaliser devient un fait important.

2.1.2  Construction de fonctions propres de spin

Une fonction propre de spin obéit à

S2 |Yñ  =  (h/2p)2 S(S+1) |Yñ

Avec les identités

S2   =  Sx2 + Sy2 + Sz2 =  S- S+ + (h/2p) Sz + Sz2 S-  =   å i  s-(i) ;   S+  =   å i  s+(i) S-S+   =   å i  s-(i) s+(i) + å i ¹ j  s-(i) s+(j) Sz |Yñ  =   1 2 (na - nb) |Yñ æ è å i  s-(i) s+(i) ö ø  |Yñ =  nb  (h/2p)2 |Yñ æ è å i ¹ j  s-(i) s+(j) ö ø  |Yñ = (h/2p)2  å P  Pab Permutation d¢un spina avec un spin b  |Yñ

nous obtenons finalement une expression pour S² :

S2  |Yñ  =   ì í î å P  Pab + 1 4 ((na - nb)2+2 (na+nb) ) ü ý þ  (h/2p)2 |Yñ

(16)

La méthode directe pour obtenir des fonctions propres de spin sera d'établir la matrice áF_I|S²|F_Jñ dans la base de déterminants simples, de diagonaliser cette matrice et de ne retenir que les fonctions propres correspondantes au spin global (valeurs propres) souhaité.

Pour n'obtenir que des combinaisons linéaires de bon spin, nous pouvons nous servir d'un opérateur de projection. L'opérateur

1 (h/2p)2 (S2 -  (h/2p)2 Sk(Sk+1))

appliqué sur une fonction de spin ne nous laisse que la partie qui ne correspond pas au spin S_k. Il suffit alors de retirer d'une fonction d'onde Y toutes les contributions qui ne correspondent pas au spin S_n souhaité par le produit

OSn  =   Õ k ¹ n  1 (h/2p)2 (S2 - (h/2p)2 Sk(Sk+1))

(17)

Si on voulait, par exemple, construire la composante du quartet (S = 3/2) du déterminant |f₁f₂[`(f)]₃ñ = |aabñ on appliquerait l'opérateur éliminant la contribution du doublet (S = 1/2) :

O3/2 |aabñ  = ì í î å P  Pab + 1 4 ((na - nb)2+2 (na+nb) ) ü ý þ |aabñ- 3 4 |aabñ =  |aabñ+ |abañ+|baañ  = |f1f2 f   3  ñ+ |f1 f   2  f3ñ+| f   1  f2f3ñ

Par symétrie nous avons trouvé alors tout de suite la seule combinaison linéaire correspondant à S = 3/2 et S_z = 1/2 à retenir :

Y = 1 Ö3 æ è |f1f2 f   3  ñ+ |f1 f   2  f3ñ+ | f   1  f2f3ñ ö ø

au lieu des trois déterminants séparés. Plus symboliquement pour |Yñ = a |Y_3/2ñ+ b |Y_1/2ñ et O_3/2 = 1/(^h/_2p)² (S²-³/₄(^h/_2p)²) :

O3/2|Yñ  =  a  3 2 · 5 2  |Y3/2ñ- 3 4  |Y3/2ñ+ b  3 4  |Y1/2ñ- b  3 4  |Y1/2ñ =  a  12 4 |Y3/2ñ = 3 a |Y3/2ñ

La composante |Y_1/2ñ a disparu, comme désiré.

2.1.3  En considérant la symétrie spatiale

Pour établir la symétrie spatiale d'un déterminant il nous faut maintenant un tableau de multiplication de représentations irréductibles pour le groupe spatiale auquel notre molécule appartient. L'affaire est simple pour des représentations irréductibles unidimensionnelles - les seuls groupes avec des représentations irréductibles unidimensionnelles sont les sous-groupes du groupe D_2h, ce qui explique pourquoi un certain nombre de programmes de chimie théorique ne peuvent traiter que ces opérations de symétrie. Sinon, la logique interne des programmes, et le maintien de tableaux de multiplication de représentations irréductibles multidimensionnelles et la décomposition des produits en nouvelles représentations irréductibles, devient rapidement très complexe.

Après analyse de la symétrie spatiale, la matrice H devient une matrice décomposée en blocs pour chaque symétrie spatiale différente, et, par conséquent, sa diagonalisation peut se faire bloc par bloc.

Par exemple (voir les exercices) le groupe C_2v a les représentations irréductibles a₁, a₂, b₁ et b₂ avec le tableau de multiplication suivant

C2v a1 a2 b1 b2 a1 a1 a2 b1 b2 a2 a1 b2 b1 b1 a1 a2 b2 a1

La partie non-nulle de la matrice de l'hamiltonien se réduit fortement puisqu'il n'y a que les produits directs a₁ Äa₁ = a₂ Äa₂ = b₁ Äb₁ = b₂ Äb₂ qui donnent a₁ : Figure

2.2  Exemple H₂

Pour combiner les deux réductions possibles nous regardons une seule molécule d'hydrogène H₂, dans un base minimale. Les deux orbitales moléculaires à former s et s^* sont de symétrie ``gerade'' (pair) et ``ungerade'' (impair). A partir de la référence nous pouvons construire 6 déterminants différents

Figure

pour lesquelles symétrie spatiale et valeur de S² sont données. Nous voyons que le seul déterminant ayant les mêmes symétries et de spin et d'espace est la diexcitation g[`g]® u[`u], donc il suffira de construire une matrice 2 × 2 pour un IC complet en base minimale.

La fonction d'onde est Y = F_{1[`1]</sub> + c<sub>2[`2]} F_{2[`2]</sub> avec le seul paramètre c<sub>2[`2]} à déterminer. Rappellons-nous les quantités nécessaires pour constriure la matrice H :

EHF   =  áF1[`1]|H|F1[`1]ñ = 2 h11+(11|11) E2[`2]  =  áF2[`2]|H|F2[`2]ñ = 2 h22+(22|22) áF1[`1]|H|F2[`2]ñ = (12|12) = K12

Le système à valeurs propres est donc

æ ç è EHF K12 K12 E2[`2] ö ÷ ø æ ç è 1 c2[`2] ö ÷ ø   =  E æ ç è 1 c2[`2] ö ÷ ø

(18)

avec les valeurs propres

E±   =  EHF + E2[`2]-EHF 2 D  ±   æ  ú Ö æ ç è E2[`2]-EHF 2 ö ÷ ø 2   +K212   =  EHF + D± Ö D2+K122 ECorr  (19)