Méthodes de simulation

La simulation numérique a pris son essor au début des années 50, lorsque les ordinateurs développés pendant la seconde guerre mondiale pour le décryptage de messages codés et la mise au point d’armes nucléaires, furent disponibles pour des applications non militaires. Depuis, l’utilisation de la simulation numérique pour l’étude et la compréhension des phénomènes physiques et chimiques est devenue une discipline de recherche à part entière. La simulation moléculaire constitue une approche complémentaire aux études expérimentales et aux modèles analytiques et théoriques. En effet, les observables macroscopiques mesurées expérimentalement sont, dans de nombreux cas, difficiles à interpréter au niveau microscopique à cause de la complexité des phénomènes mis en jeu. D’autre part, les approches théoriques se basent généralement sur des approximations qui peuvent être fortes rendant parfois le système théorique peu représentatif du système réel. Les prédictions théoriques sont alors difficiles à valider expérimentalement car le système modèle utilisé est trop éloigné de la réalité expérimentale. D’une manière générale, les simulations numériques peuvent être utilisées en complément des expériences pour l’interprétation des données expérimentales. En effet, elles permettent de découpler l’effet des différents paramètres expérimentaux afin d’étudier et de comprendre l’influence de chacun d’entre eux sur les phénomènes étudiés. D’autre part, elles permettent de valider ou d’affiner un modèle analytique permettant ainsi la rationalisation à l’échelle microscopique des tendances observées expérimentalement. Mais la force de la simulation moléculaire réside sans doute dans sa valeur prédictive des phénomènes physico-chimiques. En effet, les simulations numériques peuvent être utilisées pour faire des prédictions sur des phénomènes qui interviennent dans des conditions inaccessibles expérimentalement. Il est alors possible de réaliser des expériences in silico dans des conditions qui peuvent être extrêmes (température et pression) ou dangereuses, par exemple pour l’étude des réactions intervenant dans les réacteurs nucléaires, ou tout simplement pour étudier des phénomènes ultra-rapides ou très lents difficilement mesurables expérimentalement. De plus, elles donnent une vision microscopique de la matière, permettant l’étude à l’échelle atomique des phénomènes physico-chimiques (par exemple, identification des mécanismes réactionnels). En outre, les différentes méthodes de simulation numérique offrent plusieurs niveaux de description du système, de l’échelle microscopique à l’échelle mésoscopique voire même macroscopique, permettant ainsi d’aller vers une description multi-échelle du système étudié. Elles permettent également de mesurer un grand nombre d’observables microscopiques, qui peuvent ensuite être reliées aux grandeurs macroscopiques par les lois de la physique statistique.

La simulation moléculaire

 La simulation moléculaire regroupe plusieurs méthodes différentes telles que la méthode MonteCarlo ou la dynamique moléculaire. Néanmoins, pour toutes ces méthodes de simulation, le système étudié est décrit de manière discrète comme un ensemble de particules de taille plus ou moins grande (atome, groupe d’atomes, molécules ou macro-molécules). Les différentes méthodes de simulation se distinguent alors par le niveau de description du système, qui est défini en fonction de la nature des informations recherchées et de la puissance de calcul disponible. D’autre part, toutes les méthodes de simulation moléculaire sont basées sur les concepts et les lois de la physique statistique, qui permettent de relier les grandeurs microscopiques calculées (énergie du système, positions des particules) aux observables macroscopiques. La simulation moléculaire consiste alors à explorer l’espace des phases du système simulé pour en déterminer les propriétés. 

La physique statistique : du microscopique au macroscopique

 Du point de vue microscopique et à une température suffisamment éloignée du zéro absolu, un fluide ou un solide, peut être considéré comme un système de N particules ponctuelles. A un instant donné, chaque particule d’indice i est caractérisée par sa position q𝑖 et sa quantité de mouvement p𝑖 . L’espace à 6𝑁 dimensions que l’on appelle espace des phases, regroupe alors l’ensemble des états microscopiques accessibles au système, chaque point de cet espace caractérisant un microétat du système. La valeur d’une observable A macroscopique correspondra alors à la moyenne de la grandeur microscopique associée a sur l’ensemble des microétats 𝑗 pondérés par leur probabilité d’occupation appelée probabilité de Boltzmann et notée 𝑃𝑗 : ⟨𝐴⟩ = ∑︁ 𝑗 𝑎𝑗𝑃𝑗 (2.1) La thermodynamique statistique introduite par Gibbs (1902), permet l’étude des systèmes physico-chimiques à travers la définition d’ensembles statistiques. Un ensemble statistique regroupe une infinité de copies d’un même système dans un état microscopique différent. La notion d’ensemble statistique fait alors le lien entre les microétats du sytème et les observables macroscopiques. Chaque ensemble statistique est défini par des conditions thermodynamiques différentes (i.e. conditions expérimentales). L’expression de la probabilité d’occupation d’un microétat du système sera alors différente en fonction de l’ensemble statistique dans lequel on se place pour l’étude du système. L’ensemble canonique (𝑁, 𝑉, 𝑇) décrit par exemple un système fermé, de volume constant en équilibre avec un thermostat. 

L’échantillonnage de l’espace des phases

 Les méthodes de simulation moléculaire sont basées sur la thermodynamique statistique et permettent de calculer les propriétés, observables physiques, d’un système à partir des données microscopiques sur les différentes particules constituant ce système. En effet, selon les principes de la thermodynamique statistique, la mesure d’une observable macroscopique correspond à la moyenne des valeurs prises par son équivalent microscopique sur un grand nombre d’états du système (comme explicité précédemment). Chaque état du système correspondant à un point particulier de l’espace des phases du système. Le calcul de la valeur d’une observable macroscopique revient donc à échantillonner l’espace des phases du système de manière à déterminer un ensemble discret de points qui soient représentatifs du système sur lequel effectuer les moyennes d’ensemble. Pour s’assurer que cet ensemble de points de l’espace des phases correspond à un échantillon représentatif de l’ensemble statistique considéré, une probabilité d’occurrence (i.e. probabilité de Boltzmann associée à l’ensemble statistique de travail) est affectée à chaque point de l’espace des phases de cet échantillon. Toute la difficulté de la simulation moléculaire réside donc dans la génération de cet échantillon de l’espace des phases du système, afin qu’il soit suffisamment représentatif pour bien décrire l’espace des phases du système, tout en gardant un temps de calcul raisonnable. En pratique, il existe deux grandes méthodes d’échantillonnage de l’espace des phases d’un système : la dynamique moléculaire et la méthode Monte-Carlo (figure 2.1). La première méthode est la plus intuitive pour générer les configurations d’un système puisqu’elle consiste à suivre l’évolution temporelle des particules de ce système. Dans ce cas, les grandeurs macroscopiques seront calculées à partir de moyennes temporelles effectuées sur l’ensemble de la trajectoire de la dynamique. La méthode Monte-Carlo est quant à elle une méthode stochastique d’exploration de l’espace des configurations * du système. Des moyennes d’ensemble permettront alors de calculer les valeurs des observables macroscopiques du système. D’après l’hypothèse ergodique, la moyenne d’ensemble est identique à la moyenne temporelle. L’algorithme permettant l’échantillonnage de l’espace des phases est donc différent en fonction de la méthode de simulation choisie pour l’étude du système, mais conduit aux mêmes valeurs moyennes des observables macroscopiques. 

La description des interactions intermoléculaires du système

 La prise en compte des interactions intermoléculaires du système étudié peut être effectuée à un niveau ab initio ou classique. Cependant, la description de ces interactions à un niveau classique ne permet pas de rendre compte de la rupture et de la création de liaisons chimiques. Au contraire, les méthodes quantiques qui considèrent explicitement les électrons, sont capables de décrire les réactions chimiques, mais sont en revanche plus coûteuses en temps de calcul. Je décrirai dans la suite brièvement ces deux approches.

Approche quantique 

La description quantique des interactions entre les particules d’un système contenant 𝐾 noyaux et 𝑁 électrons fait intervenir l’équation de Schrödinger. Pour une configuration (R𝑖) donnée des noyaux et si le hamiltonien du système ne dépend pas explicitement du temps, l’équation de Schrödinger se simplifie et prend la forme suivante : ⎛ ⎝− 1 2 ∑︁ 𝑖 𝛻2 𝑖 + 𝑉^ ext(r1, …, r𝑁 ) +∑︁ 𝑖̸=𝑗 1 |r𝑖 − r𝑗 | ⎞ ⎠ 𝜓(r1, …, r𝑁 ) = 𝐸él 𝜓(r1, …, r𝑁 ) (2.4) où l’opérateur 𝑉^ ext contient le potentiel d’interaction entre les électrons et les noyaux et les potentiels extérieurs (par exemple, un champ électrique ou magnétique extérieur appliqué au — 30 — 2.1 — La simulation moléculaire système). Un moyen de résoudre cette équation consiste à décomposer la fonction d’onde polyélectronique 𝜓(r1, …, r𝑁 ) dans une base de fonctions mono-électroniques {𝜙𝑖} mettant alors en jeu de nombreuses intégrales du type ⟨𝜙𝑗 | 𝐴^ |𝜙𝑖⟩, où 𝐴^ est un opérateur quelconque. Si la résolution analytique d’un tel système est impossible et sa résolution numérique difficile du fait du grand nombre de dimensions impliquées, il existe néanmoins trois approches généralement utilisées : les méthodes de Monte-Carlo quantique (QMC), les méthodes de chimie quantiques (Hartree-Fock et post-Hartree-Fock) et la théorie de la fonctionnelle de la densité. Les méthodes de Monte-Carlo quantique sont basées sur des algorithmes stochastiques pour la résolution de l’équation de Schrödinger et le calcul des intégrales poly-électroniques.Les méthodes traditionnelles de chimie quantique reposent quant à elles sur le formalisme des déterminants de Slater. La corrélation électronique est négligée dans la méthode Hartree-Fock mais est réintroduite dans les méthodes dites post-Hartree-Fock telles que la méthode Moller-Plesset (MP2).[98] Je ne donnerai pas plus d’informations sur les méthodes précédentes et invite le lecteur intéressé à se référer au livre de A. Szabo et N. S. Ostlund.[99] La dernière approche, que je présenterai plus en détail dans la suite du chapitre, est basée sur la théorie de la fonctionnelle de la densité (DFT). 

Approche classique

 La modélisation des interactions intermoléculaires du système étudié peut également se faire en utilisant des potentiels d’interaction analytiques, qui permettent d’approximer l’énergie potentielle 𝑈 d’interaction entre deux atomes du système sans faire intervenir la structure électronique. Par conséquent, le coût de calcul est beaucoup plus faible que dans le cas des méthodes quantiques, et permet alors l’étude de systèmes de plus grande taille, de l’ordre de plusieurs milliers atomes. Néanmoins, ces méthodes ne permettent pas de décrire des systèmes où la structure électronique varie notablement au cours du phénomène étudié comme dans le cas des réactions chimiques (rupture et formation de liaisons, transfert électronique, variation du degré d’oxydation). Les potentiels d’interaction utilisés dépendent de la nature du système étudié et de la grandeur que l’on souhaite calculer. Ils sont généralement choisis additifs, leur forme analytique est la somme de différents termes décrivant chaque type d’interaction (𝐸élec énergie électrostatique, 𝐸ind énergie d’induction, 𝐸disp énergie de dispersion et 𝐸rép l’énergie de répulsion) : 𝑈 = 𝑈élec + 𝑈ind + 𝑈disp + 𝑈rép (2.5) Les paramètres qui interviennent dans ces différents termes d’interaction peuvent être issus de calculs quantiques sur la structure électronique ou être ajustés de manière empirique à partir d’un jeu de propriétés mesurées expérimentalement. La difficulté réside dans l’obtention de paramètres de potentiel qui soient relativement transférables d’un système à un autre. La description du système peut se faire à différents niveaux, on peut ainsi considérer tous les atomes, un centre de force sera alors attribué à chaque atome du système ou bien, attribué un centre de force à un groupement fonctionnel tel que –CH3 ou –CH2 pour se placer à un niveau de description supérieur (dit d’atomes unifiés).