Kazushi Fujimoto

MODYLAS: A Highly Parallelized General-Purpose Molecular Dynamics Simulation Program for Large-Scale Systems with Long-Range Forces Calculated by Fast Multipole Method (FMM) and Highly Scalable Fine-Grained New Parallel Processing Algorithms

Our new molecular dynamics (MD) simulation program, MODYLAS, is a general-purpose program appropriate for very large physical, chemical, and biological systems. It is equipped with most standard MD techniques. Long-range forces are evaluated rigorously by the fast multipole method (FMM) without using the fast Fourier transform (FFT). Several new methods have also been developed for extremely fine-grained parallelism of the MD calculation. The virtually buffering-free methods for communications and arithmetic operations, the minimal communication latency algorithm, and the parallel bucket-relay communication algorithm for the upper-level multipole moments in the FMM realize excellent scalability. The methods for blockwise arithmetic operations avoid data reload, attaining very small cache miss rates. Benchmark tests for MODYLAS using 65 536 nodes of the K-computer showed that the overall calculation time per MD step including communications is as short as about 5 ms for a 10 million-atom system; that is, 35 ns of simulation time can be computed per day. The program enables investigations of large-scale real systems such as viruses, liposomes, assemblies of proteins and micelles, and polymers.

All-atom molecular dynamics calculation study of entire poliovirus empty capsids in solution.

Small viruses that belong, for example, to the Picornaviridae, such as poliovirus and foot-and-mouth disease virus, consist simply of capsid proteins and a single-stranded RNA (ssRNA) genome. The capsids are quite stable in solution to protect the genome from the environment. Here, based on long-time and large-scale 6.5 × 10(6) all-atom molecular dynamics calculations for the Mahoney strain of poliovirus, we show microscopic properties of the viral capsids at a molecular level. First, we found equilibrium rapid exchange of water molecules across the capsid. The exchange rate is so high that all water molecules inside the capsid (about 200,000) can leave the capsid and be replaced by water molecules from the outside in about 25 μs. This explains the capsids tolerance to high pressures and deactivation by exsiccation. In contrast, the capsid did not exchange ions, at least within the present simulation time of 200 ns. This implies that the capsid can function, in principle, as a semipermeable membrane. We also found that, similar to the xylem of trees, the pressure of the solution inside the capsid without the genome was negative. This is caused by coulombic interaction of the solution inside the capsid with the capsid excess charges. The negative pressure may be compensated by positive osmotic pressure by the solution-soluble ssRNA and the counter ions introduced into it.

Free energy profiles for penetration of methane and water molecules into spherical sodium dodecyl sulfate micelles obtained using the thermodynamic integration method combined with molecular dynamics calculations

The free energy profiles, ΔG(r), for penetration of methane and water molecules into sodium dodecyl sulfate (SDS) micelles have been calculated as a function of distance r from the SDS micelle to the methane and water molecules, using the thermodynamic integration method combined with molecular dynamics calculations. The calculations showed that methane is about 6-12 kJ mol(-1) more stable in the SDS micelle than in the water phase, and no ΔG(r) barrier is observed in the vicinity of the sulfate ions of the SDS micelle, implying that methane is easily drawn into the SDS micelle. Based on analysis of the contributions from hydrophobic groups, sulfate ions, sodium ions, and solvent water to ΔG(r), it is clear that methane in the SDS micelle is about 25 kJ mol(-1) more stable than it is in the water phase because of the contribution from the solvent water itself. This can be understood by the hydrophobic effect. In contrast, methane is destabilized by 5-15 kJ mol(-1) by the contribution from the hydrophobic groups of the SDS micelle because of the repulsive interactions between the methane and the crowded hydrophobic groups of the SDS. The large stabilizing effect of the solvent water is higher than the repulsion by the hydrophobic groups, driving methane to become solubilized into the SDS micelle. A good correlation was found between the distribution of cavities and the distribution of methane molecules in the micelle. The methane may move about in the SDS micelle by diffusing between cavities. In contrast, with respect to the water, ΔG(r) has a large positive value of 24-35 kJ mol(-1), so water is not stabilized in the micelle. Analysis showed that the contributions change in complex ways as a function of r and cancel each other out. Reference calculations of the mean forces on a penetrating water molecule into a dodecane droplet clearly showed the same free energy behavior. The common feature is that water is less stable in the hydrophobic core than in the water phase because of the energetic disadvantage of breaking hydrogen bonds formed in the water phase. The difference between the behaviors of the SDS micelles and the dodecane droplets is found just at the interface; this is caused by the strong surface dipole moment formed by sulfate ions and sodium ions in the SDS micelles.

Why Is Poly(oxyethylene) Soluble in Water? Evidence from the Thermodynamic Profile of the Conformational Equilibria of 1,2-Dimethoxyethane and Dimethoxymethane Revealed by Raman Spectroscopy

The origin of high solubility of poly(oxyethylene) in water has been an open question. Although it is thought that the high solubility of poly(oxyethylene) arises from an increase of the trans-gauche-trans (tgt) conformer in water, the relationship between the increase of the tgt conformer and the solubility is unclear. In this study, we have investigated the conformational equilibria of 1,2-dimethoxyethane, which is a model molecule for poly(oxyethylene), by using Raman spectroscopy, and determined the change in the population and the free energy of each conformer with the aid of density functional theory calculations. The free energy of transfer of the tgt conformer from the pure liquid to the water phase is -6.1 ± 0.2 kJ mol(-1). Furthermore, the fraction of the tgt conformer increases from 0.37 to 0.78. Thus, the net contribution of the tgt conformer is -4.8 ± 0.2 kJ mol(-1), which is 79% of the total free energy of transfer (-6.07 kJ mol(-1)). This demonstrates that the high solubility of 1,2-dimethoxyethane originates from the lowest free energy and the highest fraction of the tgt conformer in water. We also successfully explain the thermodynamic mechanism of the low solubility of dimethoxymethane, which is the model molecule for poly(oxymethylene).

Molecular dynamics study of solubilization of immiscible solutes by a micelle: Free energy of transfer of alkanes from water to the micelle core by thermodynamic integration method

Free energy of transfer, DeltaG(w-->m), from water phase to a sodium dodecyl sulfate (SDS) micelle core has been calculated for a series of hydrophobic solutes originally immiscible with water by thermodynamic integration method combined with molecular dynamics calculations. The calculated free energy of transfer is in good correspondence to the experiment as well as the theoretical free energy of transfer. The calculated DeltaG(w-->m)s are all negative, implying that the alkane molecules are more stable in the micelle than in the water phase. It decreases almost linearly as a function of the number of carbon atoms of the alkanes longer than methane with a decrement of 3.3 kJ mol(-1) per one methylene group. The calculated free energy of transfer indicates that, for example, at the micelle concentration of 50 CMC (critical micelle concentration), about only 1 of 6 micelles or 1 of 32 000 micelles does not contain a solute methane or n-octane molecule, respectively.

Molecular dynamics study of the potential of mean force of SDS aggregates

In our previous study, all-atomistic molecular dynamics (MD) calculations have been carried out for the aggregation of ionic sodium dodecyl sulfate in water [S. Kawada et al., Chem. Phys. Lett. 646, 36 (2016)]. Aggregates of 20-30 dodecyl sulfate ions were formed within a short MD run for 10 ns. However, further aggregation did not occur despite a long MD calculation for more than 100 ns. This suggests that strong electrostatic repulsive interactions between the aggregates prevent the fusion of the aggregates. In the present study, mean force and potential of mean force acting between two aggregates with aggregation number N = 30 were evaluated as a function of their separation by MD calculations. The repulsive force becomes strong with decreasing distance between the two aggregates before they merge into one. An origin of the repulsive force is an electric double layer formed by the sulfate group and counter sodium ions. Strength of the repulsive force is in good agreement with the theoretical value given by the Derjaguin-Landau-Verwey-Overbeek (DLVO) theory. Once the aggregates establish contact, the force between them turns to be a large attractive force that can be explained by the interfacial tension. In order to form a single micelle from the two aggregates, it is necessary for them to climb over a free energy barrier of 23 kJ/mol. Once, the barrier is overcome, the micelle is stabilized by ∼200 kJ/mol. The time constant of aggregation evaluated from the calculated free energy barrier was about 28 μs at the concentration in our previous study.

Detailed Structural Analysis of a Self-Assembled Vesicular Amphiphilic NCN-Pincer Palladium Complex by Using Wide-Angle X-Ray Scattering and Molecular Dynamics Calculations.

Wide-angle X-ray scattering experiments and all-atomistic molecular dynamics calculations were performed to elucidate the detailed structure of bilayer vesicles constructed by self-assembly of an amphiphilic palladium NCN-pincer complex. We found an excellent agreement between the experimental and calculated X-ray spectra, and between the membrane thickness determined from a TEM image and that calculated from an electron-density profile, which indicated that the calculated structure was highly reliable. The analysis of the simulated bilayer structure showed that in general the membrane was softer than other phospholipid bilayer membranes. In this bilayer assemblage, the degree of alignment of complex molecules in the bilayer membrane was quite low. An analysis of the electron-density profile shows that the bilayer assemblage contains a space through which organic molecules can exit. Furthermore, the catalytically active center is near this space and is easily accessible by organic molecules, which permits the bilayer membrane to act as a nanoreactor. The free energy of permeation of water through the bilayer membrane of the amphiphilic complex was 12 kJ mol-1 , which is much lower than that for phospholipid bilayer membranes in general. Organic molecules are expected to pass though the bilayer membrane. The self-assembled vesicles were shown to be catalytically active in a Miyaura-Michael reaction in water.

A molecular dynamics study of the breathing and deforming modes of the spherical ionic SDS and nonionic C12E8 micelles

In order to investigate shape of the micelles and its thermal fluctuations, molecular dynamics calculations have been performed for spherical ionic sodium dodecyl sulfate (SDS) and nonionic octaethyleneglycol monododecyl ether (C12E8) micelles. New statistical functions suitable for extracting the fluctuations of the shape of the spherical micelles were defined using spherical harmonics and Legendre polynomials. The breathing and deforming modes of the SDS and C12E8 micelles were analyzed in detail based on these new functions. The elastic nature of the micelle core was also discussed. The present analysis gives a new molecular picture that the micelle shape is a superposition of the various kinds of breathing and deforming modes, and each mode has a specific relaxation time of the shape fluctuation.

Molecular dynamics study of free energy of transfer of alcohol and amine from water phase to the micelle by thermodynamic integration method

Free energy of transfer of methylamine, octylamine, methanol, and octanol from water phase to sodium dodecyl sulfate (SDS) micelle has been calculated using thermodynamic integration method combined with molecular dynamics calculations. Together with the results for alkanes obtained in our previous study [K. Fujimoto, N. Yoshii, and S. Okazaki, J. Chem. Phys. 133, 074511 (2010)], the effect of polar group on the partition of hydrophilic solutes between water phase and the micelle has been investigated in detail at a molecular level. The calculations showed that the molecules with octyl group are more stable in the SDS micelle than in the water phase due to their hydrophobicity of long alkyl chain. In contrast, methanol and methylamine are stable in the water phase as well as in the micelle because of their high hydrophilicity. The spatial distribution of methylamine, octylamine, methanol, and octanol has also been evaluated as a function of the distance, R, from the center of mass of SDS micelle to the solutes. The distribution shows that the methylamine molecule is adsorbed on the SDS micelle surface, while the methanol molecule is delocalized among the whole system, i.e., in the water phase, on the surface of the micelle, and in the hydrophobic core of the micelle. The octylamine and octanol molecules are solubilized in the SDS micelle with palisade layer structure and are not found in the water phase.

G-Protein/β-Arrestin-Linked Fluctuating Network of G-Protein-Coupled Receptors for Predicting Drug Efficacy and Bias Using Short-Term Molecular Dynamics Simulation.

The efficacy and bias of signal transduction induced by a drug at a target protein are closely associated with the benefits and side effects of the drug. In particular, partial agonist activity and G-protein/β-arrestin-biased agonist activity for the G-protein-coupled receptor (GPCR) family, the family with the most target proteins of launched drugs, are key issues in drug discovery. However, designing GPCR drugs with appropriate efficacy and bias is challenging because the dynamic mechanism of signal transduction induced by ligand—receptor interactions is complicated. Here, we identified the G-protein/β-arrestin-linked fluctuating network, which initiates large-scale conformational changes, using sub-microsecond molecular dynamics (MD) simulations of the β2-adrenergic receptor (β2AR) with a diverse collection of ligands and correlation analysis of their G protein/β-arrestin efficacy. The G-protein-linked fluctuating network extends from the ligand-binding site to the G-protein-binding site through the connector region, and the β-arrestin-linked fluctuating network consists of the NPxxY motif and adjacent regions. We confirmed that the averaged values of fluctuation in the fluctuating network detected are good quantitative indexes for explaining G protein/β-arrestin efficacy. These results indicate that short-term MD simulation is a practical method to predict the efficacy and bias of any compound for GPCRs.