On the role of inherent structures in glass-forming materials: I. The vitrification process

In this work, we investigate the role of inherent structures in the vitrification process of glass-forming materials by using a two-component Lennard-Jones mixture. We start with a simplified model that describes the dynamics of the atomistic system as a Poisson process consisting of a series of transitions from one potential energy minimum (inherent structure) to another, the rate of individual transitions being described by a first-order kinetic law. We investigate the validity of this model by comparing the mean square displacement resulting from atomistic molecular dynamics (MD) trajectories with the corresponding mean square displacement based on inherent structures. Furthermore, in the case of vitrification via stepwise cooling, we identify the role of the potential energy landscape in determining the properties of the resulting glass. Interestingly, the cooling rate is not sufficient to define the resulting glass in a stepwise cooling process, because the time spent by the system at different temperatures (length of the steps) has a highly nonlinear impact on the properties of the resulting glass. In contrast to previous investigations of supercooled liquids, we focus on a range of temperatures close to and below the glass transition temperature, where the use of MD is incapable of producing equilibrated samples of the metastable supercooled state. Our aim is to develop a methodology that enables mapping the dynamics under these conditions to a coarse-grained first-order kinetic model based on the Poisson process approximation. This model can be used in order to extend our dynamical sampling ability to much broader time scales and therefore allow us to create computer glasses with cooling rates closer to those used experimentally. In a continuation to this work, we provide the mathematical formulation for lifting the coarse-grained Poisson process model and reproducing the full dynamics of the atomistic system.     

Insights into the dynamics of HIV-1 protease: a kinetic network model constructed from atomistic simulations

The conformational dynamics in the flaps of HIV-1 protease plays a crucial role in the mechanism of substrate binding. We develop a kinetic network model, constructed from detailed atomistic simulations, to determine the kinetic mechanisms of the conformational transitions in HIV-1 PR. To overcome the time scale limitation of conventional molecular dynamics (MD) simulations, our method combines replica exchange MD with transition path theory (TPT) to study the diversity and temperature dependence of the pathways connecting functionally important states of the protease. At low temperatures the large-scale flap opening is dominated by a small number of paths; at elevated temperatures the transition occurs through many structurally heterogeneous routes. The expanded conformation in the crystal structure 1TW7 is found to closely mimic a key intermediate in the flap-opening pathways at low temperature. We investigated the different transition mechanisms between the semi-open and closed forms. The calculated relaxation times reveal fast semi-open ? closed transitions, and infrequently the flaps fully open. The ligand binding rate predicted from this kinetic model increases by 38-fold from 285 to 309 K, which is in general agreement with experiments. To our knowledge, this is the first application of a network model constructed from atomistic simulations together with TPT to analyze conformational changes between different functional states of a natively folded protein.     

Lumping analysis for the prediction of long-time dynamics: from monomolecular reaction systems to inherent structure dynamics of glassy materials

In this work we develop, test, and implement a methodology that is able to perform, in an automated manner, "lumping" of a high-dimensional, discrete dynamical system onto a lower-dimensional space. Our aim is to develop an algorithm which, without any assumption about the nature of the system's slow dynamics, is able to reproduce accurately the long-time dynamics with minimal loss of information. Both the original and the lumped systems conform to master equations, related via the "lumping" analysis introduced by Wei and Kuo [Ind. Eng. Chem. Fundam. 8, 114 (1969)], and have the same limiting equilibrium probability distribution. The proposed method can be used in a variety of processes that can be modeled via a first order kinetic reaction scheme. Lumping affords great savings in the computational cost and reveals the characteristic times governing the slow dynamics of the system. Our goal is to approach the best lumping scheme with respect to three criteria, in order for the lumped system to be able to fully describe the long-time dynamics of the original system. The criteria used are: (a) the lumping error arising from the reduction process; (b) a measure of the magnitude of singular values associated with long-time evolution of the lumped system; and (c) the size of the lumped system. The search for the optimum lumping proceeds via Monte Carlo simulation based on the Wang-Landau scheme, which enables us to overcome entrapment in local minima in the above criteria and therefore increases the probability of encountering the global optimum. The developed algorithm is implemented to reproduce the long-time dynamics of a glassy binary Lennard-Jones mixture based on the idea of "inherent structures," where the rate constants for transitions between inherent structures have been evaluated via hazard plot analysis of a properly designed ensemble of molecular dynamics trajectories.     

Transferability of coarse-grained force fields: the polymer case

A key question for all coarse-graining methodologies is the degree of transferability of the resulting force field between various systems and thermodynamic conditions. Here we present a detailed study of the transferability over different thermodynamic states of a coarse-grained (CG) force field developed using the iterative Boltzmann inversion method. The force field is optimized against distribution functions obtained from atomistic simulations. We analyze the polymer case by investigating the bulk of polystyrene and polyamide-6,6 whose coarse-grained models differ in the chain length and in the number of atoms lumped in one bead. The effect of temperature and pressure on static, dynamic, and thermodynamic properties is tested by comparing systematically the coarse-grain results with the atomistic ones. We find that the CG model describing the polystyrene is transferable only in a narrow range of temperature and it fails in describing the change of the bulk density when temperature is 80 K lower than the optimization one. Moreover the calculation of the self-diffusion coefficient shows that the CG model is characterized by a faster dynamics than the atomistic one and that it overestimates the isothermal compressibility. On the contrary, the polyamide-6,6 CG model turns out to be fully transferable between different thermodynamic conditions. The transferability is checked by changing either the temperature or the pressure of the simulation. We find that, in this case, the CG model is able to follow all the intra- and interstructural rearrangements caused by the temperature changes. In addition, while at low temperature the difference between the CG and atomistic dynamics is remarkable due to the presence of hydrogen bonds in the atomistic systems, for high temperatures, the speedup of the CG dynamics is strongly reduced, leading to a CG diffusion coefficient only six times bigger than the atomistic one. Moreover, the isothermal compressibility calculated at different temperatures agrees very well with the experimental one. We find that the polymer chain length does not affect the transferability of the force field and we attribute such transferability mainly to the finer model used in describing the polyamide-6,6 than the polystyrene.     

Modelling and simulation of light-activated membrane proteins: dynamical transitions in bacteriorhodopsin

Many of the functions of membranes are carried out by proteins associated with them. A knowledge of atomic-detail membrane protein structures and dynamics is required for a full understanding of these functions. We briefly discuss recent progress in this field using modelling and simulation. One of the best characterised membrane proteins, bacteriorhodopsin, undergoes dynamical transitions with temperature. Here we present preliminary results of molecular dynamics simulation of this protein as a function of temperature, indicating the presence of dynamical transitions at approximately the temperatures seen experimentally.     

Quantum dynamics of vibration–vibration energy transfer for vibrationally excited HF colliding with H2

The rate constants for H₂–HF energy transfer processes, especially for those in vibrationally excited states, are very demanding in astrophysics and chemical laser engineering, especially for those in vibrationally excited states. Based on our recent potential energy surface, we used the coupled-states approximation including the nearest neighboring Coriolis couplings with energy-based local basis set to perform dynamics calculation for the H₂–HF energy transfer system. Rate constants for vibrational transitions (1; 3) → (0; 4), (1; 3) → (2; 2), and (0; 3) → (1; 2) were obtained. For state-to-state rate constants, transitions that have no internal angular momentum gap dominate at high temperatures. The vibrational-resolved rate constant for (1; 3) → (0; 4) initially decreases and then increases with the temperature, while those for (1; 3) → (2; 2), and (0; 3) → (1; 2) transitions monotonically increase. The calculated rate constants are in good agreement with the available experimental results. These dynamical data can be further applied to the numerical simulation of hydrogen fluoride chemical laser. © 2018 Wiley Periodicals, Inc.     

First order melting transitions of highly ordered dipalmitoyl phosphatidylcholine gel phase membranes in molecular dynamics simulations with atomistic detail

Molecular dynamics simulations with atomistic detail of the gel phase and melting transitions of dipalmitoyl phosphatidylcholine bilayers in water reveal the dependency of many thermodynamic and structural parameters on the initial system ordering. We quantitatively compare different methods to create a gel phase system and we observe that a very high ordering of the gel phase starting system is necessary to observe behavior which reproduces experimental data. We performed heating scans with speeds down to 0.5 K/ns and could observe sharp first order phase transitions. Also, we investigated the transition enthalpy as the natural intrinsic parameter of first order phase transitions, and obtained a quantitative match with experimental values. Furthermore, we performed systematic investigations of the statistical distribution and heating rate dependency of the microscopic phase transition temperature.     

The origin of persistent shear stress in supercooled liquids

The persistence of shear stress fluctuations in viscous liquids is a direct consequence of the non-zero shear stress of the local potential minima which couples stress relaxation to transitions between inherent structures. In simulations of 2D and 3D glass forming mixtures, we calculate the distribution of this inherent shear stress and demonstrate that the variance is independent of temperature and obeys a power law in density. The inherent stress is shown to involve only long wavelength fluctuations, evidence of the central role of the static boundary conditions in determining the residual stress left after the minimization of the potential energy. A temperature T(η) is defined to characterise the crossover from stress relaxation governed by binary collisions at high temperatures to low temperature relaxation dominated by the relaxation of the inherent stress. T(η) is found to coincide with the breakdown of the Stokes-Einstein scaling of diffusion and viscosity.     

A molecular dynamics simulation of the melting points and glass transition temperatures of myo- and neo-inositol

The heat of sublimation, density, melting point, and glass transition temperature are calculated for myo- and neo-inositol, using the condensed-phase optimized molecular potentials for atomistic simulation studies (COMPASS) force field and molecular dynamics techniques. Our results show that the calculated heats of sublimation and density are very close to the experimental values for both compounds. Furthermore, our simulated melting temperatures for myo- and neo-inositol also compare very well to the experimentally obtained data. The glass transition temperatures for myo- and neo-inositol have been calculated to be ca. 494 K and ca. 518 K, respectively, and the shape of the volume versus temperature plots produced are typical for a glass transition. As a result, it is our view that the COMPASS force field suitably describes these two compounds in molecular simulations and that molecular dynamics techniques, combined with this force field, can be used to simulate the melt and glass transitions for such molecules.     

Dynamical integration of a Markovian web: a first passage time approach

In this work we address the dynamics of Markovian systems by tracking the evolution of the probability distribution, utilizing mean first passage time theory to augment the set of states considered. The method is validated on a lattice system and is applied, in conjunction with landscape analysis (saddle point searches) and multidimensional transition-state theory, to an atomistic model of glassy atactic polystyrene, in order to follow its time evolution over more than ten orders of magnitude on the time scale, from less than 10(-15) up to 10(-5) s. Frequencies extracted from the eigenvalues of the rate constant matrix are in favorable agreement with experimental measurements of subglass relaxation transitions at 250 K.     

A unified theory for charge-carrier transport in organic crystals

To characterize the crossover from bandlike transport to hopping transport in molecular crystals, we study a microscopic model that treats electron-phonon interactions explicitly. A finite-temperature variational method combining Merrifield's transformation with Bogoliubov's theorem is developed to obtain the optimal basis for an interacting electron-phonon system, which is then used to calculate the bandlike and hopping mobilities for charge carriers. Our calculations on the one dimensional (1D) Holstein model at T=0 K and finite temperatures show that the variational basis gives results that compared favorably to other analytical methods. We also study the structures of polaron states at a broad range of parameters including different temperatures. Furthermore, we calculate the bandlike and hopping mobilities of the 1D Holstein model in different parameters and show that our theory predicts universal power-law decay at low temperatures and an almost temperature independent behavior at higher temperatures, in agreement with experimental observations. In addition, we show that as the temperature increases, hopping transport can become dominant even before the polaron state changes its character. Thus, our result indicates that the self-trapping transition studied in conventional polaron theories does not necessarily correspond to the bandlike to hopping transition in the transport properties in organic molecular crystals. Finally, a comparison of our 1D results with experiments on ultrapure naphthalene crystals suggests that the theory can describe the charge-carrier mobilities quantitatively across the whole experimental temperature range.     

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.     

Ab initio nonadiabatic dynamics study of ultrafast radiationless decay over conical intersections illustrated on the Na3F cluster

We present a theoretical approach for the ultrafast nonadiabatic dynamics based on the ab initio molecular dynamics carried out "on the fly" in the framework of the configuration interaction method combined with Tully's surface hopping algorithm for nonadiabatic transitions. This approach combined with our Wigner distribution approach allows us to perform accurate simulations of femtosecond pump-probe spectra in the systems where radiationless transitions among electronic states take place. In this paper we illustrate this by theoretical simulation of ultrafast processes and nonradiative relaxation in the Na(3)F cluster, involving three excited states and the ground electronic state. Furthermore, we show that our accurate simulation of the photoionization pump-probe spectrum is in full agreement with the experimental signal. Based on the nonadiabatic dynamics at high level of accuracy and taking into account all degrees of freedom, the nonradiative lifetime for the 1 (1)B(1) excited state of Na(3)F has been determined to be approximately 900 fs.     

Optical properties of alkaline-earth metal oxides from first principles

This study reports the results of an ab initio electronic and optical calculation of alkaline-earth metal oxides (MgO, CaO, SrO and BaO) in the NaCl crystal structure using the full potential linearized augmented plane wave (FP-LAPW) method within the density functional theory. The exchange-correlation potential is treated by the generalized gradient approximation within the Perdew et al scheme. Moreover, the Engel–Vosko GGA formalism is applied so as to optimize the corresponding potential for band structure calculations. The real and imaginary parts of the dielectric function ?(ω), the optical absorption coefficient I(ω), the reflectivity R(ω) and the energy loss function are calculated by random phase approximation (RPA). The calculated results show a qualitative agreement with the available experimental results in the sense that we can recognize some peaks qualitatively, those due to single particle transitions. Furthermore the interband transitions responsible for the structures in the spectra are specified. It is shown that the oxygen 2p states and metal d states play the major role in optical transitions as initial and final states respectively. The effect of the spin–orbit coupling on the optical properties is also investigated and found to be quite small, especially in the low energy region. The dielectric constants are calculated and compared with the available theoretical and experimental results.     

Local chain dynamics of a model polycarbonate near glass transition temperature: A molecular dynamics simulation

Constant pressure constant temperature molecular dynamics method is employed to investigate the atomistic scale dynamics of a model Bisphenol A polycarbonate in the vicinity of its glass transition temperature. First, the glass transition temperature and the thermal expansion coefficients of the polymer are predicted by performing simulations at different temperatures. To explore the significance of different modes of motion, various types of time correlation functions are utilized in analyzing the trajectories. In these nanosecond scale simulations, the motion of the chain segments is found to be highly localized with little reorientation of the vectors representing these segments. Detailed analysis of trajectories and the correlation functions of the backbone dihedrals and side methyl groups indicates that they exhibit numerous conformational transitions. The activation energies of the conformational transitions obtained from the simulation are generally larger than the potential barriers for the rotations of these dihedrals, however, both show the same trend. We also have estimated the phenylene ring flip activation energy as 12.6 kcal/mol and the flip frequency as 0.77 MHz at 300 K. These values fall either fall within the range determined by various NMR spectroscopy experiments or slightly out of the range. The study shows that the conformational transitions between the adjacent dihedrals are strongly correlated. Three basic cooperative modes are identified from the simulation. They are: a positive synchronous rotation of two phenylene rings, a negative synchronous rotation of two phenylene rings, and a carbonate group rotation. Above the glass transition temperature, the large scale cooperative motions become much more significant.     

Critical dynamics of dimers: implications for the glass transition

The Adam-Gibbs view of the glass transition relates the relaxation time to the configurational entropy, which goes continuously to zero at the so-called Kauzmann temperature. We examine this scenario in the context of a dimer model with an entropy-vanishing phase transition and stochastic loop dynamics. We propose a coarse-grained master equation for the order parameter dynamics which is used to compute the time-dependent autocorrelation function and the associated relaxation time. Using a combination of exact results, scaling arguments, and numerical diagonalizations of the master equation, we find nonexponential relaxation and a Vogel-Fulcher divergence of the relaxation time in the vicinity of the phase transition. Since in the dimer model the entropy stays finite all the way to the phase transition point and then jumps discontinuously to zero, we demonstrate a clear departure from the Adam-Gibbs scenario. Dimer coverings are the "inherent structures" of the canonical frustrated system, the triangular Ising antiferromagnet. Therefore, our results provide a new scenario for the glass transition in supercooled liquids in terms of inherent structure dynamics.     

Freezing and folding behavior in simple off-lattice heteropolymers

We have performed parallel tempering Monte Carlo simulations using a simple continuum heteropolymer model for proteins. All 10 heteropolymer sequences which we have studied have shown first-order transitions at low temperature to ordered states dominated by single chain conformations. These results are in contrast with the theoretical predictions of the random energy model for heteropolymers, from which we would expect continuous transitions to glassy behavior at low temperatures.     

Robustness of atomistic Gō models in predicting native-like folding intermediates

Go? models are exceedingly popular tools in computer simulations of protein folding. These models are native-centric, i.e., they are directly constructed from the protein's native structure. Therefore, it is important to understand up to which extent the atomistic details of the native structure dictate the folding behavior exhibited by Go? models. Here we address this challenge by performing exhaustive discrete molecular dynamics simulations of a Go? potential combined with a full atomistic protein representation. In particular, we investigate the robustness of this particular type of Go? models in predicting the existence of intermediate states in protein folding. We focus on the N47G mutational form of the Spc-SH3 folding domain (x-ray structure) and compare its folding pathway with that of alternative native structures produced in silico. Our methodological strategy comprises equilibrium folding simulations, structural clustering, and principal component analysis.     

Fluorescence probing of temperature-dependent dynamics and friction in ionic liquid local environments

The solvation dynamics and local orientational friction for a series of four ionic liquids have been probed using coumarin 153 (C153) as a function of temperature. These ionic liquids are comprised of nonaromatic organic cations paired with a common anion, bis(trifluoromethylsulfonyl)imide (NTf(2)-). The specific liquids are as follows: N-methyl-tri-N-butylammonium NTf(2)- (N(1444)+/NTf(2)-), N-hexyl-tri-N-butylammonium NTf(2)- (N(6444)+/NTf(2)-), N-methyl-N-butylpyrrolidinium NTf(2)- (Pyrr(14)+/NTf(2)-), and N-methyl-N-ethoxyethylpyrrolidinium NTf(2)- (Pyrr(1(2O2))+/NTf(2)-). The observed solvation dynamics and fluorescence depolarization dynamics occur over a broad range of time scales that can only be adequately fit by functions including three or more exponential components. Stretched exponential distributions cannot adequately fit our data. The solvation and reorientational dynamics of the C153 probe are studied over a range of temperatures from 278.2 to 353.2 K. For both the solvation dynamics and the probe reorientational dynamics, the observed temperature dependence is well fit by a Vogel-Tammann-Fulcher law. To correlate the observed microscopic dynamics with macroscopic physical properties, temperature-dependent viscosities are also measured. Differential scanning calorimetry is used to study the thermodynamics of the phase transitions from the liquid to supercooled liquid to glassy states. For the two tetraalkylammonium liquids, the observed melting transitions occur near 300 K, so we are able to study the dynamics in a clearly supercooled regime. Very long time scale orientational relaxation time constants dynamics on the order of 100 ns are observed in the C153 fluorescence anisotropy. These are interpreted to arise from long-lived local structures in the environment surrounding the C153 probe.