Discontinuous Molecular Dynamics Simulations Research Articles

Molecular simulation of protein dynamics in nanopores. I. Stability and folding

Discontinuous molecular dynamics simulations, together with the protein intermediate resolution model, an intermediate-resolution model of proteins, are used to carry out several microsecond-long simulations and study folding transition and stability of alpha-de novo-designed proteins in slit nanopores. Both attractive and repulsive interaction potentials between the proteins and the pore walls are considered. Near the folding temperature T(f) and in the presence of the attractive potential, the proteins undergo a repeating sequence of folding/partially folding/unfolding transitions, with T(f) decreasing with decreasing pore sizes. The unfolded states may even be completely adsorbed on the pore's walls with a negative potential energy. In such pores the energetic effects dominate the entropic effects. As a result, the unfolded state is stabilized, with a folding temperature T(f) which is lower than its value in the bulk and that, compared with the bulk, the folding rate decreases. The opposite is true in the presence of a repulsive interaction potential between the proteins and the walls. Moreover, for short proteins in very tight pores with attractive walls, there exists an unfolded state with only one alpha-helical hydrogen bond and an energy nearly equal to that of the folded state. The proteins have, however, high entropies, implying that they cannot fold onto their native structure, whereas in the presence of repulsive walls the proteins do attain their native structure. There is a pronounced asymmetry between the two termini of the protein with respect to their interaction with the pore walls. The effect of a variety of factors, including the pore size and the proteins' length, as well as the temperature, is studied in detail.

Transferable step potentials for perfluorinated hydrocarbons

The step potential equilibria and discontinuous molecular dynamics (SPEADMD) model is adapted for characterizing the interaction potentials of perfluorocarbons and their mixtures with n-alkanes. We seek to explain the peculiar behavior of these systems, especially with regard to the unfavorable mixing behavior. The methodology is based on discontinuous molecular dynamics (DMD) and second order thermodynamic perturbation theory (TPT). DMD simulation is applied to the repulsive part of the potential. The effects of disperse attractions and hydrogen bonding are treated by TPT, discretizing the attractive potential into four distinct wells of variable depth. This approach accelerates the molecular simulations in general and the parameterization of the transferable potentials in particular. C 3–C 8 straight chain perfluorocarbons are characterized along with perfluorobenzene, perfluorocyclobutane, and heptafluoropropane applying explicit atom models for all fluorine atoms. Each compound is simulated at 21 densities. Interpolation with density combines with TPT to give a complete equation of state. The depths of the attractive wells are optimized by iterating on their values until the vapor pressures computed by the resulting equation of state provide the minimum deviation from experimental data. Step potentials offer the prospect of setting the step depths without preconceived notions regarding the shape of the potential. In the case of perfluorocarbons, it is found that the minimum of the potential is reached at roughly 1.8 σ, with the range from σ to 1.8 σ acting as a soft “shoulder” region which is weakly attractive, where σ designates the site diameter. This shoulder region can explain the unusually low heat of vaporization of perfluorocarbons and high vapor pressures relative to compounds with similar molecular volume. Application of this approach results in average absolute deviations of 8% in vapor pressure and 3% in liquid density. Coincidentally, the soft shoulder also explains the weak mixture interactions between perfluorocarbons and n-alkanes. The Lorentz–Berthelot combining rule assumed to describe site–site interactions between corresponding steps leads to ɛ ijk ∼ ( ɛ iik ɛ jjk ) 1/2, where ɛ ijk is the depth of the potential well between sites i and j in the kth step. For steps at short range, the alkane steps are at their deepest, but the perfluorocarbon interactions are weak, making the mixed interaction weak. For steps at long range, the perfluorocarbon steps are deeper, but the alkane interactions are weak, again making the mixed interaction weak.

Temperature-Dependent Folding Pathways of Pin1 WW Domain: An All-Atom Molecular Dynamics Simulation of a Gō Model

We study the folding thermodynamics and kinetics of the Pin1 WW domain, a three-stranded β-sheet protein, by using all-atom (except nonpolar hydrogens) discontinuous molecular dynamics simulations at various temperatures with a Gō model. The protein exhibits a two-state folding kinetics near the folding transition temperature. A good agreement between our simulations and the experimental measurements by the Gruebele group has been found, and the simulation sheds new insights into the structure of transition state, which is hard to be straightforwardly captured in experiments. The simulation also reveals that the folding pathways at approximately the transition temperature and at low temperatures are much different, and an intermediate state at a low temperature is predicted. The transition state of this small β-protein at its folding transition temperature has a well-established hairpin 1 made of β1 and β2 strands while its low-temperature kinetic intermediate has a formed hairpin 2 composed of β2 and β3 strands. Theoretical results are compared with other simulation results as well as available experimental data. This study confirms that specific side-chain packing in an all-atom Gō model can yield a reasonable prediction of specific folding kinetics for a given protein. Different folding behaviors at different temperatures are interpreted in terms of the interplay of entropy and enthalpy in folding process.

Does confining the hard-sphere fluid between hard walls change its average properties?

We use grand canonical transition-matrix Monte Carlo and discontinuous molecular dynamics simulations to generate precise thermodynamic and kinetic data for the equilibrium hard-sphere fluid confined between smooth hard walls. These simulations show that the pronounced inhomogeneous structuring of the fluid normal to the confining walls, often the primary focus of density functional theory studies, has a negligible effect on many of its average properties over a surprisingly broad range of conditions. We present one consequence of this insensitivity to confinement: a simple analytical equation relating the average density of the confined fluid to that of the bulk fluid with equal activity. Nontrivial implications of confinement for average fluid properties do emerge in this system, but only when the fluid is both (i) dense and (ii) confined to a gap smaller than approximately three particle diameters. For this limited set of conditions, we find that "in-phase" oscillatory deviations in excess entropy and self-diffusivity (relative to the behavior of the bulk fluid at the same average density) occur as a function of gap size. These paired thermodynamic/kinetic deviations from bulk behavior appear to reflect the geometric packing frustration that arises when the confined space cannot naturally accommodate an integer number of particle layers.

Structure and dynamics of short chain molecules in disordered porous materials: A molecular dynamics simulation study

The static and dynamic properties of short polymer chains in disordered materials are studied using discontinuous molecular dynamics simulations. The polymers are modeled as chains of hard spheres and the matrix is a collection of fixed hard spheres. The simulations show that the chain size is a nonmonotonic function of the matrix concentration for all polymer concentrations. The dependence of polymer diffusion D on the degree of polymerization N becomes stronger as the matrix concentration is increased. At high matrix concentrations we observe a decoupling between translational and rotational diffusion, i.e., the rotational relaxation time becomes very large but the translational diffusion is not affected significantly. We attribute this to the trapping of a small number of polymers. Under these conditions the polymer chains diffuse via a hopping mechanism.

Numerical implementation of the exact dynamics of free rigid bodies

In this paper the exact analytical solution of the motion of a rigid body with arbitrary mass distribution is derived in the absence of forces or torques. The resulting expressions are cast into a form where the dependence of the motion on initial conditions is explicit and the equations governing the orientation of the body involve only real numbers. Based on these results, an efficient method to calculate the location and orientation of the rigid body at arbitrary times is presented. This implementation can be used to verify the accuracy of numerical integration schemes for rigid bodies, to serve as a building block for event-driven discontinuous molecular dynamics simulations of general rigid bodies, and for constructing symplectic integrators for rigid body dynamics.

Discontinuous molecular dynamics simulation study of polymer collapse

Discontinuous molecular dynamics simulations were used to study the coil-globule transition of a polymer in an explicit solvent. Two different versions of the model were employed, which are differentiated by the nature of monomer-solvent, solvent-solvent, and nonbonded monomer-monomer interactions. For each case, a model parameter lambda determines the degree of hydrophobicity of the monomers by controlling the degree of energy mismatch between the monomers and solvent particles. We consider a lambda-driven coil-globule transition at constant temperature. The simulations are used to calculate average static structure factors, which are then used to determine the scaling exponents of the system in order to determine the theta-point values lambda(theta) separating the coil from the globule state. For each model we construct coil-globule phase diagrams in terms of lambda and the particle density rho. Additionally, we explore for each model the effects of varying the range of the attractive interactions on the phase boundary separating the coil and globule phases. The results are analyzed in terms of a simple Flory-type theory of the collapse transition.

Dynamics of Chain Molecules in Disordered Materials

The dynamic behavior of hard chains in disordered materials composed of fixed hard spheres is studied using discontinuous molecular dynamics simulations. The matrix induces entanglements in the chain fluid, i.e., for high matrix densities the diffusion coefficient D scales with the chain length N as D is approximately N(-2). At high matrix densities the rotational relaxation time becomes very large but the translational diffusion is not affected significantly; i.e., the chains display a dynamic heterogeneity reminiscent of probe diffusion in supercooled liquids and glasses. We show that this is because some chains are trapped, and move via a hopping mechanism. There are no signatures of this dynamic heterogeneity in the matrix static structure, however, which is identical to that of a hard-sphere liquid.

Spontaneous fibril formation by polyalanines; discontinuous molecular dynamics simulations.

Fibrillary protein aggregates rich in beta-sheet structure have been implicated in the pathology of several neurodegenerative diseases. In this work, we investigate the formation of fibrils by performing discontinuous molecular dynamics simulations on systems containing 12 to 96 model Ac-KA(14)K-NH(2) peptides using our newly developed off-lattice, implicit-solvent, intermediate-resolution model, PRIME. We find that, at a low concentration, random-coil peptides assemble into alpha-helices at low temperatures. At intermediate concentrations, random-coil peptides assemble into alpha-helices at low temperatures and large beta-sheet structures at high temperatures. At high concentrations, the system forms beta-sheets over a wide range of temperatures. These assemble into fibrils above a critical temperature which decreases with concentration and exceeds the isolated peptide's folding temperature. At very high temperatures and all concentrations, the system is in a random-coil state. All of these results are in good qualitative agreement with those by Blondelle and co-workers on Ac-KA(14)K-NH(2) peptides. The fibrils observed in our simulations mimic the structural characteristics observed in experiments in terms of the number of sheets formed, the values of the intra- and intersheet separations, and the parallel peptide arrangement within each beta-sheet. Finally, we find that when the strength of the hydrophobic interaction between nonpolar side chains is high compared to the strength of hydrogen bonding, amorphous aggregates, rather than fibrillar aggregates, are formed.

Transferable step potentials for amines, amides, acetates, and ketones

The vapor pressure and liquid density of 11 amines, 5 primary amides, 9 acetates, and 15 ketones are studied to characterize transferable step potential models. Discontinuous molecular dynamics simulations are performed for reference potentials with details of the molecular structure and bond lengths. The simulation results are interpreted with thermodynamic perturbation theory to provide an efficient basis for molecular modeling and characterization of the induction forces. The attractive potential is represented as four steps with only the first and last steps varied independently. The two middle steps are interpolated such that each site type is characterized by three parameters: the diameter, the depth of the inner well, and the depth of the outer well. The depths of the attractive wells are optimized to fit experimental data for vapor pressure and liquid density. In general, the vapor pressure is correlated to an overall 12% average absolute deviation (%AAD) and the liquid density to 4% AAD. These deviations are larger than the errors encountered in characterizing straight chain hydrocarbons but consistent with recent results for naphthenic and aromatic compounds. The subtle variability in molecular structure caused by the consideration of isomers with heteroatoms and branching in varied positions puts a strain on the assumption of transferability. An additional challenge is encountered as the coverage of site types is extended to broader classes of families; the availability of reliable experimental data becomes much more limited. This lack of experimental data is especially a problem for the primary amides in this study and a likely problem for most future extensions.

Parametric Studies of Interaction Strengths in Polymer/CO2 Systems: Discontinuous Molecular Dynamics Simulations

Discontinuous molecular dynamics simulations are performed on homopolymer/solvent and surfactant/solvent systems. The homopolymer and surfactant molecules are modeled as freely jointed square-well chains. Solvent molecules are modeled as both hard spheres and square-well spheres. We explore how the various interaction parameters affect the types of phase behavior and micellization observed in the homopolymer/solvent and surfactant/solvent systems. Increasing the packing fraction of homopolymers in both hard-sphere solvents and square-well solvents increases the solvent's ability to dissolve homopolymers only when the segment-solvent interaction strength exceeds a critical value. Although only upper critical solution temperature (UCST) behavior is observed for homopolymers in hard-sphere solvents, both UCST and lower critical solution temperature (LCST) behavior are observed for homopolymers in square-well solvents, depending upon the interaction strengths and chain length. This indicates that it is necessary to account for the solvent-solvent attraction to model LCST behavior in supercritical CO2. Our simulation results on surfactants in hard-sphere solvents show that it is necessary to account for the interactions experienced by both the head and tail blocks in order to capture the essential features of surfactant/supercritical CO2 systems.

The dynamics of single chains within a model polymer melt

Discontinuous molecular dynamics simulations are performed on a system containing 32 hard chains of length 192 at a volume fraction of phi = 0.45 to explore the idea that localized entanglements have a significant effect on the dynamics of the individual chains within an entangled polymer melt. Anomalous behavior can still be observed when studying the dynamics of the individual chains, although increased time averaging causes the anomalous relaxation-memory-release behavior that was observed previously in the system to smooth out. First, the individual chain mean squared displacements and apparent diffusion coefficients are calculated, and a wide distribution of diffusive behavior is found. Although the apparent diffusion coefficient curve averaged over all chains displays the predicted long-time diffusive behavior, the curves for the individual chains differ both qualitatively and quantitatively. They display superdiffusive, diffusive, and subdiffusive behavior, with the largest percentage of chains exhibiting superdiffusive behavior and the smallest percentage exhibiting the predicted diffusive behavior. Next, the individual chain end-to-end vector autocorrelation functions and relaxation times are determined, and a wide distribution of stress relaxation behavior is found. The times when the end-to-end vector autocorrelation functions relax completely span almost an order of magnitude in reduced time. For some chains, the end-to-end vector autocorrelation function relaxes smoothly toward zero similar to the system average; however, for other chains the relaxation is slowed greatly, indicating the presence of additional entanglements. Almost half of the chains exhibit the anomalous behavior in the end-to-end vector autocorrelation function. Finally, the dynamic properties are displayed for a single chain exhibiting anomalous relaxation-memory-release behavior, supporting the idea that the relaxation-memory-release behavior is a single-chain property.

Computer Simulation of Block Copolymer/Nanoparticle Composites

Discontinuous molecular dynamics simulation is used to study the phase behavior of diblock copolymer/nanoparticle composites. The copolymers are modeled as chains of tangent hard-spheres with square shoulder repulsions between unlike species, while the nanoparticles are modeled as hard-spheres with a square shoulder repulsion with one of the copolymer blocks. The resulting phase diagrams are presented for composites containing nanoparticles of various sizes and interaction strengths and include lamellae, perforated lamellae, cylinders, and disordered phases. Composites containing large nanoparticles also exhibit two-phase coexistence between different copolymer phases or between a copolymer phase and a nanoparticle phase, depending on the nanoparticle interaction strength. We also present concentration profiles perpendicular to the lamellar interface for nanoparticles of different sizes and interaction strengths. Neutral nanoparticles concentrate at the interface between copolymer domains, while interacti...

Molecular modeling of isomer effects in naphthenic and aromatic hydrocarbons

The vapor pressure and liquid density of 12 aromatic and 11 naphthenic hydrocarbons are studied to characterize transferable step potential models. Discontinuous molecular dynamics simulations are combined with thermodynamic perturbation theory to provide the basis for molecular modeling. Since the attractive forces appear as perturbations in this approach, their characterization can be varied continuously, without re-simulating, to obtain optimal results in minimal time. A generalized grid search is used to optimize the attractive steps for multiple compounds and multiple site-types simultaneously.In general, the vapor pressure is correlated to 10% average absolute deviation (%AAD) and the liquid density to 4%AAD. These are slightly larger than the errors encountered in characterizing straight chain hydrocarbons, ethers, and alcohols. One complicating factor that may explain the larger errors in this instance is the subtle variability in molecular structure caused by the consideration of isomers with branching in different positions around the ring. The transferable potentials account for this variability in a configurational sense, but not in a quantum mechanical sense. A detailed study reveals the extent to which the vapor pressure of several isomers can be correlated qualitatively and quantitatively. A preliminary study suggests one manner in which quantum mechanical corrections can be implemented.

Effects of interactions among surfactants, water and oil on equilibrium configuration of surfactant-water-oil systems

The distribution and configuration of surfactants at interface in surfactant-water-oil systems have been investigated using discontinuous molecular dynamic simulations. There exists a certain equilibrium concentration of surfactants at interface for the systems with certain interactions among surfactant, water and oil. The interface length and equilibrium morphology of the systems are dependent on the equilibrium concentration of surfactants at interface and the total amount of surfactants. The interaction strengths among surfactant, water and oil determine the equilibrium concentration of surfactants at interface. Three typical configurations of surfactants at interface have been observed: ailene1 surfactant 0 molecules are perpendicular to the interface and arranged closely; ailene2 perpendicular to the interface and arranged at interval 0 of two particles; ailene3 lie down in the interface partly. 0

Phase behavior in model homopolymer/CO2 and surfactant/CO2 systems: discontinuous molecular dynamics simulations.

Discontinuous molecular dynamics simulations are performed on surfactant (HmTn)/solvent systems modeled as a mixture of single-sphere solvent molecules and freely jointed surfactant chains composed of m slightly solvent-philic head spheres (H) and n solvent-philic tail spheres (T), all of the same size. We use a square-well potential to account for the head-head, head-solvent, tail-tail, and tail-solvent interactions and a hard-sphere potential for the head-tail and solvent-solvent interactions. We first simulate homopolymer/supercritical CO2 (scCO2) systems to establish the appropriate interaction parameters for a surfactant/scCO2 system. Next, we simulate surfactant/scCO2 systems and explore the effect of the surfactant volume fraction, packing fraction, and temperature on the phase behavior. The transition from the two-phase region to the one-phase region is located by monitoring the contrast structure factor of the equilibrated surfactant/scCO2 system, and the micelle to unimer transition is located by monitoring the aggregate size distribution of the equilibrated surfactant/scCO2 system. We find a two-phase region, a micelle phase, and a unimer phase with increasing packing fraction at fixed temperature or with increasing temperature at fixed packing fraction. The phase diagram for the surfactant/scCO2 system in the surfactant volume fraction-packing fraction plane and the density dependence of the critical micelle concentration are in qualitative agreement with experimental observations. The phase behavior of a surfactant/scCO2 system can be directly related to the solubilities of the corresponding homopolymers that serve as the head and tail blocks for the surfactant. The influence of surfactant structure (head and tail lengths) on the phase transitions is explored.

Transferable Step Potentials for the Straight-Chain Alkanes, Alkenes, Alkynes, Ethers, and Alcohols

Discontinuous molecular dynamics simulation and thermodynamic perturbation theory have been used to study thermodynamic and transport properties for hydrocarbons and oxygenated compounds. The fundamental basis of the method relies on a stepwise characterization of the disperse interactions and Wertheim potentials for the hydrogen bonding. The study on the transferable step potential energies, which are used to correlate the physical properties, the vapor pressure and liquid density, from available experimental data in the literature for straight-chain alkanes is extended to straight-chain alkenes, alkynes, ethers, and alcohols. The transferable multistep potentials for the new training set of pseudo-atoms governed by two transferable square-well depths and diameters are introduced, being transferable for the selected model compounds in each family. Those transferable parameters are applied to a set of validation compounds, and the results are analyzed. Overall, the vapor pressures and liquid densities are...

Effect of the interaction strength and range on the wetting behavior in confined binary fluid

The effects of interaction strength and range on the wetting behavior in two-dimensional binary fluid (which is composed of water-like and water–oil particles) confined between two parallel plane walls have been investigated by using the free energy analysis and discontinuous molecular dynamic simulations. Some rules are found. When the concentration of water-like particles c w is comparable to that of oil-like particles c o, the concave-water-lens may exist. As the strength or range of the preferential attraction between the walls and water-like particles increase, the system undergoes the transition from the partial wetting to complete wetting, and the occupied rate of water-like particles on the wall η varies gradually from a low fractional value to unity. When c w is very smaller or larger than c o, the concave-water-lens exist no long, and η increases abruptly from a low fractional value to unity.

Assembly and Kinetic Folding Pathways of a Tetrameric β-Sheet Complex: Molecular Dynamics Simulations on Simplified Off-Lattice Protein Models

We have performed discontinuous molecular dynamic simulations of the assembly and folding kinetics of a tetrameric β-sheet complex that contains four identical four-stranded antiparallel β-sheet peptides. The potential used in the simulation is a hybrid Go-type potential characterized by the bias gap parameter g, an artificial measure of a model protein's preference for its native state, and the intermolecular contact parameter η, which measures the ratio of intermolecular to intramolecular native attractions. The formation of the β-sheet complex and its equilibrium properties strongly depend on the size of the intermolecular contact parameter η. The ordered β-sheet complex in the folded state and nonaligned β-sheets or tangled chains in the misfolded state are distinguished by measuring the squared radius of gyration R g 2 and the fraction of native contacts Q. The folding yield for the folded state is high at intermediate values of η, but is low at both small and large values of η. The folded state at small η is liquid-like, but is solid-like at both intermediate and large η. The misfolded state at small η contains nonaligned β-sheets and tangled chains with poor secondary structure at large η. Various folding pathways via dimeric and trimeric intermediates are observed, depending on η. Comparison with experimental results on protein aggregation indicates that intermediate η values are most appropriate for modeling fibril formation and small η values are most appropriate for modeling the formation of amorphous aggregates.

Thermodynamics and stability of a beta-sheet complex: molecular dynamics simulations on simplified off-lattice protein models.

We have performed discontinuous molecular dynamics simulations of the thermodynamics and stability of a tetrameric beta-sheet complex that contains four identical four-stranded antiparallel beta-sheet peptides. The potential used in the simulation is a hybrid Go-type potential characterized by the bias gap parameter g, an artificial measure of the preference of a model protein for its native state, and the intermolecular contact parameter eta, which measures the ratio of intermolecular to intramolecular native attractions. Despite the simplicity of the model, a complex set of thermodynamic transitions for the beta-sheet complex is revealed that shows there are three distinct oligomer (partially ordered, ordered, and highly ordered beta-sheet complex) states and four noninteracting monomers phases. The thermodynamic properties of the three oligomer states strongly depend on both the size of the intermolecular contact parameter eta and the temperature. The partially ordered beta-sheet complex is made up of four ordered globules and is observed at intermediate to large eta at high temperatures. The ordered beta-sheet complex contains four native beta-sheets and is located at small to intermediate eta at low temperatures in the phase diagram. The highly ordered beta-sheet complex has fully-stiff beta-sheet strands, the same as the global energy minimum structure, and is observed for all eta at low temperatures.