Simulation de la dynamique moléculaire

Qu'est-ce que la dynamique moléculaire ?

La dynamique moléculaire est une technique de simulation numérique permettant de modéliser le mouvement d'un système de particules (noyaux d'atomes, électrons, ion) au cours du temps dans un certain environnement (température, pression...). Dans le cadre de la mécanique newtonienne, la particule (assimilée à une masse ponctuelle) est soumise à une force (dérivant d'un potentiel d'énergie) qui produit une accélération sur celle-ci. Les équations sont intégrées numériquement en utilisant un pas de temps infinitésimal, garantissant la conservation de l’énergie du système, typiquement de une à deux femtosecondes. La simulation consistera alors à calculer à chaque instant la position, la vitesse et l'accélération de chacune des particules, ainsi que les forces entre particules, en utilisant les résultats obtenus à l'instant précédent. Il est à noter que les échelles de taille et de temps accessibles par la dynamique moléculaire sont de l’ordre de quelques dizaines de nanomètres, et de quelques picosecondes.

Les différents types d'interaction entre particules

Les principales interactions entre particules sont les forces électrostatiques de Coulomb, les forces de van der Waals (qui sont prises en compte dans le potentiel de Lennard-Jones (incluant l'interaction attractive de Van der Waals, prépondérante à longue distance, et la répulsion électronique de Pauli, dominante à courte distance, et empêchant l'interpénétration mutuelle des nuages électroniques de deux atomes)), la liaison hydrogène, ou encore, la liaison halogène, ainsi que les interactions intra-moléculaires (liaisons covalentes, liaisons ioniques, liaisons métalliques).

Différents types d'interactions moléculaires

Deux méthodes de simulation de la dynamique moléculaire : classique, et quantique

Le calcul des forces d'interaction peut se faire selon deux méthodes : la dynamique moléculaire ab initio, où les forces sont calculées à partir des premiers principes de la mécanique quantique (et notamment par l'équation de Schrödinger), et la dynamique moléculaire classique (régie par les lois de la mécanique classique, ou newtonienne), où les forces dérivent d'un potentiel interatomique fixé empiriquement. D'une part, les méthodes dites ab initio procèdent au calcul la structure électronique des particules, et prennent en compte le caractère quantique des interactions mises en jeu. Elles n'incluent aucun paramètre empirique ou semi-empirique dans leurs équations et dérivent directement des principes théoriques, sans inclusion de données expérimentales. D'autre part, les méthodes dites classiques font intervenir des potentiels empiriques, ou champs de forces.

Paramétrisation de champs de forces

Un champ de forces décrit la structure de l'énergie potentielle d'un système de particules, à partir de résultats expérimentaux et de calculs précis en mécanique quantique, par la somme d'un ensemble de potentiels énergétiques paramétrés.

La forme fonctionnelle de base d'un champ de force comprend les termes de liaisons relatifs aux atomes liés par des liaisons covalentes, et les termes d'interactions (dits non-liés ou non-covalents) qui décrivent les interactions à longue portée (électrostatique, Van der Waals, interactions à N corps (entre paires, triplets, quadruplets d'atomes)). Un ensemble typique de paramètres comprend des valeurs pour la masse atomique, le rayon de van der Waals, la charge partielle d'atomes individuels, et des valeurs d'équilibre de longueurs de liaison, de mesures d'angles plans et dièdres pour des ensembles d'atomes liés (paires, triplets, quadruplets), et les valeurs de constante de force pour chaque potentiel.

La communauté scientifique a par la suite développé une multitude de champs de forces, dont certains incluent des modèles pour la polarisabilité dans lesquels les charges des particules sont influencées par celles de leurs voisines (champs de forces polarisables basés sur des dipôles induits, des charges ponctuelles, des multipôles distribués, la densité, la théorie de la polarisation de liaison), mais aussi des champs de forces réactifs, des champs de forces à gros grains...

Tous les champs de forces sont fondés sur de nombreuses approximations et construits à partir de différents types de données expérimentales. En conséquence, ils sont considérés comme empiriques. Il est à souligner que le niveau de complexité de la fonctionnelle d'énergie détermine la précision des champs de forces, ceux-ci étant à leur tour injectés dans l'équation du mouvement de Newton pour le calcul des trajectoires des particules. Ainsi, la paramétrisation de la fonction d'énergie potentielle permet de calculer de manière plus ou moins précise et coûteuse (en termes de ressources de calcul) les grandeurs thermodynamiques fondamentales d'un système moléculaire, comme la densité, l'enthalpie de vaporisation, voire le coefficient de diffusion d’un liquide moléculaire.

La mécanique quantique pour la prédiction précise de la structure électronique

L'équation de Schrödinger décrit l'évolution dans le temps d'une particule massive non relativiste en mécanique quantique. Elle lie la fonction d'onde, l'énergie cinétique des électrons, le potentiel électrostatique entre électrons et noyaux, et le potentiel moyen créé par les autres électrons. La fonction d'onde correspond à la représentation de l'état quantique : la probabilité de trouver une particule au voisinage d'une certaine position à tel instant est proportionnelle au carré du module de la fonction d'onde, aussi appelé densité de probabilité de présence.

La méthode de Hartree–Fock est une méthode de résolution approchée de l'équation de Schrödinger, prenant en compte l'impact de la densité électronique, mais ignorant les contributions d'interaction entre plusieurs corps. La méthode du cluster couplé, basée sur la méthode d'orbitale moléculaire Hartree-Fock, ajoute un terme de correction afin de prendre en compte la corrélation électronique (c'est-à-dire l'interaction entre électrons).

Dans sa forme dépendante du temps associée au hamiltonien de Coulomb, l'équation de Schrödinger prédit la plupart des propriétés d'une molécule, y compris sa forme, c'est-à-dire sa structure tridimensionnelle, et notamment sa structure électronique. Or de nombreuses preuves spectroscopiques ont mis en évidence que l'hamiltonien coulombien ne prenait pas en compte certains termes, qui, bien que plus petits que les énergies cinétiques et coulombiennes, n'en sont pas moins non négligeables. Elles ont ainsi conduit à l'introduction d'un nouveau degré de liberté pour les électrons et noyaux, le spin, souvent assimilé au moment cinétique (qui correspond au moment du vecteur relatif à la quantité de mouvement par rapport à un point). L'hamiltonien moléculaire est l'opérateur hamiltonien de l'énergie des électrons et des noyaux d'une molécule, qui, associé à l'équation de Schrödinger, est à la base du calcul des propriétés des molécules, comme la conductivité, les propriétés optiques et magnétiques, ou encore la réactivité chimique. Il est la somme de cinq classes de termes représentant l'énergie totale du système : les opérateurs d'énergie cinétique pour chaque noyau et chaque électron du système, les opérateurs d'énergie potentielle entre électrons et noyaux (attraction coulombienne) et de répulsions coulombiennes entre électrons, ainsi que les opérateurs d'énergie potentielle de répulsions coulombiennes entre noyaux (énergie de répulsion nucléaire).

Pour quelles applications ?

Issue de la physique du solide et de la chimie quantique, la dynamique moléculaire trouve désormais de nombreuses applications : biologie, biochimie, chimie organique, physique des particules, physique de la matière condensée, sciences des matériaux, mécanique... Elle permet de simuler des systèmes biologiques (peptides, protéines, lipides, anticorps, vaccins) qu'ils soient en solution ou sous forme de cristaux ou de membranes, la réactivité chimique à haute précision en solution pour la chimie organique et bio-organique, le comportement des matériaux, des métaux lourds en solution, des liquides ioniques (catalyse, carburant aérospatial…), des systèmes chargés...

La première application concerne la pharmaceutique avec la mesure de l’activité d’un médicament grâce au calcul de l'énergie libre (ΔG), qui permet de qualifier l'affinité d'une molécule pour une cible thérapeutique. La simulation de la dynamique moléculaire est ainsi utilisée pour la modélisation du pharmacophore (ensemble des groupements fonctionnels disposés selon un arrangement spatial adéquat, assurant la fixation du médicament sur le récepteur), la modélisation de la protéine (modélisation d'homologie, annotation de domaines protéiques, qualité de la structure protéique), l'identification de site de liaison protéique, la prédiction de modes multiples de liaison, l'ingénierie protéique (prédiction d'agrégation de protéine, estimation des surfaces accessibles au solvant et de l'hydropathie), la toxicologie prédictive, la modélisation d'anticorps.

En second lieu, la simulation de la dynamique moléculaire est utile au secteur de l'aérospatial grâce au calcul des propriétés thermodynamiques (stabilité thermique, compressibilité, diffusion...) permettant de caractériser le comportement du carburant sous forme de liquide ionique.

Une troisième application concerne la science des matériaux avec la mise en relation entre la structure atomique d'un matériau et ses propriétés macroscopiques. Les calculs de dynamique moléculaire permettent alors de prédire les propriétés de miscibilité, de mélange, de cohésion, de mouillabilité, de diffusion, d'adhésion aux surfaces, optiques, structurelles, conformationnelles, vibratoires, thermophysiques... Par là-même, ils rendent possible la conception de divers matériaux, polymères (peintures, revêtements, lubrifiants, emballages, gels, adhésifs), catalyseurs, composites, métaux et alliages, céramiques, semi-conducteurs, batteries et piles à combustibles, nanomatériaux...

L'apport de nos chercheurs

Tinker-HP, la simulation moléculaire haute performance

Nos chercheurs ont développé Tinker-HP, un logiciel de dynamique moléculaire intégrant une méthode de résolution massivement parallèle, qui utilise des champs de forces polarisables de nouvelle génération incluant des multipôles électrostatiques distribués. Issu d'une collaboration scientifique pluridisciplinaire entre chimie, informatique et mathématiques, Tinker-HP permet de faire des simulations de dynamique moléculaire en temps long (de l'ordre de la microseconde) sur de grands systèmes moléculaires ayant jusqu'à plusieurs millions d'atomes.

En résumé, Tinker-HP intègre une méthodologie plus rapide que la mécanique quantique (des millions de fois) mais aussi précise, grâce aux champs de forces polarisables de nouvelle génération formulés sur la chimie quantique (multipôles et polarisabilités distribués) incluant la physique à N-corps, ainsi qu'à la résolution massivement parallèle sur des milliers de cœurs (CPU) ou sur des cartes graphiques (GPU). Au-delà de l'approximation par paires, Tinker-HP inclut les effets explicites d'induction multi-corps, comme la polarisation et le transfert de charge, ce qui permet un calcul plus précis des positions, vitesses, accélérations et forces à l’échelle atomique. D'autre part, le couplage multi-échelles quantique/classique (QM/MM) permet le calcul explicite de la structure électronique de sous-systèmes complexes de particules, donnant accès à des phénomènes à l'échelle subatomique. De manière concrète, Tinker-HP rend accessible la réactivité chimique d’un système moléculaire, autant qu'il rend possible la prédiction quantitative des propriétés thermodynamiques des systèmes moléculaires.

Les principaux domaines d'application de la simulation moléculaire

Tinker-HP peut traiter des systèmes biologiques, des métaux lourds en solution (actinides, lanthanides), des liquides ioniques (catalyse, carburant aérospatial…), et des systèmes chargés.

Pharmaceutique

Mesure de l’activité d’un médicament : prédiction quantitative de propriétés thermodynamiques de molécules par le calcul d’affinité d’un ligand pour une cible thérapeutique (énergie libre ΔG)

Spatial

Calcul des propriétés thermodynamiques pour la caractérisation du comportement du carburant sous forme de liquide ionique (stabilité thermique, compressibilité, diffusion...)

Matériaux

Mise en relation entre la structure atomique et les propriétés macroscopiques (miscibilité, cohésion, mouillabilité, diffusion, adhésion aux surfaces, optiques, structurelles, conformationnelles, vibratoires, thermophysiques...)

Conception de polymères (peintures, revêtements, lubrifiants, emballages alimentaires, gels, adhésifs), composites, métaux et alliages, céramiques, semi-conducteurs, batteries et piles à combustibles, nanomatériaux...

Tinker-HP, au cœur de la technique

Cette pluridisciplinarité a notamment mené au développement de nouvelles méthodes de résolution des équations physiques adaptées à une utilisation massivement parallèle des infrastructures de calcul, ce qui a permis une réduction drastique des temps de restitution pour l’utilisateur. Des techniques de décomposition spatiale en 3D et de distribution de mémoire appliquées aux modèles de champs de force polarisables permettent ainsi de traiter les interactions à courte portée. Les interactions électrostatiques et de polarisation à longue portée sont prises en compte par le couplage de la méthode d'Ewald avec des solveurs de polarisation itératifs et non-itératifs de Krylov. Des facteurs d'accélération d'entre quatre et cinq ordres de grandeur en temps ont été constatés sur de grands systèmes par rapport à une algorithmique standard séquentielle. Enfin, il est à noter que la conception du code permet l'utilisation d'une large variété de systèmes de calcul, allant des stations de travail de laboratoire aux superordinateurs de type petascale composés de milliers de cœurs. Tinker-HP propose ainsi le premier environnement de calcul à haute performance pour le développement de la prochaine génération de champs de force polarisables.

Quelques exemples quant au potentiel de Tinker-HP

A titre d'exemple, Tinker-HP a simulé (en 2017) une boite d'eau avec 23,3 millions d'atomes sur un superordinateur composé de 12 288 cœurs sur 512 nœuds. Bien entendu, la productivité de la simulation a été moindre par rapport à des systèmes de tailles plus raisonnables, mais Tinker-HP a été capable de simuler une douzaine de picosecondes par jour, ce qui serait suffisant pour calculer des propriétés telles que les potentiels électrostatiques, ou même pour calculer de petites évolutions de dynamique moléculaire de l'ordre de la nanoseconde. Un autre exemple concerne la simulation d'une protéine solvatée de 200 000 atomes, dont le calcul durait entre 4 et 5 mois avec les approches antérieures alors qu'il ne dure plus que quelques heures avec Tinker-HP (en 2017).