Determination of Electrostatic Parameters for a Polarizable Force Field Based on the Classical Drude Oscillator

Abstract

A procedure to determine the electrostatic parameters has been developed for a polarizable empirical force field based on the classical Drude oscillator model. Atomic charges and polarizabilities for a given molecule of interest were derived from restrained fitting to quantum-mechanical electrostatic potentials (ESP) calculated at the B3LYP/ cc-pVDZ or B3LYP/aug-cc-pVDZ levels on grid points located on concentric Connolly surfaces. The determination of the atomic polarizabilities requires a series of perturbed ESP maps, each one representing the electronic response of the molecule in the presence of a background charge placed on Connolly surfaces primarily along chemical bonds and lone pairs. Reference values for the partial atomic charges were taken from the CHARMM27 additive all-atom force field, and those for the polarizabilities were based on adjusted Miller's ahp atomic polarizability values. The fitted values of atomic polarizabilities were scaled to reflect the reduced polarization expected for the condensed media and/or to correct for the systematic underestimation of experimental molecular polarizabilities by B3LYP calculations. Following correction of the polarizabilities, the atomic charges were adjusted to reproduce gas-phase dipole moments. The developed scheme has been tested on a set of small molecules representing functional moieties of nucleic acids. The derived electrostatic parameters have been successfully applied in a preliminary polarizable molecular dynamics simulation of a DNA octamer in a box of water with sodium counterions. Thus, this study confirms the feasibility of the use of a polarizable force field based on a classical Drude model for simulations of biomolecules in the condensed phase.

Empirical force fields for biological macromolecules: Overview and issues

“…The first reported study of macromolecular simulations involved proteins in the absence of solvent, with the duration of the simulations being only 2 ps. 84 More recently, simulations of several small proteins using a fully polarizable force field, including solvent and for time durations into the ns range 85 and a simulation of DNA in solution, including counterions, 86 have been reported. Thus, progress is being made towards the development of polarizable force fields for biomolecules.…”

Kirkwood-Buff analysis of aqueous N-methylacetamide and acetamide solutions modeled by the CHARMM additive and Drude polarizable force fields

“…Electrostatic parameters, including the partial atomic charges, atomic polarizabilities, and the atom-atom Thole scale factors, were developed for each model compound following the methodology described before. 23 Initial values of the partial atomic charges were either taken from the C22 additive all-atom force field 46 or, in the case of ACEM, the initial charges were from NMA. Starting values of the polarizabilities were adjusted Miller's atomic hybrid polarizabilities (ahp) values, 23,89,90 and the Thole factors were set to 1.3, the value originally determined for benzene.…”

“…23 Initial values of the partial atomic charges were either taken from the C22 additive all-atom force field 46 or, in the case of ACEM, the initial charges were from NMA. Starting values of the polarizabilities were adjusted Miller's atomic hybrid polarizabilities (ahp) values, 23,89,90 and the Thole factors were set to 1.3, the value originally determined for benzene. The fitted electrostatic models yielded dipole moments in good agreement with the reference values as shown in Table III.…”

“…21 Polarizable force field development in our laboratory has primarily been based on the classical Drude oscillator model. [22][23][24][25] In the model a point charge, called the Drude oscillator or particle, is attached to each non-hydrogen atom via a harmonic spring. Induced polarization is modeled by Drude particles relaxing in response to their surrounding electrostatic environment and seeking their local energy minima for any given configuration.…”

“…32 As described below the optimization was performed following standard protocols used for the Drude model in our laboratory. 23,24,[29][30][31] The final parameter set was then used in the KB analysis.…”