Convergence of the Many-Body Expansion for Energy and Forces for Classical Polarizable Models in the Condensed Phase

We analyze convergence of energies and forces for the AMOEBA classical polarizable model when evaluated as a many-body expansion (MBE) against the corresponding <i>N</i>-body parent potential in the context of a condensed-phase water simulation. This is in contrast to most MBE formulations based on quantum mechanics, which focus only on convergence of energies for gas-phase clusters. Using a single water molecule as a definition of a body, we find that truncation of the MBE at third order, 3-AMOEBA, captures direct polarization exactly and yields apparent good convergence of the mutual polarization energy. However, it renders large errors in the magnitude of polarization forces and requires at least fourth-order terms in the MBE to converge toward the parent potential gradient values. We can improve the convergence of polarization forces for 3-AMOEBA by embedding the polarization response of dimers and trimers within a complete representation of the fixed electrostatics of the entire system. We show that the electrostatic embedding formalism helps identify the specific configurations involving linear hydrogen-bonding arrangements that are poorly convergent at the 3-body level. By extending the definition of a body to be a large water cluster, we can reduce errors in forces to yield an approximate polarization model that is up to 10 times faster than the parent potential. The 3-AMOEBA model offers new ways to investigate how the properties of bulk water depend on the degree of connectivity in the liquid

A New Method for Treating Drude Polarization in Classical Molecular Simulation

With polarization becoming an increasingly common feature in classical molecular simulation, it is important to develop methods that can efficiently and accurately evaluate the many-body polarization solution. In this work, we expand the theoretical framework of our inertial extended Langrangian, self-consistent field iteration-free method (iEL/0-SCF), introduced for point induced dipoles, to the polarization model of a Drude oscillator. When applied to the polarizable simple point charge model (PSPC) for water, our iEL/0-SCF method for Drude polarization is as stable as a well-converged SCF solution and more stable than traditional extended Lagrangian (EL) approaches or EL formulations based on two temperature ensembles where Drude particles are kept “colder” than the real degrees of freedom. We show that the iEL/0-SCF method eliminates the need for mass repartitioning from parent atoms onto Drude particles, obeys system conservation of linear and angular momentum, and permits the extension of the integration time step of a basic molecular dynamics simulation to 6.0 fs for PSPC water

Evolution of the Potential Energy Landscape with Static Pulling Force for Two Model Proteins

The energy landscape is analyzed for off-lattice bead models of protein L and protein G as a function of a static pulling force. Two different pairs of attachment points (pulling directions) are compared in each case, namely, residues 1/56 and 10/32. For the terminal residue pulling direction 1/56, the distinct global minimum structures are all extended, aside from the compact geometry that correlates with zero force. The helical turns finally disappear at the highest pulling forces considered. For the 10/32 pulling direction, the changes are more complicated, with a variety of competing arrangements for beads outside the region where the force is directly applied. These alternatives produce frustrated energy landscapes, with low-lying minima separated by high barriers. The calculated folding pathways in the absence of force are in good agreement with previous work. The N-terminal hairpin folds first for protein L and the C-terminal hairpin for protein G, which exhibits an intermediate. However, for a relatively low static force, where the global minimum retains its structure, the folding mechanisms change, sometimes dramatically, depending on the protein and the attachment points. The scaling relations predicted by catastrophe theory are found to hold in the limit of short path lengths

Calculating the Bimolecular Rate of Protein–Protein Association with Interacting Crowders

We have recently introduced a method termed Poisson–Boltzmann semianalytical method (PB-SAM) for solving the linearized Poisson–Boltzmann equation for large numbers of arbitrarily shaped dielectric cavities with controlled precision. In this work we extend the applicability of the PB-SAM approach by deriving force and torque expressions that fully account for mutual polarization in both the zero- and first-order derivatives of the surface charges, that can now be embedded into a Brownian dynamics scheme to look at electrostatic-driven mesoscale assembly and kinetics. We demonstrate the capabilities of the PB-SAM approach by simulating the protein concentration effects on the bimolecular rate of association of barnase and barstar, under periodic boundary conditions and evaluated through mean first passage times. We apply PB-SAM to the pseudo-first-order reaction rate conditions in which either barnase or barstar are in great excess relative to the other protein (124:1). This can be considered a specific case in which the PB-SAM approach can be applied to crowding conditions in which crowders are not inert but can form interactions with other molecules

Assessing Ion–Water Interactions in the AMOEBA Force Field Using Energy Decomposition Analysis of Electronic Structure Calculations

AMOEBA is a molecular mechanics force field that addresses some of the shortcomings of a fixed partial charge model, by including permanent atomic point multipoles through quadrupoles, as well as many-body polarization through the use of point inducible dipoles. In this work, we investigate how well AMOEBA formulates its non-bonded interactions, and how it implicitly incorporates quantum mechanical effects such as charge penetration (CP) and charge transfer (CT), for water–water and water–ion interactions. We find that AMOEBA’s total interaction energies, as a function of distance and over angular scans for the water dimer and for a range of water-monovalent cations, agree well with an advanced density functional theory (DFT) model, whereas the water-halides and water-divalent cations show significant disagreement with the DFT result, especially in the compressed region when the two fragments overlap. We use a second-generation energy decomposition analysis (EDA) scheme based on absolutely localized molecular orbitals (ALMOs) to show that in the best cases AMOEBA relies on cancellation of errors by softening of the van der Waals (vdW) wall to balance permanent electrostatics that are too unfavorable, thereby compensating for the missing CP effect. CT, as another important stabilizing effect not explicitly taken into account in AMOEBA, is also found to be incorporated by the softened vdW interaction. For the water-halides and water-divalent cations, this compensatory approach is not as well executed by AMOEBA over all distances and angles, wherein permanent electrostatics remains too unfavorable and polarization is overdamped in the former while overestimated in the latter. We conclude that the DFT-based EDA approach can help refine a next-generation AMOEBA model that either realizes a better cancellation of errors for problematic cases like those illustrated here, or serves to guide the parametrization of explicit functional forms for short-range contributions from CP and/or CT

Mechanism of Nucleation and Growth of Aβ40 Fibrils from All-Atom and Coarse-Grained Simulations

In this work, we characterize the nucleation and elongation mechanisms of the “diseased” polymorph of the amyloid-β 40 (Aβ40) fibril using an off-lattice coarse-grained (CG) protein model. After determining the nucleation size and subsequent stable protofibrillar structure from the CG model, validated with all-atom simulations, we consider the “lock and dock” and “activated monomer” fibril elongation mechanisms for the protofibril by statistical additions of a monomer drawn from four different ensembles of the free Aβ40 peptide to grow the fibril. Our CG model shows that the dominant mechanism for fibril elongation is the lock and dock mechanism across all monomer ensembles, even when the monomer is in the activated form. Although our CG model finds no thermodynamic difference between the two fibril elongation mechanisms, the activated monomer is found to be kinetically faster by a factor of 2 for the “locking” step compared with all other structured or unstructured monomer ensembles

Assessing Ion–Water Interactions in the AMOEBA Force Field Using Energy Decomposition Analysis of Electronic Structure Calculations

AMOEBA is a molecular mechanics force field that addresses some of the shortcomings of a fixed partial charge model, by including permanent atomic point multipoles through quadrupoles, as well as many-body polarization through the use of point inducible dipoles. In this work, we investigate how well AMOEBA formulates its non-bonded interactions, and how it implicitly incorporates quantum mechanical effects such as charge penetration (CP) and charge transfer (CT), for water–water and water–ion interactions. We find that AMOEBA’s total interaction energies, as a function of distance and over angular scans for the water dimer and for a range of water-monovalent cations, agree well with an advanced density functional theory (DFT) model, whereas the water-halides and water-divalent cations show significant disagreement with the DFT result, especially in the compressed region when the two fragments overlap. We use a second-generation energy decomposition analysis (EDA) scheme based on absolutely localized molecular orbitals (ALMOs) to show that in the best cases AMOEBA relies on cancellation of errors by softening of the van der Waals (vdW) wall to balance permanent electrostatics that are too unfavorable, thereby compensating for the missing CP effect. CT, as another important stabilizing effect not explicitly taken into account in AMOEBA, is also found to be incorporated by the softened vdW interaction. For the water-halides and water-divalent cations, this compensatory approach is not as well executed by AMOEBA over all distances and angles, wherein permanent electrostatics remains too unfavorable and polarization is overdamped in the former while overestimated in the latter. We conclude that the DFT-based EDA approach can help refine a next-generation AMOEBA model that either realizes a better cancellation of errors for problematic cases like those illustrated here, or serves to guide the parametrization of explicit functional forms for short-range contributions from CP and/or CT

The Importance of the Scaffold for <i>de Novo</i> Enzymes: A Case Study with Kemp Eliminase

We report electric field values relevant to the reactant and transition states of designed Kemp eliminases KE07 and KE70 and their improved variants from laboratory directed evolution (LDE), using atomistic simulations with the AMOEBA polarizable force field. We find that the catalytic base residue contributes the most to the electric field stabilization of the transition state of the LDE variants of the KE07 and KE70 enzymes, whereas the electric fields of the remainder of the enzyme and solvent <i>disfavor</i> the catalytic reaction in both cases. By contrast, we show that the electrostatic environment plays a large and stabilizing role for the naturally occurring enzyme ketosteroid isomerase (KSI). These results suggest that LDE is ultimately a limited strategy for improving <i>de novo</i> enzymes since it is largely restricted to optimization of chemical positioning in the active site, thus yielding a ∼3 order magnitude improvement over the uncatalyzed reaction, which we suggest may be an absolute upper bound estimate based on LDE applied to comparable <i>de novo</i> Kemp eliminases and other enzymes like KSI. Instead <i>de novo</i> enzymatic reactions could more productively benefit from optimization of the electrostatics of the protein scaffold in early stages of the computational design, utilizing electric field optimization as guidance

Higher-Order Extended Lagrangian Born–Oppenheimer Molecular Dynamics for Classical Polarizable Models

Generalized extended Lagrangian Born–Oppenheimer molecular dynamics (XLBOMD) methods provide a framework for fast iteration-free simulations of models that normally require expensive electronic ground state optimizations prior to the force evaluations at every time step. XLBOMD uses dynamically driven auxiliary degrees of freedom that fluctuate about a variationally optimized ground state of an approximate “shadow” potential which approximates the true reference potential. While the requirements for such shadow potentials are well understood, constructing such potentials in practice has previously been <i>ad hoc</i>, and in this work, we present a systematic development of XLBOMD shadow potentials that match the reference potential to any order. We also introduce a framework for combining friction-like dissipation for the auxiliary degrees of freedom with general-order integration, a combination that was not previously possible. These developments are demonstrated with a simple fluctuating charge model and point induced dipole polarization models

Coexistence of Multilayered Phases of Confined Water: The Importance of Flexible Confining Surfaces

Flexible nanoscale confinement is critical to understanding the role that bending fluctuations play on biological processes where soft interfaces are ubiquitous or to exploit confinement effects in engineered systems where inherently flexible 2D materials are pervasively employed. Here, using molecular dynamics simulations, we compare the phase behavior of water confined between flexible and rigid graphene sheets as a function of the in-plane density, ρ_2D. We find that both cases show commensurate mono-, bi-, and trilayered states; however, the water phase in those states and the transitions between them are qualitatively different for the rigid and flexible cases. The rigid systems exhibit discontinuous transitions between an (<i>n</i>)-layer and an (<i>n</i>+1)-layer state at particular values of ρ_2D, whereas under flexible confinement, the graphene sheets bend to accommodate an (<i>n</i>)-layer and an (<i>n</i>+1)-layer state coexisting in equilibrium at the same density. We show that the flexible walls introduce a very different sequence of ice phases and their phase coexistence with vapor and liquid phases than that observed with rigid walls. We discuss the applicability of these results to real experimental systems to shed light on the role of flexible confinement and its interplay with commensurability effects