Presented By: Yuchun Lin

1. International Workshop on Continuum Modeling of Biomolecules September 14-16, 2009 in Beijing, China An image-based reaction field method for electrostatic interactions in molecular dynamics simulations Presented By: Yuchun Lin Department of Mathematics & Statistics Department of Physics & Optical Science University of North Carolina at Charlotte

2. Introduction & Background Molecular Dynamics Simulation Simulation of biological macromolecules is a key area of interest: • Understand the dynamic mechanisms of macromolecular function (protein folding, enzymaticcatalysis) • Predict the energetics of various biological processes (ligand association, protein stability) • Design novel molecules with particular properties (drug design, protein engineering) • It still has some issues. • Accurate simulations require the solvent to be treated carefully. • Long range interaction: Truncation of electrostatic interaction leads to artifacts.

3. Introduction & Background e Reaction Field Hybrid implicit/explicit e Explicit Implicit More Accurate & Less Efficient More Efficient & Less Accurate

4. Hybrid Solvation Models • Kirkwood expansion • --- slow convergence at boundary Friedman image expression --- approximated & less accurate accurate up to O(1/ε) • Repulsive potential applied • --- strong surface effect • Exact solution of PB in particular geometries • Arbitrary geometry Numerical Solution: D. Beglov, et al, J. Chem. Phys. 100(1994) 9050 H. Alper, et al, J. Chem. Phys.,99(1993) 9847 G. Brancato, et al, J. Chem. Phys. 122(2005) 154109 Image Approximation: P. K. Yang, et al, J. Phys. Chem. B, 106 (2002) 2973. G. Petraglio, et al, J. Chem. Phys. 123(2005) 044103 A. Wallqvist, Mol. Sim. 10(1993) 13–17. Numerical Solution: W. Im, et al., J. Chem. Phys. 114(2001) 2924 Generalized Born Model: M. S. Lee, et al, J. Comput. Chem. 25 (2004) 1967 D. Bashford, et al., Annu. Rev. Phys. Chem. 51 (2000) 129–152

5. Basic Idea Image-based method to compute reaction field Known Drawbacks Friedman expression for reaction field is approximate Surface effects are non-negligible or not removed easily Our Solutions Periodic boundary conditions for non-electrostatic Multiple image charges method Y. Lin, A. Baumketner, S. Deng, Z. Xu, D. Jacob, W. Cai, An image-based reaction field method for electrostatic interactions in molecular dynamics simulations of aqueous solutions, J. Chem. Phys., 2009, under review

6. Theory: RF in multiple-image charges approach Poisson-Boltzmann equations: H. L. Friedman, Mol. Phys. 29 (1975) 1533–1543 W. Cai, S. Deng, D. Jacobs, J. Comput. Phys., 223(2007), 846-864 S. Deng, W. Cai, Comm. Comput. Phys., 2(2007), 1007-1026

7. Theory: RF in multiple-image charges approach With Kirkwood expansion on pure solution case For using image method, let , , Where and First series is the potential of Kelvin image: Using the integral identity and rewrite second series as:

8. Theory: RF in multiple-image charges approach Now the reaction field inside the cavity is: Next, we construct discrete image charge by Gauss-Radauquadrature: Here are the Gauss-Radauquadrature weights and points. Since s1=-1 and then x1=rK, the classical Kelvin image charge and the first discrete image charge can be combined, leading to:

9. Theory: Integration of the RF model with MD Role of a buffer layer between explicit and implicit solvents : A. Wallqvist, Mol. Sim. 10(1993) 13–17 L.Wang, J. Hermans, J. Phys. Chem. 99(1995) 12001

10. Theory: Integration of the RF model with MD Choice of boundary conditions: d L Choice of box type: For Cubic Box: For L = 45Å, τ = 5Å box, Cube allows only 2 Å for d

11. Model • Three parameters: • Number of image charge (Ni) • Ni=0, 1, 2, 3 • Thickness of buffer layer (τ) • τ=2Å, 4Å, 6Å, 8Å • Box size (L) • L=30Å, 45Å, 60Å d Fast Multipole Method is applied L. Greengard, V. Rokhlin, J. Comput. Phys. 73 (1987) 325–348 For L=45Å, τ=5Å box, TO allows 17 Å for d For Truncated Octahedron:

12. Results: Relative Density L=30Å Standard Deviation L=45Å • A buffer layer of at least 6 Å is required to yield uniform density. • Large surface effect at low τ. L=60Å

13. Results: Radial Distribution Function Number of image charges (Ni) is not critical Effect of buffer layer thickness (τ) is unnoticeable Effect of box size (L), converges on L=60Å with PME

14. Results: Diffusion Coefficient Reaction field is critical for the proper description of diffusion (Unit: 10-9m2s-1)

15. Results: Dielectric Constant The convergence with the number of image charges occurs at Ni = 1 V. Ballenegger, J. P. Hansen, J. Chem. Phys. 122(2005) 114711 L=60Å, τ=4Å

16. Results: Dielectric Constant The dependence of ε on the thickness of buffer layer is week Å Å Å Å Å

17. Results: Dielectric Constant Dielectric properties require large simulation boxes and RF corrections PME: ε = 90± 10

18. Summary & Conclusions • Summary: • Large enough buffer layer is important • Large box size produces good bulk • properties of simulated water • Reaction field is essential for proper • description of dielectric permittivity • Conclusion: • A new solvation model is proposed. Static, structural and dynamic properties of water are well reproduced compared to PME. • Applications to biological system are our future work. Optimal parameters L = 60Å, τ= 6Å, Ni = 1. W. Cai, S. Deng, D. Jacobs, J. Comput. Phys., 223(2007), 846-864 S. Deng, W. Cai, Comm. Comput. Phys., 2(2007), 1007-1026 S. Deng, W. Cai, J. Comput. Phys. 227 (2007) 1246–1266. Y. Lin, A. Baumketner, S. Deng, Z. Xu, D. Jacobs, W. Cai, J. Chem. Phys., 2009, under review

19. Acknowledgement Advisors: Dr. AndrijBaumketner Dr. Wei Cai Dr. Shaozhong Deng Dr. Don Jacobs Group Members: Dr. Xia Ji Dr. Haiyan Jiang Dr. Boris Ni Dr. ZhenliXu Ms. Katherine Baker Mr. Wei Song Ms. Ming Xiang Funding by：