Hybrid Quantum Mechanics / Molecular Mechanics (QM/MM) Approaches - - PowerPoint PPT Presentation
Hybrid Quantum Mechanics / Molecular Mechanics (QM/MM) Approaches -Treatment of the electrostatic QM/MM interface - Mauro Boero Institut de Physique et Chimie des Matriaux de Strasbourg University of Strasbourg - CNRS, F-67034 Strasbourg,
1
Institut de Physique et Chimie des Matériaux de Strasbourg University of Strasbourg - CNRS, F-67034 Strasbourg, France and @Institute of Materials and Systems for Sustainability, Nagoya University - Oshiyama Group, Nagoya Japan
2
Molecules in the gas phase Solids and liquids
good properties reproduced using periodic boundary conditions (PBC)
Complex disordered systems
No periodicity Partitioning of the system: QM/MM
No environment 3
Partitioning the system: shopping list
QM methods
by a classical force field (MM)
parts
QM/MM QM/Interface
4
Partitioning of a QM system into 2 parts A and B: The non-linear (NL) correction
where we use the “nuclear density”: and are unknown.
5
QM/MM description of the subsystem B
Point charge representation of MM atoms at fix MM geometry
Because of the breakdown of the point charge representation: 2-, 3- and 4-body terms are needed:
6
QM description of A & MM description of B We collect all the approximations into the energy term
7
QM MM This term is small when
( and are small), and
no bonds crossing the QM/MM boundary: bonds crossing the QM/MM boundary:
crossing bond
We distinguish the two cases ( is small)
review in : M. Boero, Lect. Notes Phys. 795, pag. 81-98, Springer, Berlin Heidelberg 2010
8
We parameterize the potential using atom centered potentials, i.e. centered on the link atoms.
with Gaussian-type projectors, QM MM where Linking atom The parameters are determined through a fitting procedure. E CORR(RI,[]) drdr'i
*
(r)VI
OECP(r RI,r'RI )i (r')
VI
OECP(r,r')
Ylm(r)pl(r)1
ml l
pl(r')Ylm
* (r')
pl(r) rl exp(r2 /(2 2
2))
E CORR(RI,[])
O.A. von Lilienfeld, I. tavernelli, U. Rothlisberger, J. Chem. Phys. 122, 014113 (2005)
9
According to the Hohenberg-Kohn theorem
QM MM
“ the external potential is determined, within a trivial additive constant , by the electron density . This also determines the ground state wave function and all other electronic properties of the system”
In our QM/MM scheme we optimize the parameters
the QM volume (obtained using a full QM description
function:
O.A. von Lilienfeld, J. Chem. Phys. 122, 014113 (2005)
10
Linking atom
electronic density along z-axis Distance along the z-axis [Å]
z
11
HOOC1-CH3 HOOC1-Dconv (1V) HOOC1-Dopt (DCACP) ESP (CH3/D) ESP (C1) Dipole [D]
1.67 2.87 1.41 0.00
[-0.30 + 3*0.1]
0.32
0.74 0.28 0.54
Example: Acetic acid (Box size: 8 Å, gas phase, 80 Ry PW cutoff)
O.A. von Lilienfeld, I. tavernelli, U. Rothlisberger, J. Chem. Phys. 122, 014113 (2005)
12
i)},{RI}] H DFT[{i(r i)};{RI}] H int[{i(r i)},{RI}] H MM [{RI}]
13
i)},{RI}] H DFT[{i(r i)};{RI}] H int[{i(r i)},{RI}] H MM [{RI}]
i)};{RI}]: QM-part: Hartree and xc interaction O(Nel NG NG)
i)},{RI}]:
QM/MM interface O(N cl NG)
MM part: classical (ff) potential O(N cl N cl) atoms classical
number the is functions set
waves plane
number the is
cl G
N N
14
I 1 Ncl
I 1 Ncl
int 3 3 int
I I I I I I I
15
I 1 Ncl
16
RESP(,RI )
J QM
I NN
1 RI r 2
I NN
[(r)] r RI
3
2
I 1 NN
1
review in : M. Boero, Lect. Notes Phys. 795, pag. 81-98, Springer, Berlin Heidelberg 2010
17
However in all the known cases it is r1 ~ 10-12 a.u. r2 ~ 20-25 a.u.
Only NN < MM atoms in this shell
1 r 2
18
ff
I NN
To avoid incompatibilities due to QM (electronic) vs. classical (point charges) description of the electrostatic, we introduce the screened Coulomb potential
ff
I NN
where
n rn
n1 rn1
is a “covalent” radius of atom I, and n is an integer (n=3)
n
19
potential (ESP) due to the QM charge density seen by the close MM atoms
charge fluctuations have been observed in unrestrained ESP charges during dynamics. Namely, we minimize the norm
D
I QM
J NN
2
D qI H
I QM
2
ESP RESTRAIN
qI
D = qI RESP
20
D
I QM
J NN
2
D qI H
I QM
2
is minimized on the fly during the dynamics. wq = weight parameter to reduce charge fluctuations:
The potential VJ is the Coulomb potential generated on the MM atom J by the electronic density distribution
q
where u(|r - rJ|) is a Coulomb potential modified at short range to avoid spurious over-polarization effects.
21
22
D qI H
I QM
H are the so-called Hirshfeld charges. They are defined as
H
at r RI
at r RK
K
JqI D T J I
J
2
J 1/RI RJ
H
= low index running on QM
NN UQM
23
D 0
JqK D T J K
J
J 0
D
1tK K
J AK J J
J T J J
24
RESP[(r)]
D[(r)]
I QM
DFT[(r)]
The D-RESP potential v RESP(r) E RESP (r) E RESP qI
D I QM
qI
D
(r) while the (extra) forces on the atoms are
FJ rJ E RESP E RESP RJ E RESP qI
D I QM
qI
D
RJ
can be computed in a much more efficient way than the exact DFT potential The corresponding potential on the electrons becomes
25
= center of the multipolar expansion (geometric center of the QM system)
5 , , , 3 , , 3
J I I z y x J I z y x J I
26
ff
I NN
D[(r)]qJ
I QM
J MM (ESP)
J MM(MP)
x,y,z
3 J MM(MP)
r )
,x,y,z
5 J MM(MP)
r )(RI r )...
27
&QMMM COORDINATES INPUT TOPOLOGY ADD HYDROGEN AMBER ARRAYSIZES ... END ARRAYSIZES BOX TOLERANCE BOX WALLS CAPPING CAP HYDROGEN ELECTROSTATIC COUPLING [LONG RANGE] ESPWEIGHT EXCLUSION fGROMOS,LISTg FLEXIBLE WATER [ALL,BONDTYPE] FORCEMATCH ... END FORCEMATCH GROMOS HIRSHFELD [ON,OFF] MAXNN NOSPLIT RCUT NN RCUT MIX RCUT ESP RESTART TRAJECTORY SAMPLE INTERACTING [OFF,DCD] SPLIT TIMINGS UPDATE LIST VERBOSE WRITE LOCALTEMP [STEP fn ltg] &END
28
ELECTROSTATIC COUPLING [LONG RANGE] Section: &QMMM
into account for all classical atoms at a distance R ≤ RCUT_NN from any quantum atom and for all the MM atoms at a distance of RCUT_NN < r ≤ RCUT_MIX and a charge larger than 0.1e (NN atoms).
all MM-atoms with RCUT MIX < r ≤ RCUT ESP are coupled to the QM system by a ESP coupling Hamiltonian (EC atoms).
the rest of the classical atoms is explicitly kept into account via interacting with a multipole expansion for the QM-system up to quadrupolar order. A file named MULTIPOLE is produced.
NN-area and in the EC-area list via the force-field charges.
the quantum system by the force-field charges (mechanical coupling). The files INTERACTING.pdb, TRAJECTORY_INTERACTING, MOVIE_INTERACTING, TRAJ_INT.dcd, and ESP (or some of them) are created. The list of NN and EC atoms is updated every 100 MD steps. This can be changed using the keyword UPDATE LIST. The default values for the cut-offs are RCU_ NN=RCUT_MIX=RCUT_ESP=10 a.u.. These values can be changed by the keywords RCUT_NN, RCUT_MIX, and RCUT_ESP with rnn ≤ rmix ≤ resp. 29
The Redistributed Charge scheme (RC) is a scheme to improve the electrostatic at the linking atom.
It tunes the electrostatic balance between quantum and classical description at the interface.
30
The RC scheme is used in conjunction with hydrogen capping approach Nomenclature for the MM atoms:
M1 M2x M2y M2z M3 M3 M3 Q1 Q2x Q2y Q2z Q3 Q3 Q3 31
The RC scheme is used in conjunction with hydrogen capping approach The MM partial charge of atom M1 is redistributed evenly over 3 point charges q0=qM1/n, where n is the number of M1-M2 bonds. The capping hydrogen (or linking atom, HL) carries no charge. Location of the charges q0:
RM 1 Rq0 RM 1 RM 2 Cq0 0.5
32
The RC method introduces an error in the M1-M2 dipole. The dipole q0R (R=|RM1-RM2|) reduces by a factor of (1 – Cq0). For Cq0 = 0.5, the contribution reduces by 50%. In the RCD method, the values of redistributed charges q0 and of the charges
contributions to the M1−M2 bond dipoles are preserved.
, 2 , 2 q q k M RCD k M q RCD
33
In the RC and RCD schemes, the QM subsystem is polarized by the MM system, but the MM subsystem is not polarized by the QM subsystem, resulting in an unbalanced treatment of the electrostatic interactions between the QM and MM subsystems. The polarized boundary RC (PBRC) and polarized-boundary RCD (PBRCD) schemes make improvements in this regard to the RC and RCD schemes, respectively, by allowing self-consistent mutual polarizations between the QM and MM subsystems near the boundary. In the PBRC and PBRCD schemes, the polarization of the MM subsystem due to the QM electrostatic potential is accomplished by adjusting the MM atomic partial charges in the QM/MM calculations according to the principles
34
First, the MM atomic partial charges are assigned to the MM atoms, and the QM/MM calculation is performed. The electric field generated by the QM subsystem (nuclei and electronic wavefunctions) is computed and imposed on the MM, and a new set of MM atomic partial charges is determined according to the electronegativity equalization and charge conservation. The new set of partial charges replaces the old set of partial charges, and a new QM/MM calculation is performed with the updated partial charges (new external potential for the DFT calculation). convergence ?
STOP Convergence= until the variations in partial charges are smaller than preset thresholds, or until the number of iterations exceeds a preset value. no yes 35
Charge equalization method by Rappé and Goddard.
[A. K. Rappé and w. A. Goddard, W. A. (1991) J. Phys. Chem. 95, 3358-3363 (1991).]
EQ(Q
1,Q2,...,QNcl )
(EI 0 I
0QI 1
2 QI
2JII 0 ) I
1 2 QI
I J
QJJIJ QI d3r (r) r RI
I
EI 0 is the unperturbed reference charge I
0 = E
Q
I0
=1/2(IP - EA) is the first order response (electronegativity) term) Coulomb
(repulsive response
second the is EA)
= Q E =
I0 2 2
II
J JIJ
0 Coulomb interaction between unit charges at I and J (at disrance R =| RI - RJ |)
Consider the classical energy expression:
36
The charge equalization method is based on the minimization of the energy
1,Q2,...,QNcl )
i.e. under the constraints
2 1
I J J IJ J II I I N I
cl
QI
I1 Ncl
This leads to a linear system of equations CD D where
1 Qtot, Di i 0 1 0 for i 2
37
Electronegativity equalization method by Mortier et. al.
[W. J. Mortier, S. K. Ghosh, and S. Shankar, J. Am. Chem. Soc. 108, 4315- 4320 (1986) .]
The existence of a unique chemical potential everywhere in the molecule establishes the electronegativity equalization principle that Mortier write as I (Q
1,Q2,...,QNcl ) (I 0 I ) 2(I 0 I )QI
QJ RIJ
JI
cl
N
2 1
with and are the neutral atom electronegativity and hardness, respectively, QI, and QJ are the charges on atoms I and J, and RIJ, is the internuclear
electronegativity and hardness that arise as a consequence of bonding.
I I
I
38
and are obtained by calibrating through small-molecule calculations and are transferable and consistently usable for calculating charges in any molecule.
I I
39