implicit water is only included implicitly as an effective (averaged) interaction with no atomic detail vs. explicit solvation an explicit representation of water is included at the molecular level What does solvent do? implicit water is only included implicitly as an effective (averaged)

interaction with no atomic detail vs. explicit solvation an explicit representation of water is included at the molecular level What does solvent do? viscous damping force dielectric screening / polarization hydrophobic effect structure implicit or no solvent What does solvent do? viscous damping force

dielectric screening U dielectric constant, = 40 (not important) qi qj rij gas phase liquid phase ?

What does solvent do? viscous damping force dielectric screening (not important) U qi qj rij gas phase liquid phase ? This is BULK solvent screening. At short range,

no screening 80 vs. 1 distance What does solvent do? viscous damping force dielectric screening (not important) U

Simple approximation: qi qj Distance dependent dielectric rij rij or 4rij 80 1 distance What does solvent do? viscous damping force dielectric screening

U 80 (not important) Simple approximation: qi qj Distance dependent dielectric rij rij or 4rij Better? Sigmoidal dielectric

1 What does solvent do? viscous damping force dielectric screening hydrophobic effect (not important) surface area approaches Observation: Gsolvation for the saturated hydrocarbons in water is linearly related to the solvent accessible surface area the line passing

through the nonpolar Ala, Val, Leu and Phe has a slope of 22 cal/2 Creighton Proteins (FM Richards, Ann Rev Biophys Bioeng 6, 151-176 (1977)) surface area approaches Observation: Gsolvation for the saturated hydrocarbons in water is linearly related to the solvent accessible surface area exposed solvent accessible surface area (SASA) Gresidue = iAi atoms,i

free energy of interaction of a solute with water atomic solvation parameters based on free energies of transfer Problems: sensitive to is, parameterization, surface area and change in conformation in dynamics you need derivatives of SASA what about polarization effects? what is the effective solvent polarization? (solve Poisson equation) Born: isolated point charge (q) in a spherical cavity of

radius r immersed in a dielectric continuum with dielectric constant ? GBorn outside r inside q what is the effective solvent polarization? (solve Poisson equation) Born: isolated point charge (q) in a spherical cavity of radius r immersed in a dielectric continuum with

dielectric constant GBorn q2 1 1 2r in out what is the effective solvent polarization? (solve Poisson equation) Born: isolated point charge (q) in a spherical cavity of radius r immersed in a dielectric continuum with dielectric constant GBorn

q2 1 1 2r Onsager: neutral system with dipole GOnsager 1 2( 1) 3 2 2 1 r GBSA: generalized Born surface area approach

{ SASA term Gsolvation = Gcavity + Gvdw + Gpolarization ? q qi qj GBSA: generalized Born surface area approach

{ SASA term Gsolvation = Gcavity + Gvdw + Gpolarization 1 - 2 (1- e -fGB ) q q /f i j

GB fGB = (r2ij + 2ije-Dij)1/2 ij = (ij)1/2 Dij = r2ij/ 22ij : dielectric constant, qis: charges, rij: distance between atom pairs. i: Born radii; calculated numerically for each charged atom in solute; values change as calculation proceeds. : modification to incorporate salt effects at low salt via a Debye-Huckel term GBSA: generalized Born surface area approach {

SASA term Gsolvation = Gcavity + Gvdw + Gpolarization 1 - 2 (1- e -fGB ) q q /f i j GB

fGB = (r2ij + 2ije-Dij)1/2 ij = (ij)1/2 Dij = r2ij/ 22ij : dielectric constant, qis: charges, rij: distance between atom pairs. i: Born radii; calculated numerically for each charged atom in solute; values change as calculation proceeds. : modification to incorporate salt effects at low salt via a Debye-Huckel term this expression gives the Born equation for superimposed charges when rij = 0, the Onsager reaction field within 10% for a dipole in a spherical cavity when rij < 0.1 ij, and the Born plus Coulomb dielectric polarization energy within 1% for two charged spheres when rij > 2.5 ij. (Still et al., 1990). GBSA: generalized Born surface area

approach Excellent agreement to FEP for neutral small molecules with only ~2-4x overhead compared to gas phase dynamics Implemented in MacroModel, AMSOL, AMBER and CHARMM. [Recent resurgence in interest! Stable MD simulation / folding !!!] Can give stable trajectories of proteins and nucleic acids; ~2-5x cost of in-vacuo simulations. Poisson-Boltzmann electrostatics treats the solute as a low (fixed) dielectric medium containing atomic charges surrounded by a molecular surface immersed in a dielectric continuum electrostatic component from solution of PB equations; most often done numerically by finite difference hydrophobic/cavity term from surface area gives reasonable free energies of solvation can include salt effects through solution of non-linear

PB Fixed, low dielectric in interior. Charges located at atomic positions. Surface (or dielectric discontinuity) is defined by SASA or vdw surface or Dielectric continuum Poisson-Boltzmann electrostatics treats the solute as a low (fixed) dielectric medium containing atomic charges surrounded by a molecular surface immersed in a dielectric continuum electrostatic component from solution of PB equations; most often done numerically by finite difference hydrophobic/cavity term from surface area gives reasonable free energies of solvation can include salt effects through solution of non-linear PB

Issues: no microscopic detail strong dependence on radii time consuming in MD sensitivity to grid sensitivity to interior/exterior dielectric Programs: DELPHI (Honig/Sharp), MEAD (Bashford), UHBD (McCammon), GRASP (Nicholls), APBS, ... implicit METHODS: modified dielectric functions surface area methods continuum or Poisson-Boltzmann generalized Born methods [reaction field, integral equation methods] BENEFITS: Averaging is implicit; gives free energy of

solvation Conformational sampling is faster CURSE: cannot represent specific (direct) structural water interaction the computational cost is still large force field balance / METHODS: parameterization modified dielectric functions methods still under development surface areamethods (i.e. in flux) implicit

continuum or Poisson-Boltzmann generalized Born methods [reaction field, integral equation methods] BENEFITS: Averaging is implicit; gives free energy of solvation Conformational sampling is faster explicit solvation viscous damping dielectric screening hydrophobicity, , ? ISSUES: water model sampling boundary conditions (periodic or non-) long range forces, cutoff

methods/effects METHODS: [Langevin dipole] TIP3P and other models (SHAKE) explicit solvation what water model to use? what representation or boundary conditions? how to make the calculation tractable? how much water to add? Explicit water models flexible water models: few are in general use (see Levitt et al., 1995; Ferguson, 1995; Mizan et al., 1994, etc.)

rigid water models: require SHAKE/RATTLE to constrain bonds and angles TIP3P: method of choice with Cornell et al. force field [well balanced with 6-31 G* charges since prepolarized through larger fixed dipole] TIP4P/TIP5P: supported within LEaP (and PLEP) SPC/E: works as well as TIP3P (or better), not supported by default in AMBER. TIP: SPC/E: qO qH r*O O r(OH) = 0.9572 , >HOH = 104.52 r(OH) = 1.0 ,

>HOH = 109.47122 TIP3P -0.834 0.417 1.7683 0.152 TIP4P -1.04* 0.520 1.7699 0.155 SPC/E -0.8476 0.4238 1.7766 0.1554

TIP3P, TIP4P: Jorgensen et al., 1983 SPC/E: Berendsen et al., 1987 rigid water models: how do they perform? TIP3P (g/cm3) 0.982 (*) [density] 0.99 cut 0.97 PME Ei (kcal/mol) -9.86 (*)

[interaction energy] -9.7 cut -9.5 PME -11.3 cut -11.1 PME D (x10-5 cm2/s) 5.3 cut 5.1 PME 2.5 cut 2.7 PME [diffusion]

SPC/E experiment 0.997 1.0 cut 0.99 PME -9.92 2.3 at 25 reasonable density, interaction energies and 1st peak radial distribution function occupancies

3 site models underestimate compressibility TIP3P has too little structure beyond the 1st peak; less tetrahedrality than TIP4P and tends to diffuse too rapidly TIP4P has excellent density and proper oxygen structure however hydrogens are somewhat misplaced In simulations of nucleic acids, TIP3P and SPC/E lead to roughly equivalent hydration patterns explicit solvation what water model to use? what representation or boundary conditions? how to make the calculation tractable? how much water to add? explicit water boundary

conditions finite vs. infinite surrounding? (vacuum, continuum, ...) size and shape (sphere, orthorhombic, truncated octahedral) explicit water boundary conditions blob or cap truly periodic or Ewald vacuum boundary conditions blob or cap spherical shell of atoms around the site of interest

Problems: - surface tension leads to high pressure - reduced fluctuations - vacuum interface, waters can drift... vacuum boundary conditions blob or cap spherical shell of atoms around the site of interest Problems: - surface tension leads to high pressure - reduced fluctuations - vacuum interface, waters can drift... alternatives: stochastic boundary Langevin forces on near surface waters surface waters held fixed at bulk density

vacuum boundary conditions blob or cap spherical shell of atoms around the site of interest Problems: - surface tension leads to high pressure - reduced fluctuations - vacuum interface, waters can drift... alternatives: solvent boundary potentials how to develop??? + vacuum boundary conditions

blob or cap spherical shell of atoms around the site of interest Problems: - surface tension leads to high pressure - reduced fluctuations - vacuum interface, waters can drift... alternatives: reaction field or dielectric continuum What happens to waters that leave? (explicit to implicit conversion?) explicit water boundary conditions truly periodic or Ewald

explicit solvation what water model to use? what representation or boundary conditions? how to make the calculation tractable? how much water to add? nonbonded interactions Can we speed up the calculations by limiting the number of pair interactions?

CUTOFFS x x spherical minimum image reorientational motion slowed! Roberts & Schnitker, 1995 nonbonded interactions nonbonded interactions

10.0 kcal at 10 ! nonbonded interactions Can we avoid artifacts due to the cutoff of pair interactions? Avoid minimum image Dont split dipoles: charge group neutrality (AMBER) Use two cutoffs: CUT2ND>CUT x outer interactions are updated less frequently Do not ignore long ranged electrostatic interactions Ewald (Ewald, 1921), particle mesh Ewald (Essman et al., 1995; Darden & Sagui, 2000), particle-particle particle-mesh Ewald (Luty & van Gunsteren, 1996) Fast multipole method (Greengard & Rokhlin, 1989; Board et al., 1992; Lambert et

al., 1996) Cell multipole (Ding et al., 1992) improving cutoffs? shift switch shift/switch the energy or force discontinuity improving cutoffs? shift/switch the energy or force discontinuity Better to use longer cutoff than worry about details of switch/shift Best cutoff method is atom-based force shift

Switch is not-advised for electrostatics Truncation may lead to heating due to dipole reorientation near the cutoff shift switch but what about artifacts from the cutoff? Examples of the deleterious effects of the cutoff DNA falls apart (Cheatham et al., 1995) electrostatic potential (Pettitt & Smith, 1991) anomalous water transport

(Feller et al., 1996) alpha helical peptides (Schreiber & Steinhauser, 1992) attractive ion PMF (Bader & Chandler, 1992) d[CCAACGTTGG]2 overall extension of duplex, middle base pairs fray terminal base pairs fray, duplex bends

120 ps, 12 cutoff 50 ps, 16 cutoff (group based truncation) Smith & Pettitt J. Chem. Phys. 95, 8430-8441 (1991) atom switch Ewald Coulombic potential group switch

S. E. Feller, R. W. Pastor, A. Rojnuckarin, S. Bogusz, B. R. Brooks, Effect of electrostatic force truncation on interfacial and transport properties of water. J. Phys. Chem. 100, 17011-17020 (1996). density weighted water polarization profile Ewald 12 shift 12 force shift 10-12 force switch Electrostatic potential profile across vapor/ water/vapor interface Ewald (solid)

Biochemistry 31, 5856-5860 J. Phys. Chem. 96, 6423-6427 (1992) Ewald cutoff W(R): PMF of two ions Fe2+-Fe3+ CUTOFF Fe2.5+-Fe2.5+ EWALD non-integral charge to avoid dipole correction; leads to same structure and energetics not cutoff spline dependent! Artifact is > 8 kBT!!!

+ dielectric continuum + in simulation with cutoffs, two like charged ions moved closer together Periodic boundary conditions: Ewald electrostatics truly periodic or Ewald U rij

electrostatics... q q 1 i j __ UQ = + UQcorr 2 i,j=1| ri - rj | r < cutoff N ij (assume =1) U

rij electrostatics... q q 1 i j __ UQ = + UQcorr 2 i,j=1| ri - rj | r < cutoff N (assume =1) ij

boundary conditions? ...extend to truly periodic N , 1 q q __ i j UQ = 2 n i,j=1 | ri - rj + n | n: sum over lattice vectors n = nxX + nyY + nzZ , implies n = 0 & i = j omitted

, n is conditionally convergent Ewald UQ = Udirect +Ureciprocal +U self +Usurface Udirect: charges j i j periodic image

x Ewald 2 UQ = Udirect +Ureciprocal +U self +Usurface Udirect: charges j i screening counterion charge density, r 2

erfc(x) = 1 - 1/2 e-t dt 0 j periodic image (r) = 3e -2r2 3/2 r: position relative to center : width x Ewald

2 UQ = Udirect +Ureciprocal +U self +Usurface Udirect: charges j i screening counterion charge density, r i U direct= 2

erfc(x) = 1 - 1/2 e-t dt 0 j periodic image (r) = 3e -2r2 3/2 r: position relative to center : width erfc(|ri-rj|) ...extend to periodic system, qj add n and sum over all |ri-rj| charges:

x Ewald 2 UQ = Udirect +Ureciprocal +U self +Usurface charges j i Udirect: screening counterion charge density, r i

U direct= 2 erfc(x) = 1 - 1/2 e-t dt 0 j periodic image (r) = 3e -2r2 3/2 r: position relative to center : width erfc(|ri-rj|) ...extend to periodic system,

qj add n and sum over all |ri-rj| 1 __ Udirect= 2 n , charges: N q q i j

i,j=1 erfc(|ri-rj+n|) |ri-rj+n| Simplify? adjust so only terms within the cutoff are significant... Ewald: reciprocal sum Ureciprocal: ...solution to Poissons equation with periodic boundary conditions and source term(r) [expand in Fourier series, solve term by term...] Ewald: reciprocal sum Ureciprocal:

...solution to Poissons equation with periodic boundary conditions and source term(r) [expand in Fourier series, solve term by term...] 1 __ Ur = V e m=0 2 2 2 - m / m 2

N S(m) = qi e | S(m) | 2im.|ri| 2 V: volume of unit cell m: mxX*+myY*+MzZ* ? = c.F(Q) i=1

The key to the speed is the FFT on Q, the cardinal spline interpolated charge grid... particle mesh Ewald: other terms... self term: Uself = - 1/2 qi2 surface term: assume spherical boundary (dipole, surroundings, shape) tin foil boundary: = implies Usurface = 0

Usurface= - 2 3v PME extensions... r-6 Ewald (particle mesh) PME into gibbs (free energy simulation) further parallelization and single PE optimization, beyond SPMD? Full conversion to MPI and distribution of source with AMBER 4.x

PME by default in AMBER6 PME with dipoles, quadrupoles (i.e. polarization) qi.ri2 Consequences of imposing true periodicity dipole reorientation vs. induced correlations + - vs. + - + - + -

Consequences of imposing true periodicity dipole reorientation vs. induced correlations + - vs. + - + - + surface effects? tin-foil vs. vacuum boundary conditions non-ideal (effectively increased) molality anisotropic interaction energy uniform neutralization background, G force field artifacts? Consequences of imposing true

periodicity Smith & Pettitt: true periodicity is not a problem in solvents with a reasonably high dielectric! dipole rotation protein rotation (effect < kT) Cheatham & Kollman, Norberto de Souza & Ornstein no apparent dependence on box size on DNA 5 , 10 or 15 of water surrounding the DNA little salt effects .but: Consequences of imposing true periodicity 20 box

30 box 40 box artificial stabilizati on of helix poly-alanine octapeptide, 2 ns simulations Hnenberger & McCammon J. Phys. Chem. B. 104, 3668-3675 (2000) ontinuum calculations show artifacts, Confirmed in MD simulations

Cutoff or Ewald? If using a cutoff, use a good cutoff method such as atom-based force shifting. Be aware of potential artifacts (water transport properties, ) If using Ewald (i.e. the correct physics) make sure there is enough surrounding solvent or sufficient dielectric to screen image artifacts. Note: for many systems (such as systems that are not highly charged) these artifacts may not be large!