4/26/2005 CHEM 408 - Sp05L15 - 1Molecular Dynamics:MD/MC Common ElementsMD/MC Common Elements• representative phase-space trajectory via F=ma • energies & forces required• dynamical observables available• most natural in microcanonical (N, V, E) ensemble (but others possible)Monte Carlo:• non-dynamical sampling of phase-space via random walk• forces not required• no true dynamics available• canonical ensemble (N, V, T) common, but many others easily accessible (N, P, T), (µ, V, T), …• phase equilibrium more easily modeled(Leach, Ch. 6)Despite these differences, MD & MC methods share many common features4/26/2005 CHEM 408 - Sp05L15 - 21. Boundary Conditionsno boundaries (vacuum)dielectric or other cavity• some early simulations (esp. proteins) simulated free molecules in vacuum or with various types of artificial boundaries• poor choice for most condensed phase simulations because of large fraction of molecules affected by vacuum or artificial walls• consider cubic sample of N = n3molecules; if surface consists of skin of m mols. thickness then the surface fraction is:256 water molecules333)2(nmnnfsurface−−=nNf(m=2) f (m=1)5 125 99% 78%10 1000 78% 49%20 8000 49% 27%40 64000 27% 14%100 1000000 12% 6%4/26/2005 CHEM 408 - Sp05L15 - 3Periodic Boundary ConditionsFive Basic Geometries in 3d:Other BCs for Surfaces, etc.(5th is parallelepiped)Figures from Leach4/26/2005 CHEM 408 - Sp05L15 - 4• periodic BCs are almost universally applied• most often cubic BCs because of simplicity; truncated octahadron also used, esp. for large solutes• some potential problems:- inhibit fluctuations of size > box length- introduce artificial periodicity into systembut these effects are seldom of much consequence4/26/2005 CHEM 408 - Sp05L15 - 52. Potential Truncation• in most cases it is not useful or necessary to calculate the interactions between molecules in the entire periodic array defined by the periodic BCs• in the “minimum image” convention one only considers interactions of a given molecule with its N-1 closest neighborsFigure from Allen & Tildesley, Computer Simulation of Liquids (Oxford, 1987)• most often, fewer interactions are calculated by adopting a spherical cutoff, i.e. neglecting interactions beyond some Rc≤L/2• for N molecules there are N(N-1)/2 MI interactions (still costly)RcMI• reduces the N2scaling to N1for large N4/26/2005 CHEM 408 - Sp05L15 - 6• rather than checking distances each step, a (Verlet) neighbor list, which is updated much less frequently is used to keep track of which pairs of molecules should interact Figure from Allen & Tildesley, Computer Simulation of Liquids (Oxford, 1987)RcRlist• simple truncation can lead to problems in MD simulations where forces are important• one solution is to use shifted potentials:r / σ1.01.52.02.53.0v(r)/ε-1.0-0.8-0.6-0.4-0.20.00.2LJtruncatedshiftedRcmore realistic is Rc≥ 2.5σExample LJ truncation⎩⎨⎧≥<−=cccshiftedRrRrRvrvrvfor 0for )()()(4/26/2005 CHEM 408 - Sp05L15 - 7⎪⎩⎪⎨⎧≥<<+−≤=uullrrrrrrrrrrrSfor 0for )(6-)(15 )(101for 1)(5433ρρρ⎪⎪⎩⎪⎪⎨⎧≥<<−−≤=uulluulrrrrrrrrrrrrSfor 0for )()(for 1)(1⎪⎪⎩⎪⎪⎨⎧≥<<−−−≤=uullululrrrrrrrrrrrrrSfor 0for )()23(for 1)(32⎟⎟⎠⎞⎜⎜⎝⎛−−≡lulrrrrr)(ρ• a much better solution is to use a switching function S(r) to smoothly taper the potential to zero over some small range of r• some examples:(r-rl)/(ru-rl)0.00.20.40.60.81.0S(r)0.00.20.40.60.81.0#1#2#30V(r)S(r) effect4/26/2005 CHEM 408 - Sp05L15 - 8• electrostatic interactions more difficultFigures from LeachExample: Water Dimer Interaction(b) atom-based truncation(a) full V(r)@ 8 ÅV = -.27 kcal/molRc= 8Å• truncation must be done on neutral units:4/26/2005 CHEM 408 - Sp05L15 - 93. Long-Range Forces• interactions which decay as r-d(d=dimensionality) or more slowly are considered long-range interactions (net is Ûv(r)rd-1dr)• in 3d ion-ion (r-1), ion-dipole (r-2) and dipole-dipole (r-3) are considered long-range interactions• such interactions have important contributions from beyond the simulation box or potential cutoff range and, especially in the case of charged species, simple truncation can lead to incorrect results example of Cl-Cl interactiondiffusion constant?4/26/2005 CHEM 408 - Sp05L15 - 10• beyond a cutoff distance Rcthe system is treated like a homogeneous dielectric medium with dielectric constant εs εsRcE1• dipoles within cutoff region induce a polarization of the continuum region which creates a “reaction field” at molecule #1: ∑∈+−=cVjjcssREµεεrr31112)1(2• the effect is to modify interactions between dipolar molecules to: jicssdirectijRvvµµεεrr⋅+−−=3112)1(• the reaction field method is inexpensive and easy to implement; for MD simulations switching functions are needed; dynamics of continuum not included a) Reaction Field Method• useful for systems without charged speciesFigure from Allen & Tildesley, Computer Simulation of Liquids (Oxford, 1987)4/26/2005 CHEM 408 - Sp05L15 - 11b) Ewald Method• calculate all interactions with infinite lattice implied by periodic BCs using ideas developed for treating ionic crystals • consider charge version (also versions for dipoles & other multipoles) ∑∑>=−=NijNijijielectrrqqV10||41rrπε• within primary box • for periodic lattice ∑∑∑∞===+−=011'0||2141rrrrrnNjNijijielectnrrqqVπεlattice vector)0,1,2(=nrnr0r(exclude j0=i0)• Σ is only conditionally convergent (must group terms to get sensible results); consider spherical arrays of cells4/26/2005 CHEM 408 - Sp05L15 - 12• in the Ewald method the direct sum is made short ranged by adding a neutralizing charge distribution around each charge: )exp()(222/33rqriiαπαρ−−=rqiρi•ρi(r) is then subtracted in a 2nd summation performed in reciprocal “k” space −ρir-spacesumk-spacesum)()exp()( rrkrdkiirrrrrρρ⋅−=∫• the result is: r-space sumr-space sumk-space sumL = cell lengthself term and surface term (for array in