Spectral representation of Hartree‐Fock exchange operator

We report a new mathematical result: it is possible to construct a spectral representation of the two particles Coulomb potential in the form of | r ? r ^′|^?1 = ∑_λλg( r ) g_λ( r ^′). We call this formula λ‐decomposition. Two special nontrivial cases of λ‐decomposition are reported together with the numerical analysis of the convergence for one of them. It is shown how λ‐decomposition allows to construct a new fast algorithm for Hartree‐Fock exchange operator calculation, in which the calculation of electron repulsion integrals (ERIs) is completely avoided. The connection between the new method and the resolution of identity and Cholesky decomposition based approaches has been established. Finally, the accuracy of ERIs evaluation within the new approach has been studied numerically. The results demonstrate that it is possible to achieve the accuracy of 10^?10 for the ERIs in wide range of their orbital exponents with relatively small number of terms in λ‐decomposition. © 2011 Wiley Periodicals, Inc. Int J Quantum Chem, 2012     

Interactive quantum chemistry: a divide-and-conquer ASED-MO method

We present interactive quantum chemistry simulation at the atom superposition and electron delocalization molecular orbital (ASED-MO) level of theory. Our method is based on the divide-and-conquer (D&C) approach, which we show is accurate and efficient for this non-self-consistent semiempirical theory. The method has a linear complexity in the number of atoms, scales well with the number of cores, and has a small prefactor. The time cost is completely controllable, as all steps are performed with direct algorithms, i.e., no iterative schemes are used. We discuss the errors induced by the D&C approach, first empirically on a few examples, and then via a theoretical study of two toy models that can be analytically solved for any number of atoms. Thanks to the precision and speed of the D&C approach, we are able to demonstrate interactive quantum chemistry simulations for systems up to a few hundred atoms on a current multicore desktop computer. When drawing and editing molecular systems, interactive simulations provide immediate, intuitive feedback on chemical structures. As the number of cores on personal computers increases, and larger and larger systems can be dealt with, we believe such interactive simulations-even at lower levels of theory-should thus prove most useful to effectively understand, design and prototype molecules, devices and materials.     

The treatment of solvation by a generalized Born model and a self-consistent charge-density functional theory-based tight-binding method

We present a model to calculate the free energies of solvation of small organic compounds as well as large biomolecules. This model is based on a generalized Born (GB) model and a self-consistent charge-density functional theory-based tight-binding (SCC-DFTB) method with the nonelectrostatic contributions to the free energy of solvation modeled in terms of solvent-accessible surface areas (SA). The parametrization of the SCC-DFTB/GBSA model has been based on 60 neutral and six ionic molecules composed of H, C, N, O, and S, and spanning a wide range of chemical groups. Effective atomic radii as parameters have been obtained through Monte Carlo Simulated Annealing optimization in the parameter space to minimize the differences between the calculated and experimental free energies of solvation. The standard error in the free energies of solvation calculated by the final model is 1.11 kcal mol(-1). We also calculated the free energies of solvation for these molecules using a conductor-like screening model (COSMO) in combination with different levels of theory (AM1, SCC-DFTB, and B3LYP/6-31G*) and compared the results with SCC-DFTB/GBSA. To assess the efficiency of our model for large biomolecules, we calculated the free energy of solvation for a HIV protease-inhibitor complex containing 3,204 atoms using the SCC-DFTB/GBSA and the SCC-DFTB/COSMO models, separately. The computed relative free energies of solvation are comparable, while the SCC-DFTB/GBSA model is three to four times more efficient, in terms of computational cost.     

Proton-coupled electron transfer of the phenoxyl/phenol couple: effect of Hartree-Fock exchange on transition structures

Proton-coupled electron transfer (PCET) and hydrogen atom transfer (HAT) reactions of the phenoxyl/phenol couple are studied theoretically by using wave function theory (WFT) as well as DFT methods. At the complete active space self-consistent field (CASSCF) level, geometry optimization is found to give two transition states (TSs); one is the PCET type with two benzene rings being nearly coplanar, and the other is the HAT type with two benzene rings taking a stacking structure. Geometry optimization at the (semilocal) DFT level, on the other hand, is found to give only one transition state (i.e., the PCET-type one) and fail to obtain the stacking TS structure. By comparing various levels of theories (including long-range corrected DFT functionals), we demonstrate that the Hartree-Fock exchange at long range plays a critical role in obtaining the sufficient stacking stabilization of the present open-shell system, and that the sole addition of empirical dispersion correction to semilocal DFT functionals may not be adequate for describing such a stacking interaction. Next, we investigate the solvent effect on the PCET and HAT TS thus obtained using the reference interaction site model self-consistent field (RISM-SCF) method. The results suggest that the free energy barrier increases with increasing polarity of the solvent, and that the solvent effects are stronger for the PCET TS than the stacking HAT TS pathway. The reason for this is discussed based on the dipole moment of different TS structures in solution.     

Computerized implementation of higher‐order electron‐correlation methods and their linear‐scaling divide‐and‐conquer extensions

We have implemented a linear‐scaling divide‐and‐conquer (DC)‐based higher‐order coupled‐cluster (CC) and Møller–Plesset perturbation theories (MPPT) as well as their combinations automatically by means of the tensor contraction engine, which is a computerized symbolic algebra system. The DC‐based energy expressions of the standard CC and MPPT methods and the CC methods augmented with a perturbation correction were proposed for up to high excitation orders [e.g., CCSDTQ, MP4, and CCSD(2)_TQ]. The numerical assessment for hydrogen halide chains, polyene chains, and first coordination sphere (C1) model of photoactive yellow protein has revealed that the DC‐based correlation methods provide reliable correlation energies with significantly less computational cost than that of the conventional implementations. © 2017 Wiley Periodicals, Inc.     

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.     

A semiempirical long-range corrected exchange correlation functional including a short-range Gaussian attenuation (LCgau-B97)

We applied an improved long‐range correction scheme including a short‐range Gaussian attenuation (LCgau) to the Becke97 (B97) exchange correlation functional. In the optimization of LCgau‐B97 functional, the linear parameters are determined by least squares fitting. Optimizing μ parameter (0.2) that controls long‐range portion of Hartree‐Fock (HF) exchange to excitation energies of large molecules (Chai and Head‐Gordon, J Chem Phys 2008, 128, 084106) and additional short‐range Gaussian parameters (a = 0.15 and k = 0.9) that controls HF exchange inclusion ranging from short‐range to mid‐range (0.5–3 Å) to ground state properties achieved high performances of LCgau‐B97 simultaneously on both ground state and excited state properties, which is better than other tested semiempirical density functional theory (DFT) functionals, such as ωB97, ωB97X, BMK, and M0x‐family. We also found that while a small μ value (～0.2) in LC‐DFT is appropriate to the local excitation and intramolecular charge‐transfer excitation energies, a larger μ value (0.42) is desirable in the Rydberg excitation‐energy calculations. © 2011 Wiley Periodicals, Inc. J Comput Chem, 2011     

Seminumerical calculation of the Hartree-Fock exchange matrix: application to two-component procedures and efficient evaluation of local hybrid density functionals

A two-component extension of the seminumerical procedure for the calculation of the Hartree-Fock (HF) exchange matrix recently presented by Neese et al. (Chem Phys 2009, 356, 98) was implemented into the program system TURBOMOLE. It is demonstrated that this allows for efficient self-consistent treatment of spin-orbit coupling at HF and hybrid density functional theory level. One-component HF calculations were performed to study the accuracy of integration grids and the exploitation of the molecular point group symmetry. The efficiency was tested, and for one-component calculations compared to the implementation realized by Neese. It was further demonstrated that local hybrid density functionals can be evaluated with this technique. The "prototype" of this class of functionals, Lh-BLYP, was applied to an organic molecule with more than 150 atoms.     

DCMB that combines divide‐and‐conquer and mixed‐basis set methods for accurate geometry optimizations,total energies,and vibrational frequencies of large molecules

We present a method, named DCMB, for the calculations of large molecules. It is a combination of a parallel divide‐and‐conquer (DC) method and a mixed‐basis (MB) set scheme. In this approach, atomic forces, total energy and vibrational frequencies are obtained from a series of MB calculations, which are derived from the target system utilizing the DC concept. Unlike the fragmentation based methods, all DCMB calculations are performed over the whole target system and no artificial caps are introduced so that it is particularly useful for charged and/or delocalized systems. By comparing the DCMB results with those from the conventional method, we demonstrate that DCMB is capable of providing accurate prediction of molecular geometries, total energies, and vibrational frequencies of molecules of general interest. We also demonstrate that the high efficiency of the parallel DCMB code holds the promise for a routine geometry optimization of large complex systems. © 2012 Wiley Periodicals, Inc.     

Parallel implementation of efficient charge–charge interaction evaluation scheme in periodic divide‐and‐conquer density‐functional tight‐binding calculations

A low‐computational‐cost algorithm and its parallel implementation for periodic divide‐and‐conquer density‐functional tight‐binding (DC‐DFTB) calculations are presented. The developed algorithm enables rapid computation of the interaction between atomic partial charges, which is the bottleneck for applications to large systems, by means of multipole‐ and interpolation‐based approaches for long‐ and short‐range contributions. The numerical errors of energy and forces with respect to the conventional Ewald‐based technique can be under the control of the multipole expansion order, level of unit cell replication, and interpolation grid size. The parallel performance of four different evaluation schemes combining previous approaches and the proposed one are assessed using test calculations of a cubic water box on the K computer. The largest benchmark system consisted of 3,295,500 atoms. DC‐DFTB energy and forces for this system were obtained in only a few minutes when the proposed algorithm was activated and parallelized over 16,000 nodes in the K computer. The high performance using a single node workstation was also confirmed. In addition to liquid water systems, the feasibility of the present method was examined by testing solid systems such as diamond form of carbon, face‐centered cubic form of copper, and rock salt form of sodium chloride. © 2017 Wiley Periodicals, Inc.     

The Fock Matrix Analysis for Atomic Orbitals in Molecular Orbitals II. The Electronic Structure of N2 Molecules

Weinhold's natural hybrid orbitals can be chosen as the molecular adapted atomic orbitals to build the canonical molecular orbitals of N₂ molecules. The molecular Fock matrix expanded in the natural hybrid orbitals can reveal deeper insight of the electronic structure and reaction of the N₂ molecule. For example, the magnitude of F_ab can signify the bonding character of the paired electrons as well as the diradical character of the unpaired electrons for both σ‐ and π‐types. Discarding the concept of the overlap between non‐orthogonal atomic orbitals, the different orbitals for different spins in the unrestricted Hartree‐Fock wavefunction reveal that there are three pairs of opposite spin density flows between two atoms, which proceed until the bonding molecular orbitals form.     

Truncation of small matrix elements based on the Euclidean norm for blocked data structures

Methods for the removal of small symmetric matrix elements based on the Euclidean norm of the error matrix are presented in this article. In large scale Hartree-Fock and Kohn-Sham calculations it is important to be able to enforce matrix sparsity while keeping errors under control. Truncation based on some unitary-invariant norm allows for control of errors in the occupied subspace as described in (Rubensson et al. J Math Phys 49, 032103). The Euclidean norm is unitary-invariant and does not grow intrinsically with system size and is thus suitable for error control in large scale calculations. The presented truncation schemes repetitively use the Lanczos method to compute the Euclidean norms of the error matrix candidates. Ritz value convergence patterns are utilized to reduce the total number of Lanczos iterations.     

WebMTA: a web-interface for ab initio geometry optimization of large molecules using molecular tailoring approach

A web-interface for geometry optimization of large molecules using a linear scaling method, i.e., cardinality guided molecular tailoring approach (CG-MTA), is presented. CG-MTA is a cut-and-stitch, fragmentation-based method developed in our laboratory, for linear scaling of conventional ab initio techniques. This interface provides limited access to CG-MTA-enabled GAMESS. It can be used to obtain fragmentation schemes for a given spatially extended molecule depending on the maximum allowed fragment size and minimum cut radius values provided by the user. Currently, we support submission of single point or geometry optimization jobs at Hartree-Fock and density functional theory levels of theory for systems containing between 80 to 200 first row atoms and comprising up to 1000 basis functions. The graphical user interface is built using HTML and Python at the back end. The back end farms out the jobs on an in-house Linux-based cluster running on Pentium-4 Class or higher machines using an @Home-based parallelization scheme (http://chem.unipune.ernet.in/ approximately tcg/mtaweb/).     

Two-level hierarchical parallelization of second-order Møller-Plesset perturbation calculations in divide-and-conquer method

A two‐level hierarchical parallelization scheme including the second‐order Møller–Plesset perturbation (MP2) theory in the divide‐and‐conquer method is presented. The scheme is a combination of coarse‐grain parallelization assigning each subsystem to a group of processors, with fine‐grain parallelization, where the computational tasks for evaluating MP2 correlation energy of the assigned subsystem are distributed among processors in the group. Test calculations demonstrate that the present scheme shows high parallel efficiency and makes MP2 calculations practical for very large molecules. © 2011 Wiley Periodicals, Inc. J Comput Chem, 2011     

Novel method of self‐interaction corrections in density functional calculations

It is demonstrated that the commonly applied self‐interaction correction (SIC) used in density functional theory does not remove all self‐interaction. We present as an alternative a novel method that, by construction, is totally free from self‐interaction. The method has the correct asymptotic 1/r dependence. We apply the new theory to localized f electrons in praseodymium and compare with the old version of SIC, the local density approximation (LDA) and with an atomic Hartree–Fock calculation. The results show a lowering of the f level, a contraction of the f electron cloud and a lowering of the total energy by 13 eV per 4 f electron compared to LDA. The equilibrium volume of the new SIC method is close to the ones given by LDA and the older SIC method and is in good agreement with experiment. The experimental cohesive energy is in better agreement using the new SIC method, both compared to LDA and another SIC method. © 2001 John Wiley & Sons, Inc. Int J Quant Chem 81: 247–252, 2001     

Binding affinity prediction with different force fields: examination of the linear interaction energy method

A systematic study of the linear interaction energy (LIE) method and the possible dependence of its parameterization on the force field and system (receptor binding site) is reported. We have calculated the binding free energy for nine different ligands in complex with P450cam using three different force fields (Amber95, Gromos87, and OPLS-AA). The results from these LIE calculations using our earlier parameterization give relative free energies of binding that agree remarkably well with the experimental data. However, the absolute energies are too positive for all three force fields, and it is clear that an additional constant term (gamma) is required in this case. Out of five examined LIE models, the same one emerges as the best for all three force fields, and this, in fact, corresponds to our earlier one apart from the addition of the constant gamma, which is almost identical for the three force fields. Thus, the present free energy calculations clearly indicate that the coefficients of the LIE method are independent of the force field used. Their relation to solvation free energies is also demonstrated. The only free parameter of the best model is gamma, which is found to depend on the hydrophobicity of the binding site. We also attempt to quantify the binding site hydrophobicity of four different proteins which shows that the ordering of gamma's for these sites reflects the fraction of hydrophobic surface area.     

A new method for partition of interaction energy. Relation between stabilization energy and orbital mixing

A new method for partition of interaction energy is proposed. The scheme given here easily connects the calculated stabilization energy with the orbital mixing in analyzing orbital interactions of molecules. The method can reveal the relation between the change of electron distribution and stabilization energy. As an example, orbital interaction energies are estimated for diazocompounds, diazomethane and diazirine.     

A survey of the parallel performance and accuracy of Poisson solvers for electronic structure calculations

We present an analysis of different methods to calculate the classical electrostatic Hartree potential created by charge distributions. Our goal is to provide the reader with an estimation on the performance—in terms of both numerical complexity and accuracy—of popular Poisson solvers, and to give an intuitive idea on the way these solvers operate. Highly parallelizable routines have been implemented in a first‐principle simulation code (Octopus ) to be used in our tests, so that reliable conclusions about the capability of methods to tackle large systems in cluster computing can be obtained from our work. © 2013 Wiley Periodicals, Inc.     

Substituent effects on heats of formation, group interactions, and detonation properties of polyazidocubanes

We have calculated the heats of formation (HOFs) for a series of polyazidocubanes by using the density functional theory (DFT), Hartree-Fock, and MP2 methods with 6-31G* basis set as well as semiempirical methods. The cubane skeleton was chosen for a reference compound, that is, the cubane skeleton was not broken in the process of designing isodesmic reactions. There exists group additivity for the HOF with respect to the azido group. The semiempirical AM1 method also produced reliable results for the HOFs of the title compounds, but the semiempirical MINDO3 did not. The relationship between HOFs and molecular structures was discussed. It was found that the HOF increases 330-360 kJ/mol for each additional number of the azido group being added to the cubane skeleton. The distance between azido groups slightly influences the values of HOFs. The interacting energies of neighbor azido groups in polyazidocubanes are in the range of 2.3 approximately 6.6 kJ/mol, which are so small and less related to the substituent numbers. The average interaction energy between nearest neighbor --N3 groups in the most stable conformer of octaazidocubane is 2.29 kJ/mol at the B3LYP/6-31G* level. The relative stability related to the number of azido groups of the title compounds was assessed based on the calculated HOFs, the energy gaps between the frontier orbitals, and the bond orders of the C--N3 and C--C bonds. The predicted detonation velocity of hepta- and octa-derivatives is over 9 km/s, and the detonation pressure of them is ca. 40 GPa or over.