slide1
Download
Skip this Video
Download Presentation
M.F. Ruiz-L ó pez, UMR CNRS-UHP 7565 Université Henri Poincaré, Nancy

Loading in 2 Seconds...

play fullscreen
1 / 76

M.F. Ruiz-L pez, UMR CNRS-UHP 7565 Universit Henri Poincar , Nancy - PowerPoint PPT Presentation


  • 181 Views
  • Uploaded on

GDR « Agrégation, Fragmentation et Thermodynamique des Systèmes Complexes Isolés » Modélisation, Dynamique et Thermodynamique Des Systèmes Complexes Isolés Ecole thématique du CNRS Porquerolles, 22 au 27 mai 2005. Modélisation et dynamique de systèmes. complexes avec un champ de force mixte.

loader
I am the owner, or an agent authorized to act on behalf of the owner, of the copyrighted work described.
capcha
Download Presentation

PowerPoint Slideshow about 'M.F. Ruiz-L pez, UMR CNRS-UHP 7565 Universit Henri Poincar , Nancy' - albert


An Image/Link below is provided (as is) to download presentation

Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author.While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server.


- - - - - - - - - - - - - - - - - - - - - - - - - - E N D - - - - - - - - - - - - - - - - - - - - - - - - - -
Presentation Transcript
slide1

GDR « Agrégation, Fragmentation et Thermodynamique des Systèmes Complexes Isolés »

Modélisation, Dynamique et Thermodynamique Des Systèmes Complexes Isolés

Ecole thématique du CNRS

Porquerolles, 22 au 27 mai 2005

Modélisation et dynamique de systèmes

complexes avec un champ de force mixte

QM/MM ou QM/QM\'

M.F. Ruiz-López, UMR CNRS-UHP 7565

Université Henri Poincaré, Nancy

slide2

I. Historique, aspects fondamentaux, mise en œuvre

  • Problèmes posés par l’étude de gros systèmes moléculaires
  • La méthode QM/MM
  • Traitement des liaisons frontière, la méthode LSCF
  • ONIOM, approches QM/QM
  • Couplage avec la Dynamique Moléculaire
  • Couplage avec le modèle du Continuum
  • Exemple : étude d’une molécule d’eau dans l’eau
  • II. Potentiels d’interaction.
  • Nécessité d ’une paramétrisation
  • Stratégies de paramétrisation
  • Au delà du potentiel de Lennard-Jones
  • (Termes électrostatiques dans les méthodes semiempiriques)
  • III. Simulation de réactions chimiques
  • Difficultés
  • Calcul de profils d ’énergie libre
  • Calcul de trajectoires
slide3

Articles de revue

  • Gao J.,Methods and Applications of Combined Quantum Mechanical and Molecular Mechanical Potentials, Lipkowitz K.B; Boyd D.B., Eds.; Reviews in Computational Chemistry Vol. 7; VCH, 1996; 119-185
  • P. Amara and M. Field Combined Quantum Mechanics and Molecular Mechanics Potentials, Encyclopedia of Computational Chemistry, Vol. 1, Ed. P.v.R. Schleyer, Wiley & sons, 1998, Pg 431-437
  • M. F. Ruiz-López and J.L. Rivail, Combined Quantum Mechanics and Molecular Mechanics Approaches to Chemical and Biochemical Reactivity., Encyclopedia of Computational Chemistry, Vol. 1, Ed. P.v.R. Schleyer, Wiley & sons, 1998, Pg 437-448
  • G. Monard and K.M. Merz Jr. , Combined Quantum Mechanics and Molecular Mechanics Methodologies Applied toBiomolecular Systems, Acc. Chem. Resd. 1999, 32, 904-911
  • THEOCHEM Special Issue « Combined QM/MM calculations in Chemistry and Biochemistry », Vol. 632 (2003)
slide4

Problèmes posés par l’étude de gros systèmes moléculaires

  • Objectifs :
  • Décrire les propriétés électroniques d’un système « complexe »: formation ou rupture des liaisons chimiques, spectroscopies,….
  • Difficultés :
  • Grand nombre d’électrons
  • Grand nombre de degrés de liberté nucléaires
  • Approches « traditionnelles » :
  • Calculs poussés sur un système modèle de petite taille
  • Calculs approchés (semiempiriques) sur un modèle « réaliste »
  • Approches actuelles :
  • Méthodes quantiques à croissance linéaire
  • Dynamique Moléculaire « Ab initio », Car-Parrinello, ….
  • Méthodes mixtes QM/MM, ONIOM, QM/QM’…
slide5

Problèmes posés par l’étude de gros systèmes moléculaires

Comparaison du temps de calcul en fonction de la méthode

(les calculs ab initio et DFT utilisent une base 6-31G*)

slide6

La méthode QM/MM

  • Principe :
  • Processus chimique (plus ou moins) local
  • Il est donc possible de diviser le système complet en deux parties :
  • Sous-système chimique : il est décrit par une approche de mécanique quantique: ab
  • initio, DFT, semiempirique
  • Environnement : il est décrit par un champ de force de mécanique moléculaire

QM System

QM System

MM system

MM system

slide7

La méthode QM/MM

Un modèle quantique/classique pionnier: le modèle du continuum

  • le soluté est décrit par une méthode de Chimie Quantique
  • le solvant est représenté par un milieu continu diélectrique (polarisable)
  • interaction électrostatique prise en compte dans l’hamiltonien du soluté
  • interaction non-électrostatique évaluée par des formules paramétrées
  • D. Rinaldi & J. L. Rivail , Theoret. Chim. Acta, 32 (1973) 57
  • O. Tapia & O. Goscinski, Mol. Phys., 29 (1975) 1653
  • S. Miertus, E. Scrocco, & J. Tomasi, Chem. Phys., 117 (1981) 117

Diélectrique polarisable

Molécule « quantique »

H = H° + Vsoluté-solvant

slide8

La méthode QM/MM

Premiers modèles QM/MM

Niveau semiempirique

Warshel & Levitt J. Mol. Biol., 103 (1976) 227

P. A. Bash, M. J. Field & M. Karplus, JACS, 109 (1987) 8092

M. J. Field, P. A. Bash & M. Karplus, J. Comp. Chem., 11 (1990) 700

Ab initio /MM

Alagona, Desmeules, Ghio, Kollman (JACS 1984)

Singh, Kollman (JCC,1986)

Stanton, Little, Merz (JPC 1995)

Freindorf,Gao (JCC 1996)

DFT/MM

Stanton, Hartsough, Merz (JPC, 1993, 1995), ions en solution

Tuñón, Martins-Costa, Millot, Ruiz-López (JMM, CPL, 1995) molécules en solution

Parrinello et col. (JCP, 1999), ondes planes

Local SCF method

Théry, Rinaldi, Rivail, Maigret, Ferenczy (semiempirical, JCC 1994)

Assfeld, Rivail (ab initio, CPL 1996)

slide9

La méthode QM/MM

HTOTAL= HQM+ HQM / MM+ HMM

HQM Hamiltonien du système QM

HQM Hamiltonien du système MM

MM

QM

HQM / MMHamiltonien d’interaction

On suppose qu’il n’y a pas de liaison chimique entre les sous-systèmes

QM et MM, pour le moment

slide10

La méthode QM/MM

Termes électrostatiques

On suppose que le sous-système MM contient des charges qs et des moments dipolaires permanents °s

s, t MM centres

n QM noyaux

i QM électrons

Les champs électriques EMM et EQM sont définis par

avec

slide11

La méthode QM/MM

Termes de polarisation (cas d’un champ de force MM polarisable)

Supposons maintenant que le sous-système MM est polarisable avec

inds = moments dipolaires induits :

s, t MM centres

n QM noyaux

i QM électrons

(s = polarisabilité isotrope)

Énergie nécessaire à la création du dipole induit (théorie de la réponse linéaire)

La somme des 2 termes conduit à :

slide12

La méthode QM/MM

Termes de van der Waals

Interactions non-électrostatiques. En général = Lennard-Jones. Excellent rapport qualité/prix dans les simulations standard.

with

répulsion

dispersion

Principales limitations du potentiel de Lennard-Jones :

- pas très flexible (2 paramètres seulement)

- indépendant des coordonnées électroniques => correction « post-SCF »

slide13

La méthode QM/MM

Equations Hartree-Fock (Kohn-Sham)

l’opérateur de Fock doit contenir tous les termes dépendant explicitement des coordonnées électroniques.

Modèle non-polarisable

Modèle polarisable

L ’opérateur doit inclure en plus la contribution à l ’énergie de polarisation totale (MM + QM/MM) qui contient de termes quadratiques de la matrice densité. La contribution est :

La résolution des équations SCF et le calcul des moments induits dans la région MM doivent être faits par un processus itératif auto-cohérent. Coût = environ 5 fois celui d ’un modèle non-polarisable !

slide14

La méthode QM/MM

Energie totale

 fonction d ’onde perturbée du sous-système QM

Dérivées premières

Calcul analytique possible (comme dans la méthode QM)

=> L’optimisation de géométrie du système QM/MM est possible

=> forces sur les atomes, dynamique moléculaire

slide15

QM

MM

MM

Traitement des liaisons frontière, la méthode LSCF

Recent reviews

Reuter, N.; Dejaegere, A.; Maigret, B.; Karplus, M. J Chem Phys, 2000, 104, 1720

N. Ferré, X. Assfeld, J.L. Rivail, J. Comp. Chem., 2002, 23, 610

G Monard, X. Prat-Resina, A. González-Lafont, J.M. Lluch, Int. J. Quantum Chem., 2003, 93, 229

slide16

Traitement des liaisons frontière, la méthode LSCF

  • Approche “Link-atom”
  • Singh, Kollman (JCC,1986)
  • M. J. Field, P. A. Bash & M. Karplus (JCC, 1990)
  • On ajoute des atomes H pour saturer la valence des atomes QM frontière.
  • Ces atomes n’interagissent pas avec les atomes MM (pas toujours)
  • Ils sont inclus dans le calcul QM
  • Avantages :
  • Méthode intuitive, simple à mettre en œuvre, réaliste si les atomes frontière sont loin.
  • Inconvénients :
  • L’expression de l ’énergie total et celle de ses dérivées ne sont pas bien définies
  • Rôle très important des charges MM proches
slide17

Traitement des liaisons frontière, la méthode LSCF

Méthodes dérivées du « Link-atom »

Pseudo-halogen

(Hyperchem, 1994)

Cummins, Gready (JPC B 2000)

Pseudo-bond approach

Zhang, Lee, Yang (JCP 1999)

Connexion Atoms

Antes, Thiel (JPC A, 1999)

H  Y (7 electrons)

Parameterized atom

Effective Core potentials

slide18

Traitement des liaisons frontière, la méthode LSCF

Méthodes utilisant des orbitales gelées

LSCF approach (Local SCF)

Théry, Rinaldi, Rivail, Maigret, Ferenczy (semiempirical, JCC 1994)

Assfeld, Rivail (ab initio, CPL 1996)

Ferré, Assfeld, Rivail (JCC 2002)

Fragment SCF

Náray–Szabó, Surján (CPL 1983)

Náray–Szabó (Comp Chem 2000)

GHO approach (Generalized Hybrid Orbital)

Gao, Amara, Alhambra, Field (JPC, 1998)

Frozen orbital QM/MM

Philipp, Friesner (JCC 1999)

slide19

Traitement des liaisons frontière, la méthode LSCF

Localized SCF (LSCF)

Théry, Rinaldi, Rivail, Maigret, Ferenczy (semiempirical, JCC 1994)

Assfeld, Rivail (ab initio, CPL 1996)

Ferré, Assfeld, Rivail, J. Comput Chem., 23, 610 (2002)

  • Principe de la méthode :
  • Décrire les liaisons frontière par des orbitales hybrides strictement
  • localisées et gelées pendant le calcul SCF
  • (SLBOs (= Strictly localized bond orbitals)
  • Pour ce faire, on doit transformer la base d’AOs
  • Les SLBO sont obtenues dans un calcul de référence

Orbitales hybrides SLBO

Base d’orbitales atomiques

slide20

Traitement des liaisons frontière, la méthode LSCF

  • Méthode :
  • 1) Orthogonalisation des SLBO L orbitales gelées orthonormales :
  • 2) Soustraire la projection des orbitales atomiques sur les SLBOs
  • 3) Éliminer les L dépendances linéaires (orthogonalisation canonique)
  • B est rectangulaire Nx(N-L) (équivalente à S-1/2 dans le SCF standard)
  • On obtient (N-L) fonctions m ’ orthogonales entre elles et aux SLBOs
slide21

Traitement des liaisons frontière, la méthode LSCF

Related method : the GHO approach (Generalized Hybrid Orbital)

Gao, Amara, Alhambra, Field (JPC, 1998)

- Semiempirical level

- Use of 4 hybrid orbitals at the MM frontier atom (“hybrid atom”)

- Auxiliary orbitals are not included in the SCF calculation

- Auxiliary orbitals + nucleus charge act as an Effective Core Potential (ECP)

- The ECP is optimized modifying the original semiempirical parameters of

the boundary atom

slide22

ONIOM et méthodes QM/QM

  • Principe de la méthode ONIOM
  • Choix d’un sous-système, le « modèle »
  • 3 calculs:
  • système modèle au niveau de calcul HAUT => E1
  • système modèle au niveau de calcul BAS => E2
  • système complet au niveau de calcul BAS => E3
  • L’énergie du système complet au niveau HAUT, E4,
  • est estimée à l’aide de l’équation :
  • E4 = E3 + (E1 - E2)
  • Liaison frontière : link-atom

Low

level

« modèle »

High

level

IMOMM

F. Maseras and K. Morokuma, J. Comp. Chem. 16 (1995) 1170-1179.

IMOMO

S. Humbel, S. Sieber and K. Morokuma, J. Chem. Phys. 105 (1996) 1959-1967.

ONIOM

M. Svensson, S. Humbel, R. D. J. Froese, T. Matsubara, S. Sieber, and K. Morokuma, J. Phys. Chem., 100, 19357-19363 (1996).

slide23

ONIOM et méthodes QM/QM

Réaction PM3 B3LYP

PhOH  PhO¯ + H+ 331.1 350.9

DMP + H+ DMPH+ -194.7 -238.5

PhOH + DMP  PhO¯ + DMPH+136.4 112.4

E(B3LYP) = E (PM3) - 24 kcal/mol

Équation valable pour d ’autres systèmes proches ?

slide24

ONIOM et méthodes QM/QM

E(B3LYP) = 25.2 + 0.989 E(PM3)

R= 0.999

Correlation between PM3 and B3LYP results for the energies of (p-R-PhO¯ + HDMP+) and (p-R-PhO¯··HDMP+) (relative energies with respect to the separated neutral molecules).

slide26

ONIOM et méthodes QM/QM

Méthode QM/SE-I

Cui, Guo & Karplus, J. Chem. Phys., 117, 5617 (2003))

Système SE (semiempirique)

  • Principe:
  • Approximation NDDO pour toutes les interactions SE et QM/SE
  • Matrice densité QM/SE => Pmn = 0
  • Seules les interactions de Coulomb QM/SE sont donc prises en compte au niveau de l’opérateur de Fock
  • Les interaction non-électrostatiques sont estimées par un potentiel de van der Waals (comme dans QM/MM)
  • Les interactions coulombiennes sont calculées par l’une de ces 2 méthodes :
  • base auxiliaire de gaussiennes pour décrire  dans SE
  • distribution multipolaire pour représenter  dans SE (charges de Mulliken dans la pratique)

Système QM (ab initio ou DFT)

slide27

Core b

Core a

Buffer a

Buffer b

ONIOM et méthodes QM/QM

Méthode basée sur l ’approche Divide & Conquer

V. Gogonea, L.W. Westerhoff and K.M. Merz Jr., J. Chem. Phys., 113, 5604 (2000)

Principe:

T valeur fictive (1000 K, valeur typique)

eF determiné itérativement, Nb total é constant

slide28

ONIOM et méthodes QM/QM

Méthodes basées sur la DFT

Gordon-Kim, JCP, 1972

Principe:

Supposons que l’on connaît la densité électronique des systèmes A et B.

L ’énergie du système AB, EAB est calculée à partir de la densité AB = A +B en utilisant la fonctionnelle de la densité :

EAB = F [AB] = T [AB] + Vne [AB] + Vee [AB]

T kinetic term

Vne Coulomb Z - e-

Vee Coulomb e- - e- + exchange-correlation

Système A

Système B

slide29

ONIOM et méthodes QM/QM

Implémentations :

Frozen-density: Wesolowski, Warshel JPC, 1993

Méthode itérative: Jeanvoine, Ruiz-Lopez, 1994

Extensions possibles pour un calcul QM/QM général

Carter et col., JCP 1999.

B

A

Calcul SCF de A en présence du potentiel crée par B (charges ponctuelles, développement multipolaire,…). Calcul itératif nécessaire si B n’est pas gelée.

Principale difficulté : T[]

Thomas-Fermi CF ∫5/3 (r) dr (LDA)

von Weizsacker (correction 2ème ordre)

slide30

ONIOM et méthodes QM/QM

Energie d’interaction (HCl)2

Y. Jeanvoine, D.E.A., Nancy 1994

Calcul DFT utilisant F[a+b]

Calcul DFT standard

slide31

Couplage avec la dynamique moléculaire

Approximations de base:

Traitement « classique » des noyaux

Born-Oppenheimer

Algorithmes possibles:

1) Calculs QM/MM sur un ensemble de configurations issues d’une simulation de dynamique moléculaire classique

1 simulation MD classique + quelques centaines de calculs QM/MM

Limitations : celles du champ de force (réactions chimiques)

2) Calcul auto-cohérent QM/MM en présence d’un champ électrostatique « moyen » obtenu par un simulation de dynamique classique. (Méthode ASEP = average solvent electrostatic potential)

plusieurs simulations MD classiques + quelques calculs QM

Limitations : pas de dynamique explicite pour le soluté

3) Simulations avec forces QM/MM calculées « on the fly »

1 simulation QM/MM impliquant plusieurs milliers (millions) de calculs QM/MM

slide32

Couplage avec la dynamique moléculaire

  • Simulations QM/MM « on the fly » :
  • Pour une configuration initiale donnée, il faut résoudre les équations SCF du système QM en interaction avec le système MM
  • La fonction d’onde est stockée (« guess » de l’itération suivante)
  • Calcul des forces sur les noyaux QM et MM
  • Résoudre les équations du mouvement, déduire les nouvelles positions des noyaux au temps t + Dt.
  • Itérer le processus pour la nouvelle configuration.
  • Dt : suffisamment petit pour un échantillonnage correct du système, fonction des modes de vibration les plus rapides => typiquement les liaisons XH…. (Dt = 0.2-1.0 fs )
slide33

Couplage avec la dynamique moléculaire

  • Conditions initiales
  • Configuration initiale : peut être déduite d’une simulation purement classique. Cependant, le système doit être re-équilibré car le champ de force change !!!
  • Fonction  d ’essai : à chaque pas de la dynamique, la fonction d ’onde est calculée : le choix d ’une bonne fonction d ’essai est donc très important afin de réduite au minimum le nombre d ’itérations SCF :
    • Option standard : utiliser  calculée dans l’itération précédente
    • Plus efficace : utiliser une technique d’extrapolation (matrice densité, densité électronique)

Percentage of SCF calculations requiring N iterations. Comparison between two simulations of 300 fs using the standard density guess (full line) and an extrapolated guess (dashed line)

slide34

Continuum diélectrique polarisable, constante diélectrique o

MM system

QM system

Couplage avec le modèle du continuum

S. Chalmet, D. Rinaldi and M.F. Ruiz-López Int. J. Quantum Chem., 84, 559 (2001)

Continuum (eo)

Continuum (eo)

Etotal = EQM + EQM/MM + EMM + Esolv

slide35

Exemple : étude d’une molécule d’eau dans l’eau

Première étude DFT/MM d’une molécule en solution

Tuñón et al, J. Mol. Model., 1, 196 (1995)

Modèle

1 molécule QM + 128 molécules MM

(molécule QM rigide, géométrie =TIP3P)

Niveau de calcul

QM : DFT(BP) O(7111/411/1) H(41/1)

MM : champ de force TIP3P

Interactions QM/MM vdW = TIP3P

Simulations MD

Conditions périodiques

Intégration 1 fs

Durée : 70 ps après thermalisation

ensemble NVE

slide37

Exemple : étude d’une molécule d’eau dans l’eau

Internal energy in kcal/mol Ei = 1/2(EQM+QM/MM-E°QM)

Experimental -9.92

Theoretical

Standard TIP3P MD -9.86 0.03

DFT/MM MD -9.72 2.1

Induced dipole moment

« Experimental »= 0.75 D (lower limit)

Theoretical

AIMD Car-Parrinello = 1.08 D

DFT/MM MD = 0.82 D

Chalmet et al, JCP 2001

slide38

Exemple : étude d’une molécule d’eau dans l’eau

Moment dipolaire instantané dans la simulation QM/MM MD

slide39

Exemple : étude d’une molécule d’eau dans l’eau

Comparaison de différents modèles

Modèle 1: Continuum

1 molécule QM + continuum e=80

Modèle 2: QM/MM/Continuum

1 molécule QM + 4 molécules MM TIP3P + continuum

Modèle 3: QM/MM MD

1 molécule QM + 215 molécules MM TIP3P + conditions périodiques + simulation MD 25 ps.

slide40

Exemple : étude d’une molécule d’eau dans l’eau

Potentiel électrostatique créé par le solvant

QM/MM MD

QM/MM/Continuum

Continuum

Moment dipolaire induit Dm

« Expérimental » 0.75 D

Calculs

AIMD Car-Parrinello : 1.08 D

Continuum 0.60 D

QM/MM/Continuum 0.80 D

QM/MM MD 0.82 D

slide41

II. Potentiels d’interaction.

  • Nécessité d ’une paramétrisation
  • Stratégies de paramétrisation
  • Au delà du potentiel de Lennard-Jones
  • Termes électrostatiques dans les méthodes semiempiriques
slide42

Nécessité d ’une paramétrisation

Pour un niveau QM et un champ de force MM donnés, seul l’hamiltonien

non électrostatique HvdWQM/MM n’est pas complètement défini.

s MM sites

n QM nuclei

n , n ?

slide43

Nécessité d ’une paramétrisation

Molécule d ’eau dans l ’eau : simulation DFT/TIP3P MD utilisant des

paramètres de Lennard-Jones standard TIP3P pour les interactions QM/MM

gOO(R)

slide44

Nécessité d ’une paramétrisation

Molécule d ’eau dans l ’eau : simulation AM1/TIP3P MD utilisant des

paramètres de Lennard-Jones standard TIP3P pour les interactions QM/MM

gOO(R)

expérimentale

AM1/TIP3P

slide45

Nécessité d ’une paramétrisation

Donc :

Une paramétrisation est en général nécessaire

Les paramètres n , n devraientdépendre de la méthode QM utilisée

Pour les approches QM ab initio ou DFT, les paramètres peuvent dépendre de la base d ’orbitales

Dans un traitement QM semiempirique, les paramètres peuvent dépendre de la méthode (AM1, PM3, ....)

slide46

Stratégies de paramétrisation

  • Paramétrisation :
  • Que veut-on reproduire ?
  • On peut envisager plusieurs cas :
  • Type de données : reproduire des calculs ab initio poussés, des calculs au même niveau QM ou encore des données expérimentales
  • Type de propriété : reproduire la structure (ou d’autres propriétés, ex surface d’énergie potentielle) d’un complexe, d’un liquide, d’une solution, etc
slide47

Stratégies de paramétrisation

  • Paramétrisation à partir de calculs « ab initio » sur des complexes de référence
  • Recette :
  • calculer la surface d’énergie potentielle ab initio du complexe AB
  • calculer l’énergie d’interaction QM/MM électrostatique
  • la différence représente l’interaction non-électrostatique
  • faire un « fit » de cette énergie à partir de l’expression de Lennard-Jones, c’est-à-dire, minimiser la quantité :
slide48

Stratégies de paramétrisation

Paramétrisation des interactions non électrostatiques NH3 (DFT) - H2O (TIP3P)

QMMM

QMMM

slide49

QM/MM, vdW TIP3P

classical

experimental

QM/MM

vdWdimer optimized

Stratégies de paramétrisation

Paramétrisation pour des simulations en phase condensée.

la stratégie précédente n’est pas valable : effets non-additifs !

Exemple

Tu &Laaksonen (JCP 1999)

Eau liquide

Calculs QM/MM au niveau

HF/6-311G(d,p)/TIP3P

=> Résultats très mauvais lorsque la simulation est faite en utilisant les paramètres LJ optimisés pour le dimère de l’eau.

Moment dipolaire induit :

0.35 D (exp 0.75 D)

Energie interne

-6.74 kcal/mol (exp -9.92 kcal/mol)

slide50

Stratégies de paramétrisation

  • Paramétrisation pour des simulations en phase condensée
  • Martín et al, Chem. Phys. 284 (2001) 607
  • La paramétrisation doit se faire en ajustant une propriété du liquide (ou de la solution)
  • => processus itératif (plusieurs simulations QM/MM !)
  • On peut remplacer les simulations QM/MM par des simulations ASEP, moins coûteuses, (ou par une simulation classique avec des charges atomiques appropriées)
  • Les différences avec les paramètres du champ de force classique peuvent être significatives

O(solute)-O(solvent) radial

distribution function of H2O

Experimental (full line)

ASEP/MD (dashed line)

DFT/MM/MD (dotted line).

slide51

Au delà du potentiel de Lennard-Jones

Limitations du potentiel de Lennard-Jones

1) peu flexible : 2 paramètres

répulsion

dispersion

slide52

Au delà du potentiel de Lennard-Jones

Paramétrisation des interactions NH3 (DFT) - H2O (TIP3P)

Solid : Ab initio, reference calculation

Dot : DFT/MM electrostatic energy

Dashed : Non electrostatic term, to be fitted

Solid : Ab initio, reference calculation

Dot : DFT/MM with optimized Lennard-Jones potential

slide53

Au delà du potentiel de Lennard-Jones

Amélioration en utilisant un potentiel du type : exp-6

Plus flexible ( 3 paramètres)

Meilleur comportement à courte distance (répulsion)

Mais

Incorrect lorsqueR0 (sans incidence dans la pratique, le terme électrostatique n’est pas correct non plus !)

repulsion

dispersion

slide54

Au delà du potentiel de Lennard-Jones

Fitting the interaction energy of the NH3 ···H2O complex

Solid : Ab initio, reference calculation

Dot : DFT/MM with optimized non-electrostatic potential

exp-6

Lennard-Jones

slide55

Au delà du potentiel de Lennard-Jones

Limitations du potentiel de Lennard-Jones

2) Artefacts liés à la polarisation électronique, en particulier avec des bases très étendues

slide56

Au delà du potentiel de Lennard-Jones

Limitations du potentiel de Lennard-Jones

2) Artefacts liés à la polarisation électronique, en particulier avec des bases très étendues

Potentiel pseudo-6

Le terme de répulsion est défini par une somme de pseudopotentiels centrés sur les atomes classiques

Très flexible (2*k paramètres par atome MM)

Pas de paramètres pour les atomes QM ! Transférable !!!

(r est une coordonnée électronique)

slide57

III. Simulation de réactions chimiques

  • Difficultés
  • Calcul de profils d ’énergie libre
  • Calcul de trajectoires
slide58

G

Coordonnée de réaction

Difficultés

  • Etude théorique de réactions chimiques:
  • décrire les espèces chimiques impliquées (intermédiaires, TS) => mécanisme
  • estimer les variations d’énergie libre mises en jeu => constantescinétiques
  • aspects dynamiques : coefficient de transmission => vitesse

Théorie de l ’état de transition

Coefficient de transmission

slide59

Difficultés

  • Systèmes complexes, particularités :
  • grand nombre de degrés de liberté => moyennes statistiques
  • réaction => processus plus ou moins local => justifie la séparation QM/MM
  • considérer les échelles de temps (réaction chimique vs relaxation du milieu)
slide60

Calcul de profils d ’énergie libre

  • Approche statique :
  • Localisation de points stationnaires dans la surface d’énergie potentielle (sous-système chimique)
  • L ’optimisation du système complet est parfois possible. Sinon (notamment pour localiser les TS, points selle), 2 possibilités :
    • Environnement gelé ou
    • Environnement relaxé => micro - itérations
  • Calcul d’énergie libre comparable au
  • calcul standard pour un système isolé
  • - ZPE + contribution thermiques
  • - Intrinsic Reaction Coordinate

Lluch et col, Theochem, Vol 632 (2003)

slide61

Calcul de profils d ’énergie libre

  • References
  • GRACE, A. J. Turner, V. Moliner, and I. H. Williams. Phys Chem Chem Phys, 1,1323, 1999.
  • Lluch et col, Theochem, Vol 632 (2003)
  • Monard et al, Int. J. Quantum Chem., 93, 229 (2000).
  • I. H. Hillier, (Theochem) 463, 45 (1999)
  • Zhang, Liu, and Yang JCP 112, 3483 (2000)
  • Billeter, Turner, Thiel, Phys Chem Chem Phys, 2, 2177 (2000)
  • Frisch et col, J comput Chem, 24, 760 (2003)
slide62

Calcul de profils d ’énergie libre

  • Approche dynamique :
  • Intégration de l’énergie libre (simulations) Le profil d’énergie libre se calcule en utilisant les techniques de simulation habituelles et le champ de force QM/MM : « FEP » (free energy perturbation), umbrella sampling, …
  • En général, on fait une hypothèse sur le chemin de réaction (1 coordonnée interne ou une combinaison de coordonnées ….) => Choix critique !
  • Le chemin pour le système « isolé » (en absence d ’environnement) est un mauvais choix en général ! On doit plutôt envisager :
    • une approche de type« champ moyen »(continuum, ASEP, …)
    • ou la méthode de micro - itérations
slide63

Exemple de modification du mécanisme d’une réaction par l’effet du solvant

Synthesis of b-Lactams from Fluoroketenes and Imines: Ab Initio Potential Energy Surfaces in Gas-phase and in Solution. R. López, et al J. Phys. Chem., 100, 10600 (1996)

Calcul de l’IRC ( Intrinsic Reaction Coordinate ) dans un continuum

Phase gazeuse

Solution de dichlorométhane (e=9)

Product

Product

slide64

Calcul de trajectoires

  • 1er cas : La barrière d’activation du processus est de l’ordre de kT
  • Les trajectoires peuvent être obtenues dans une simulation directe partant des réactifs. Plusieurs dizaines ou centaines de trajectoires sont nécessaires pour un calcul des moyennes statistiques.
  • 2ème cas : La barrière d’activation du processus >> kT
  • La probabilité d’atteindre le TS et les produits est très faible : la simulation
  • ne peut pas démarrer au niveau des réactifs. Une technique adaptée à ce type de problèmes (rare event) doit être employée.
slide65

Calcul de trajectoires

Exemple : transfert de proton dans le système modèleH2O ···OH¯

Simulation DFT/MM dans l’eau liquide (216 molécules classiques TIP3P)

dOO = 2.9 Å

Boîte 18.8 Å

T=25°C

40 ps + 6 ps

G ≈ 1 kcal/mol

PT ≈ 1 ps

Evolution temporaire de q et distribution de probabilité (PT environ 20-30 fs, fréquence environ 1 ps)

Tuñón et al JCP 106 (1997) 2633

slide66

Calcul de trajectoires

Exemple : transfert de proton dans le système modèleH2O ···OH¯

Simulation DFT/MM dans l’eau liquide (216 molécules classiques TIP3P)

dOO = 2.9 Å

Boîte 18.8 Å

T=25°C

40 ps + 6 ps

DG ≈ 1 kcal/mol

PT ≈ 1 ps

Comparaison de q(t) et du champ électrique crée par le solvant E(t) (délai du solvant environ 50 fs)

Tuñón et al JCP 106 (1997) 2633

slide67

Calcul de trajectoires

Exemple : transfert de proton dans le système modèleNH4+ ···NH3

Simulation AM1/MM dans l’eau liquide (216 molécules classiques TIP3P)

dNN = 2.62 Å

Boîte 18.68 Å

T=25°C

100 ps + 500 ps

DG ≈ 5 kcal/mol

PT ≈ 10-20 ps

Comparaison de q(t) et du champ électrique crée par le solvant E(t)

(Valeurs moyennes, rôle des fluctuations du solvant)

q(t)

E(t)

Li et al, Chem. Phys. 240 (1999) 93

slide68

Calcul de trajectoires

  • 2ème cas : La barrière d’activation du processus >> kT => Rare event approach
  • Le principe consiste à à démarrer les trajectoires à partir du TS :
  • Déterminer la structure de transition TS à l’aide d’une méthode de type “champ moyen”
  • (ASEP, QM/MM/Continuum, Continuum) raffinée éventuellement par simulation.
  • Equilibrer ce TS en solution (la coordonnée de réaction doit être « gelée », les autres coordonnées du système chimique et celles du solvant sont libres; une simulation de type umbrella sampling est aussi envisageable)
  • A partir de cette simulation, sélectionner un ensemble de configurations soluté-solvant indépendantes (Dt > 500 fs, vitesses selon la coordonnée de réaction aléatoires)
  • Chaque configuration sert de point initial d’une dynamique moléculaire sans contrainte.
  • Intégrer les équations du mouvement pour t > 0 puis t < 0.
  • Les 2 demi-trajectoires sont assemblées en une seule trajectoire allant de - t à + t
  • Calculer les propriétés moyennes en utilisant ces trajectoires.
slide69

Calcul de trajectoires

  • Chaque demi-trajectoire peut conduire aux réactifs ou aux produits (pas d’hypothèse de départ)
  • Chaque demi-trajectoire peut éventuellement rebrousser chemin et croiser le TS une ou plusieurs fois
  • Les simulations conduisent donc à des trajectoires réactives (reliant les réactifs aux produits) et non-réactives (reliant les réactifs aux réactifs ou les produits aux produits).
  • les trajectoires non réactives traduisent la présence des effets dynamiques de l ’environnement et contribuent à diminuer le coefficient de transmission du processus
slide70

Calcul de trajectoires

Illustration 1. Simulation DFT/MM - MD de la réaction de bromation de l ’éthylène en solution aqueuse CH2=CH2 + Br2  CH2Br–CH2Br

TS solvaté: projection sur le plan de Br-

slide71

Calcul de trajectoires

Illustration 1.

CH2=CH2 + Br2  CH2Br–CH2Br

Trajectoire réactive « typique »

slide73

Calcul de trajectoires

Illustration 2. Simulation DFT/MM MD de l ’hydrolyse de la liaison amide dans l ’eau.

R2CONHR1 + H2O  R1NH2 + R2COOH

TS for non-assisted concerted mechanism

slide74

Calcul de trajectoires

Haloalkane déshalogénase de Xanthobacter Autotrophicus Réaction SN2 entre le substrat dichloroéthane

et le groupe carboxylate de Asp124.

4851 atomes + eau = total number 17109

Soriano et al, JACS, 2004

slide75

Calcul de trajectoires

Haloalkane déshalogénase de Xanthobacter Autotrophicus Réaction SN2 entre le substrat dichloroéthane

et le groupe carboxylate de Asp124.

4851 atomes + eau = total number 17109

Soriano et al, JACS, 2004

slide76

Conclusion, perspectives

  • Approches QM/MM : outil très intéressant pour l ’étude statique et dynamique de très gros systèmes mais nécessitant, selon les cas, un certain nombre d’adaptations spécifiques (paramétrisation).
  • Approches QM/QM : pas encore des méthodes performantes, sauf ONIOM, très pratique dans son domaine d ’utilisation
  • Avenir proche :
  • Méthodes mixtes « à la carte », méthodes à plusieurs « couches »
  • Méthodes QM/MM et QM/QM avec « Linear Scaling »
  • Méthodes de dynamique moléculaire avec des approches QM « rapides »
  • (type AIMD ou Car-Parrinello)
  • => Des méthodes QM « bon marché » sont nécessaires pour pouvoir aborder l ’étude des systèmes très complexes et/ou et des temps de simulations longs :
  • Nouvelle génération des méthodes semiempiriques
  • Approches « tight binding » basés sur la DFT,
  • (Suhai et col., Phys Rev, B 58 (1998) 7260)
ad