A1263195019WYPOz

1
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
2

• 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 dune molécule deau dans leau
• II. Potentiels dinteraction.
• 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

3

• 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)

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

• Objectifs
• Décrire les propriétés électroniques dun 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

5
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)
6
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
7
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
  lhamiltonien 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
8
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)
9
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 / MM Hamiltonien dinteraction
On suppose quil ny a pas de liaison chimique
entre les sous-systèmes QM et MM, pour le moment
10
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
11
La méthode QM/MM
Termes de polarisation (cas dun 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 à
12
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 gt correction
 post-SCF 
13
La méthode QM/MM
Equations Hartree-Fock (Kohn-Sham) lopé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 !
14
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) gt Loptimisation de géométrie du système
QM/MM est possible gt forces sur les atomes,
dynamique moléculaire
15
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
16
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 ninteragissent 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
• Lexpression 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

17
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
18
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áraySzabó, Surján (CPL 1983) NáraySzabó
(Comp Chem 2000) GHO approach (Generalized
Hybrid Orbital) Gao, Amara, Alhambra, Field
(JPC, 1998) Frozen orbital QM/MM Philipp,
Friesner (JCC 1999)
19
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
  dAOs
• Les SLBO sont obtenues dans un calcul de
  référence

Orbitales hybrides SLBO
Base dorbitales atomiques
20
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

21
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
22
ONIOM et méthodes QM/QM

• Principe de la méthode ONIOM
• Choix dun sous-système, le  modèle 
• 3 calculs
• système modèle au niveau de calcul HAUT gt E1
• système modèle au niveau de calcul BAS gt E2
• système complet au niveau de calcul BAS gt E3
• Lénergie du système complet au niveau HAUT,
  E4,
• est estimée à laide 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).
23
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 ?
24
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-PhOHDMP) (relative energies with respect
to the separated neutral molecules).
25
ONIOM et méthodes QM/QM
26
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 gt Pmn 0
• Seules les interactions de Coulomb QM/SE sont
  donc prises en compte au niveau de lopé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 lune 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)
27
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
28
ONIOM et méthodes QM/QM
Méthodes basées sur la DFT Gordon-Kim, JCP, 1972
Principe Supposons que lon 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
29
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
nest pas gelée.
Principale difficulté T? Thomas-Fermi CF
? ? 5/3 (r) dr (LDA) von Weizsacker
(correction 2ème ordre)
30
ONIOM et méthodes QM/QM
Energie dinteraction (HCl)2
Y. Jeanvoine, D.E.A., Nancy 1994
Calcul DFT utilisant F?a?b
Calcul DFT standard
31
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 dune 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 dun 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
32
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 donde est stockée ( guess  de
  lité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 gt typiquement les
  liaisons XH. (Dt 0.2-1.0 fs )

33
Couplage avec la dynamique moléculaire

• Conditions initiales
• Configuration initiale peut être déduite dune
  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
  litération précédente
• Plus efficace utiliser une technique
  dextrapolation (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)
34
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
35
Exemple étude dune molécule deau dans leau
Première étude DFT/MM dune 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
36
Exemple étude dune molécule deau dans leau
(c) O(QM) - H (MM) (d) Q(MM) - H (QM)
37
Exemple étude dune molécule deau dans leau
Internal energy in kcal/mol Ei
1/2(EQMQM/MM-EQM) Experimental -9.92 Theoretic
al 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
38
Exemple étude dune molécule deau dans leau
Moment dipolaire instantané dans la simulation
QM/MM MD
39
Exemple étude dune molécule deau dans leau
Comparaison de différents modèles
Modèle 1 Continuum 1 molécule QM continuum
e80
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.
40
Exemple étude dune molécule deau dans leau
Potentiel électrostatique créé par le solvant
QM/MM/Continuum
QM/MM MD
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
41

• II. Potentiels dinteraction.
• 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

42
Nécessité d une paramétrisation
Pour un niveau QM et un champ de force MM donnés,
seul lhamiltonien non électrostatique HvdWQM/MM
nest pas complètement défini.
s MM sites n QM nuclei
?n , ?n ?
43
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)
44
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
45
Nécessité d une paramétrisation
Donc Une paramétrisation est en général
nécessaire Les paramètres ?n , ?n devraient
dé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, ....)
46
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
  dautres propriétés, ex surface dénergie
  potentielle) dun complexe, dun liquide, dune
  solution, etc

47
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 dinteraction QM/MM
  électrostatique
• la différence représente linteraction
  non-électrostatique
• faire un  fit  de cette énergie à partir de
  lexpression de Lennard-Jones, cest-à-dire,
  minimiser la quantité

48
Stratégies de paramétrisation
Paramétrisation des interactions non
électrostatiques NH3 (DFT) - H2O (TIP3P)
QM MM
QM MM
49
Stratégies de paramétrisation
Paramétrisation pour des simulations en phase
condensée. la stratégie précédente nest pas
valable effets non-additifs !
Exemple Tu Laaksonen (JCP 1999) Eau
liquide Calculs QM/MM au niveau HF/6-311G(d,p)/TIP
3P gt Résultats très mauvais lorsque la
simulation est faite en utilisant les paramètres
LJ optimisés pour le dimère de leau. Moment
dipolaire induit 0.35 D (exp 0.75 D) Energie
interne -6.74 kcal/mol (exp -9.92 kcal/mol)
50
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)
• gt 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).
51
Au delà du potentiel de Lennard-Jones
Limitations du potentiel de Lennard-Jones 1) peu
flexible 2 paramètres
répulsion
dispersion
52
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
53
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 In
correct lorsque R?0 (sans incidence dans la
pratique, le terme électrostatique nest pas
correct non plus !)
repulsion
dispersion
54
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
55
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
56
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
(2k paramètres par atome MM) Pas de paramètres
pour les atomes QM ! Transférable !!!
(r est une coordonnée électronique)
57

• III. Simulation de réactions chimiques
• Difficultés
• Calcul de profils d énergie libre
• Calcul de trajectoires

58
Difficultés

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

Théorie de l état de transition
Coefficient de transmission
59
Difficultés

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

60
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é gt 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)
61
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)

62
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 .) gt 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

63
Exemple de modification du mécanisme dune
réaction par leffet 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 lIRC ( Intrinsic Reaction Coordinate )
dans un continuum
Phase gazeuse
Solution de dichlorométhane (e9)
Product
Product
64
Calcul de trajectoires

• 1er cas La barrière dactivation du processus
  est de lordre 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 dactivation du processus
  gtgt kT
• La probabilité datteindre 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.

65
Calcul de trajectoires
Exemple transfert de proton dans le système
modèle H2O OH Simulation DFT/MM dans leau
liquide (216 molécules classiques TIP3P)
dOO 2.9 Å Boîte 18.8 Å T25C 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
66
Calcul de trajectoires
Exemple transfert de proton dans le système
modèle H2O OH Simulation DFT/MM dans leau
liquide (216 molécules classiques TIP3P)
dOO 2.9 Å Boîte 18.8 Å T25C 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
67
Calcul de trajectoires
Exemple transfert de proton dans le système
modèle NH4 NH3 Simulation AM1/MM dans leau
liquide (216 molécules classiques TIP3P)
dNN 2.62 Å Boîte 18.68 Å T25C 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
68
Calcul de trajectoires

• 2ème cas La barrière dactivation du processus
  gtgt kT gt Rare event approach
• Le principe consiste à à démarrer les
  trajectoires à partir du TS
• Déterminer la structure de transition TS à
  laide dune 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 gt 500 fs, vitesses selon la
  coordonnée de réaction aléatoires)
• Chaque configuration sert de point initial dune
  dynamique moléculaire sans contrainte.
• Intégrer les équations du mouvement pour t gt 0
  puis t lt 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.

69
Calcul de trajectoires

• Chaque demi-trajectoire peut conduire aux
  réactifs ou aux produits (pas dhypothè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

70
Calcul de trajectoires
Illustration 1. Simulation DFT/MM - MD de la
réaction de bromation de l éthylène en solution
aqueuse CH2CH2 Br2 ? CH2BrCH2Br
TS solvaté projection sur le plan de Br-
71
Calcul de trajectoires
Illustration 1. CH2CH2 Br2 ?
CH2BrCH2Br Trajectoire réactive  typique 
72
Calcul de trajectoires
73
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
74
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
75
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
76

• 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 dadaptations 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)
• gt 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)