Three pillars for achieving quantum mechanical molecular dynamics simulations of huge systems: Divide‐and‐conquer,density‐functional tight‐binding,and massively parallel computation

The linear‐scaling divide‐and‐conquer (DC) quantum chemical methodology is applied to the density‐functional tight‐binding (DFTB) theory to develop a massively parallel program that achieves on‐the‐fly molecular reaction dynamics simulations of huge systems from scratch. The functions to perform large scale geometry optimization and molecular dynamics with DC‐DFTB potential energy surface are implemented to the program called DC‐DFTB‐K. A novel interpolation‐based algorithm is developed for parallelizing the determination of the Fermi level in the DC method. The performance of the DC‐DFTB‐K program is assessed using a laboratory computer and the K computer. Numerical tests show the high efficiency of the DC‐DFTB‐K program, a single‐point energy gradient calculation of a one‐million‐atom system is completed within 60 s using 7290 nodes of the K computer. © 2016 Wiley Periodicals, Inc.     

Possible mechanism of BN fullerene formation from a boron cluster: Density‐functional tight‐binding molecular dynamics simulations

We simulate the formation of a BN fullerene from an amorphous B cluster at 2000 K by quantum mechanical molecular dynamics based on the density‐functional tight‐binding method. We run 30 trajectories 200 ps in length, where N atoms are supplied around the target cluster, which is initially an amorphous B₃₆ cluster. Most of the incident N atoms are promptly incorporated into the target cluster to form B‐N‐B bridges or NB₃ pyramidal local substructures. BN fullerene formation is initiated by alternating BN ring condensation. Spontaneous atomic rearrangement and N₂ dissociation lead to the construction of an sp² single‐shelled structure, during which the BN cluster undergoes a transition from a liquid‐like to a solid‐like state. Continual atomic rearrangement and sporadic N₂ dissociation decrease the number of defective rings in the BN cluster and increase the number of six‐membered rings, forming a more regular shell structure. The number of four‐membered rings tends to remain constant, and contributes to more ordered isolated‐tetragon‐rule ring placement. © 2016 Wiley Periodicals, Inc.     

A theoretical study on cyclacenes: Analytical tight‐binding approach

We present a theoretical study of cyclacene molecules performed at tight‐binding level. The orbital energies and eigenvectors have been analytically computed, and exact expressions for the axial component of the total position spread and polarizability tensors have been obtained. In absence of dimerization, the system has a D_nh symmetry, where n is the number of hexagonal units. The energy bands present no gap at the Fermi level, and to this fact it corresponds a diverging (per‐electron) polarizability for in the direction of the system symmetry axis. The two (degenerate) components of the polarizability on the σ_h symmetry plane, conversely, remain finite for . The total position spread tensor presents a qualitatively different behavior, since all the three components of the position spread per electron remain finite for . The results are analyzed and discussed for both axial and planar components separately as these are affected differently with respect to the increasing system size. Both dipole polarizability and total position spread have been computed using an ab initio approach for the smallest systems, to compare the analytical tight‐binding expressions with a higher‐level theory.     

Three‐body expansion of the fragment molecular orbital method combined with density‐functional tight‐binding

The three‐body fragment molecular orbital (FMO3) method is formulated for density‐functional tight‐binding (DFTB). The energy, analytic gradient, and Hessian are derived in the gas phase, and the energy and analytic gradient are also derived for polarizable continuum model. The accuracy of FMO3‐DFTB is evaluated for five proteins, sodium cation in explicit solvent, and three isomers of polyalanine. It is shown that FMO3‐DFTB is considerably more accurate than FMO2‐DFTB. Molecular dynamics simulations for sodium cation in water are performed for 100 ps, yielding radial distribution functions and coordination numbers. © 2017 Wiley Periodicals, Inc.     

Atomistic insight into chondroitin‐6‐sulfate glycosaminoglycan chain through quantum mechanics calculations and molecular dynamics simulation

Chondroitin‐6‐sulfate (C6S) is a glycosaminoglycan (GAG) constituent in the extracellular matrix, which participates actively in crucial biological processes, as well as in various pathological conditions, such as atherosclerosis and cancer. Molecular interactions involving the C6S chain are therefore of considerable interest. A computational model for atomistic simulation was built. This work describes the design and validation of a force field for a C6S dodecasaccharide chain. The results of an extensive molecular dynamics simulation performed with the new force field provide a novel insight into the structure and dynamics of the C6S chain. The intramolecular H‐bonds in the disaccharide linkage region are suggested to play a major role in determining the chain structural dynamics. Moreover, the unravelling of an additional H‐bond involving the sulfate groups in C6S is interesting as changes in sulfation have been claimed to be an important factor in several diseases. The force field will prove useful for future studies of crucial interactions between C6S and various nanoassemblies. It can also be used as a basis for modeling of other GAGs. © 2010 Wiley Periodicals, Inc. J Comput Chem, 2010     

Deformation of single‐walled carbon nanotubes by interaction with graphene: A first‐principles study

The interaction between single‐walled carbon nanotubes (SWNTs) and graphene were studied with first‐principles calculations. Both SWNTs and single‐layer graphene (SLG) or double‐layer graphene (DLG) display more remarkable deformations with the increase of SWNT diameter, which implies a stronger interaction between SWNTs and graphene. Besides, in DLG, deformation of the upper‐layer graphene is less than in SLG. Zigzag SWNTs show stronger interactions with SLG than armchair SWNTs, whereas the order is reversed for DLG, which can be interpreted by the mechanical properties of SWNTs and graphene. Density of states and band structures were also studied, and it was found that the interaction between a SWNT and graphene is not strong enough to bring about obvious influence on the electronic structures of SWNTs. © 2015 Wiley Periodicals, Inc.     

Quantum chemistry study on the mechanism of the reaction between ozone and 2,3,7,8‐TCDD

Despite of its fundamental importance, the mechanism of the reaction between ozone and dioxins are still lack detailed investigation so far. It is well-known that quantum chemical calculation is a well-established method for investigating chemical reactions. In this article, quantum chemical calculation was employed to investigate the mechanism of the reaction between ozone and dioxins, as exemplified by 2,3,7,8-TCDD. The theoretical study showed that, 2,3,7,8-TCDD was gradually destructed by ozone via six cleavages of the CC bonds. All the six cleavages of the CC bonds were calculated and discussed in detail based on the theoretical calculations by the UB3LYP/6-31G(d) method. At the same time, the energies of stationary points along the reaction process were calculated by the UMP2/6-311g(d,p)//UB3LYP/6-31G(d) method and the activation energy was obtained. The obtained activation energy was 12.25 kcal/mol, which was lower than that of the reaction between benzene and O₃(16.64 kcal/mol). This indicated that, by comparison with benzene, 2,3,7,8-TCDD could be more efficiently destructed by O₃. The reason for this result was also discussed. © 2010 Wiley Periodicals, Inc. Int J Quantum Chem, 2011     

Rapid QM/MM approach for biomolecular systems under periodic boundary conditions: Combination of the density‐functional tight‐binding theory and particle mesh Ewald method

A quantum mechanical/molecular mechanical (QM/MM) approach based on the density‐functional tight‐binding (DFTB) theory is a useful tool for analyzing chemical reaction systems in detail. In this study, an efficient QM/MM method is developed by the combination of the DFTB/MM and particle mesh Ewald (PME) methods. Because the Fock matrix, which is required in the DFTB calculation, is analytically obtained by the PME method, the Coulomb energy is accurately and rapidly computed. For assessing the performance of this method, DFTB/MM calculations and molecular dynamics simulation are conducted for a system consisting of two amyloid‐β(1‐16) peptides and a zinc ion in explicit water under periodic boundary conditions. As compared with that of the conventional Ewald summation method, the computational cost of the Coulomb energy by utilizing the present approach is drastically reduced, i.e., 166.5 times faster. Furthermore, the deviation of the electronic energy is less than . © 2016 Wiley Periodicals, Inc.     

Conformational properties of penicillins: quantum chemical calculations and molecular dynamics simulations of benzylpenicillin

Herein, we present theoretical results on the conformational properties of benzylpenicillin, which are characterized by means of quantum chemical calculations (MP2/6-31G* and B3LYP/6-31G*) and classical molecular dynamics simulations (5 ns) both in the gas phase and in aqueous solution. In the gas phase, the benzylpenicillin conformer in which the thiazolidine ring has the carboxylate group oriented axially is the most favored one. Both intramolecular CH. O and dispersion interactions contribute to stabilize the axial conformer with respect to the equatorial one. In aqueous solution, a molecular dynamics simulation predicts a relative population of the axial:equatorial conformers of 0.70:0.30 in consonance with NMR experimental data. Overall, the quantum chemical calculations as well as the simulations give insight into substituent effects, the conformational dynamics of benzylpenicillin, the frequency of ring-puckering motions, and the correlation of side chain and ring-puckering motions.     

A single‐molecule reaction cascade: First‐principles molecular dynamics simulation

The success of mechanochemistry is continued with the targeted organic synthesis of functional nano‐scale devices. In the present theoretical study, first principles molecular dynamics simulations are performed for a recently synthesized three‐ring system that functions as a molecular reaction cascade. Mechanochemical and photochemical reaction conditions are investigated. The system was designed in a way that three bonds would break consecutively in mechanochemical and sonochemical setups. We succeeded to simulate the reaction mechanisms with first‐principles molecular dynamics simulations and discuss the stereochemistry.     

BLUF Hydrogen network dynamics and UV/Vis spectra: A combined molecular dynamics and quantum chemical study

Blue light sensing using flavin (BLUF) protein photoreceptor domains change their hydrogen bond network after photoexcitation. To explore this phenomenon, BLUF domains from R. sphaeroides were simulated using Amber99 molecular dynamics (MD). Five starting configurations were considered, to study different BLUF proteins (AppA/BlrB), Trp conformations (“W_in”/“W_out”), structure determination (X‐ray/NMR), and finally, His protonation states. We found dependencies of the hydrogen bonds on almost all parameters. Our data show an especially strong correlation of the Trp position and hydrogen bonds involving Gln63. The latter is in some contradiction to earlier results (Obanayama et al., Photochem. Photobiol. 2008, 84 10031010). Possible origins and implications are discussed. Our calculations support conjectures that Gln63 is more flexible with Trp104 in W_in position. Using snapshots from MD and time‐dependent density functional theory, UV/vis spectra for the chromophore were determined, which account for molecular motion of the protein under ambient conditions. In accord with experiment, it is found that the UV/vis spectra of BLUF bound flavin are red‐shifted and thermally broadened for all calculated π → π* transitions, relative to gas phase flavin at T = 0 K. However, differences in the spectra between the various BLUF configurations cannot be resolved with the present approach. © 2012 Wiley Periodicals, Inc.     

Interface of ZnS single sheet and substrates: a first‐principles study

Two‐dimensional materials have aroused great interests because of their unique properties not seen in the bulk counterparts. The interface of the ZnS single sheet and substrates are studied in this paper. Different from isolated ZnS single sheet, here in this study, the ZnS single sheet has a remarkable corrugation feature because of the interaction between the ZnS single sheet and the substrate. The top‐site Zn means an attractive reaction with the substrate while the top‐site S means repulsive. For ZnS single sheet/Si(111) interface, the symmetry of the interface does not decrease after geometry optimization because the two layers have a good lattice matching. For ZnS single sheet/Ag(111) interface, an unbalanced interaction (attractive or repulsive) between the ZnS single sheet (Zn or S atom) and Ag surface leads to remarkable corrugation of the ZnS single sheet and the symmetry of the interface decreases. Copyright © 2017 John Wiley & Sons, Ltd.     

Interactions between hydrogenated diamond surface and adsorbates with different concentration of NO2 molecules: A first‐principles study

To elucidate the effects of NO₂ and H₂O molecules on the surface conductivity of hydrogenated diamond film, models of various adsorbates containing different molecular ratio of NO₂ and H₂O on hydrogenated diamond (100) surfaces were constructed. The adsorption energies, equilibrium geometries of adsorbates, density of states, and atomic Mulliken populations were studied by using first‐principles method. The results showed that H₂O molecule in the adsorbate could weaken the interactions between the adsorbates and hydrogenated diamond surface significantly. Compared with H₂O molecule, NO₂ molecule relaxes more dramatically when adsorbed on hydrogenated diamond surface. In addition, density of states for C(100):H–2NO₂, C(100):H–NO₂, and C(100):H–NO₂ + H₂O systems are very similar to each other, which indicates an obvious peak at valence band maximum level for all the three samples. It can be attributed to mainly single occupied molecule orbital of NO₂ molecule and slightly C–H bond of C(100):H substrate. When the adsorbates contain one NO₂ and two H₂O molecules, the peak shifts slightly into valence band, but its intensity increases significantly. All the samples exhibit p‐type surface conductivity when adsorbed with pure NO₂ molecules, and the surface conductivity remains as H₂O molecules added into the NO₂ adsorbate layer. However, for oxygenated diamond surface, very week interactions generate between diamond surface and various adsorbates. All the oxygenated diamond (100) surfaces with various adsorbates containing different NO₂ and H₂O molecules on it exhibit an insulating property.     

Structural,electronic, and magnetic properties of 3d transition metal doped GaN nanosheet: A first‐principles study

We have performed the first‐principles calculations on the structural, electronic, and magnetic properties of 3d transition‐metal? (Cr, Mn, Fe, Co, and Ni) atoms doped 2D GaN nanosheet. The results show that 3d TM atom substituting one Ga leads to a structural reconstruction around the 3d TM impurity compared to the pristine GaN nanosheet. The doping of TM atom can induce magnetic moments, which are mainly located on the 3d TM atom and its nearest‐neighbor N atoms. It is found that Mn‐ and Ni‐doped GaN nanosheet with 100% spin polarization characters seem to be good candidates for spintronic applications. When two Ga atoms are substituted by two TM dopants, the ferromagnetic (FM) ordering becomes energetically more favorable for Cr‐, Mn‐, and Ni‐doped GaN nanosheet with different distances of two TM atoms. On the contrary, the antiferromagnetic (AFM) ordering is energetically more favorable for Fe‐doped GaN nanosheet. In addition, our GGA + U calculations show the similar results with GGA calculations. © 2016 Wiley Periodicals, Inc.     

Coarse‐grained molecular dynamics simulations of protein–ligand binding

Coarse‐grained molecular dynamics (CGMD) simulations with the MARTINI force field were performed to reproduce the protein–ligand binding processes. We chose two protein–ligand systems, the levansucrase–sugar (glucose or sucrose), and LinB–1,2‐dichloroethane systems, as target systems that differ in terms of the size and shape of the ligand‐binding pocket and the physicochemical properties of the pocket and the ligand. Spatial distributions of the Coarse‐grained (CG) ligand molecules revealed potential ligand‐binding sites on the protein surfaces other than the real ligand‐binding sites. The ligands bound most strongly to the real ligand‐binding sites. The binding and unbinding rate constants obtained from the CGMD simulation of the levansucrase–sucrose system were approximately 10 times greater than the experimental values; this is mainly due to faster diffusion of the CG ligand in the CG water model. We could obtain dissociation constants close to the experimental values for both systems. Analysis of the ligand fluxes demonstrated that the CG ligand molecules entered the ligand‐binding pockets through specific pathways. The ligands tended to move through grooves on the protein surface. Thus, the CGMD simulations produced reasonable results for the two different systems overall and are useful for studying the protein–ligand binding processes. © 2014 Wiley Periodicals, Inc.     

An automated approach for the parameterization of accurate intermolecular force‐fields: Pyridine as a case study

An automated protocol is proposed and validated, which integrates accurate quantum mechanical calculations with classical numerical simulations. Intermolecular force fields, (FF) suitable for molecular dynamics (MD) and Monte Carlo simulations, are parameterized through a novel iterative approach, fully based on quantum mechanical data, which has been automated and coded into the PICKY software, here presented. The whole procedure is tested and validated for pyridine, whose bulk phase, described through MD simulations performed with the specifically parameterized FF, is characterized by computing several of its thermodynamic, structural, and transport properties, comparing them with their experimental counterparts. © 2011 Wiley Periodicals, Inc.     

Electronic structure and reactivity of guanylthiourea: A quantum chemical study

Electronic structure analysis of guanylthiourea (GTU) and its isomers has been carried out using quantum chemical methods. Two major tautomeric classes (thione and thiol) have been identified on the potential energy (PE) surface. In both the cases conjugation of pi‐electrons and intramolecular H‐bonds have been found to play a stabilizing role. Various isomers of GTU on its PE surface have been analyzed in two different groups (thione and thiol). The interconversion from the most stable thione conformer ( GTU‐1 ) to the most stable thiol conformer ( GTU‐t1 ) was found to take place via bimolecular process which involves protonation at sulfur atom of GTU‐1 followed by subsequent C? N bond rotation and deprotonation. The detailed analysis of the protonation has been carried out in gas phase and aqueous phase (using CPMC model). Sulfur atom (S1) was found to be the preferred protonation site (over N4) in GTU‐1 in gas phase whereas N4 was found to be the preferred site of protonation in aqueous medium. The mechanism of S‐alkylation reaction in GTU has also been studied. The formation of alkylated analogs of thiol isomers (alkylated guanylthiourea) is believed to take place via bimolecular process which involves alkyl cation attack at S atom followed by C? N bond rotation and deprotonation. The reactive intermediate RS(NH₂)C? N? C(NH₂)₂⁺ belongs to the newly identified ^⊕N(←L)₂ class of species and provides the necessary dynamism for easy conversion of thione to thiol. © 2009 Wiley Periodicals, Inc. J Comput Chem, 2010     

Complexes of a Zn‐metalloenzyme binding site with hydroxamate‐containing ligands. A case for detailed benchmarkings of polarizable molecular mechanics/dynamics potentials when the experimental binding structure is unknown

Zn‐metalloproteins are a major class of targets for drug design. They constitute a demanding testing ground for polarizable molecular mechanics/dynamics aimed at extending the realm of quantum chemistry (QC) to very long‐duration molecular dynamics (MD). The reliability of such procedures needs to be demonstrated upon comparing the relative stabilities of competing candidate complexes of inhibitors with the recognition site stabilized in the course of MD. This could be necessary when no information is available regarding the experimental structure of the inhibitor–protein complex. Thus, this study bears on the phosphomannose isomerase (PMI) enzyme, considered as a potential therapeutic target for the treatment of several bacterial and parasitic diseases. We consider its complexes with 5‐phospho‐d ‐arabinonohydroxamate and three analog ligands differing by the number and location of their hydroxyl groups. We evaluate the energy accuracy expectable from a polarizable molecular mechanics procedure, SIBFA. This is done by comparisons with ab initio quantum‐chemistry (QC) calculations in the following cases: (a) the complexes of the four ligands in three distinct structures extracted from the entire PMI‐ligand energy‐minimized structures, and totaling up to 264 atoms; (b) the solvation energies of several energy‐minimized complexes of each ligand with a shell of 64 water molecules; (c) the conformational energy differences of each ligand in different conformations characterized in the course of energy‐minimizations; and (d) the continuum solvation energies of the ligands in different conformations. The agreements with the QC results appear convincing. On these bases, we discuss the prospects of applying the procedure to ligand‐macromolecule recognition problems. © 2016 Wiley Periodicals, Inc.     

Structural fluctuation and concerted motions in F1‐ATPase: A molecular dynamics study

F₁‐ATPase is an adenosine tri‐phosphate (ATP)‐driven rotary motor enzyme. We investigated the structural fluctuations and concerted motions of subunits in F₁‐ATPase using molecular dynamics (MD) simulations. An MD simulation for the α₃β₃γ complex was carried out for 30 ns. Although large fluctuations of the N‐terminal domain observed in simulations of the isolated β_E subunit were suppressed in the complex simulation, the magnitude of fluctuations in the C‐terminal domain was clearly different among the three β subunits (β_E, β_TP, and β_DP). Despite fairly similar conformations of the β_TP and β_DP subunits, the β_DP subunit exhibits smaller fluctuations in the C‐terminal domain than the β_TP subunit due to their dissimilar interface configurations. Compared with the β_TP subunit, the β_DP subunit stably interacts with both the adjacent α_DP and α_E subunits. This sandwiched configuration in the β_DP subunit leads to strongly correlated motions between the β_DP and adjacent α subunits. The β_DP subunit exhibits an extensive network of highly correlated motions with bound ATP and the γ subunit, as well as with the adjacent α subunits, suggesting that the structural changes occurring in the catalytically active β_DP subunit can effectively induce movements of the γ subunit. © 2010 Wiley Periodicals, Inc. J Comput Chem 2010