Molecular cluster building algorithm: electrostatic guidelines and molecular tailoring approach

Development and testing of an algorithm for efficient MP2/CCSD(T) energy estimation of molecular clusters with the 2-body approach

This work reports the development and testing of an automated algorithm for estimating the energies of weakly bound molecular clusters employing correlated theory. Firstly, the monomers and dimers of (homo/hetero) clusters are identified, and the sum of one-body and two-body contributions to correlation energy is calculated. The addition of this contribution to the Hartree-Fock full calculation (FC) energies provides a good estimate of the total energies at Møller-Plesset second-order perturbation theory (MP2)/coupled-cluster method with singles and doubles (CCSD) (T)-level theory using augmented Dunning basis sets. The estimated energies for several test clusters show an excellent agreement with their FC counterparts, with a substantial wall-clock time saving employing off-the-shelf hardware. Furthermore, the complete basis set (CBS) limit for MP2 energy computed using the two-body approach also agrees with its CBS energy with its FC counterpart.

Breaking covalent bonds in the context of the many-body expansion (MBE). I. The purported "first row anomaly" in XHn (X = C, Si, Ge, Sn; n = 1-4)

We present a new, novel implementation of the Many-Body Expansion (MBE) to account for the breaking of covalent bonds, thus extending the range of applications from its previous popular usage in the breaking of hydrogen bonds in clusters to molecules. A central concept of the new implementation is the in situ atomic electronic state of an atom in a molecule that casts the one-body term as the energy required to promote it to that state from its ground state. The rest of the terms correspond to the individual diatomic, triatomic, etc., fragments. Its application to the atomization energies of the XHn series, X = C, Si, Ge, Sn and n = 1-4, suggests that the (negative, stabilizing) 2-B is by far the largest term in the MBE with the higher order terms oscillating between positive and negative values and decreasing dramatically in size with increasing rank of the expansion. The analysis offers an alternative explanation for the purported "first row anomaly" in the incremental Hn-1X-H bond energies seen when these energies are evaluated with respect to the lowest energy among the states of the XHn molecules. Due to the "flipping" of the ground/first excited state between CH2 (3B1 ground state, 1A1 first excited state) and XH2, X = Si, Ge, Sn (1A1 ground state, 3B1 first excited state), the overall picture does not exhibit a "first row anomaly" when the incremental bond energies are evaluated with respect to the molecular states having the same in situ atomic states.

Unusually Large Hydrogen-Bond Cooperativity in Hydrogen Fluoride Clusters, (HF)n, n = 3 to 8, Revealed by the Molecular Tailoring Approach

In this work, our recently proposed molecular tailoring approach (MTA)-based method is employed for the evaluation of individual hydrogen-bond (HB) energies in linear (L) and cyclic (C) hydrogen fluoride clusters, (HF)_n (n = 3 to 8). The estimated individual HB energies calculated at the MP2(full)/aug-cc-pVTZ level for the L-(HF)_n are between 6.2 to 9.5 kcal/mol and those in the C-(HF)_n lie between 7.9 to 11.4 kcal/mol. The zero-point energy corrections and basis set superposition corrections are found to be very small (less than 0.6 and 1.2 kcal/mol, respectively). The cooperativity contribution toward individual HBs is seen to fall between 1.0 to 4.8 kcal/mol and 3.2 to 6.9 kcal/mol for linear and cyclic clusters, respectively. Interestingly, the HB energies in dimers, cleaved from these clusters, lie in a narrow range (4.4 to 5.2 kcal/mol) suggesting that the large HB strength in (HF)_n clusters is mainly due to the large cooperativity contribution, especially for n ≥ 5 (50 to 62% of the HBs energy). Furthermore, the HB energies in these clusters show a good qualitative correlation with geometrical parameters (H···F distance and F-H···F angles), stretching frequencies of F-H bonds, and electron density values at the (3, -1) bond critical points.

Appraisal of individual hydrogen bond strengths and cooperativity in ammonia clusters via a molecular tailoring approach

In this work, we propose and test a method, based on the molecular tailoring approach (MTA), for the evaluation of individual hydrogen bond (HB) energies in ammonia (NH₃)_n clusters. This methodology was tested, in our earlier work, on water clusters. Liquid ammonia being a universal, non-aqueous ionizing solvent, such information of individual HB strength is indispensable in many studies. The estimated HB energies by an MTA-based method, in (NH₃)_n for n = 3-8, were calculated to be in the range of 0.65 to 5.54 kcal mol^-1 with the cooperativity contribution falling between -0.54 and 1.88 kcal mol^-1 both calculated at the MP2(full)/aug-cc-pVTZ level of theory. It is seen that the strong HBs in (NH₃)_n clusters were additionally strengthened by the large contribution of HB cooperativity. The accuracy of these estimated HB energies was validated by approximately estimating the molecular energy of a given cluster by adding the sum of HB energies to the sum of monomer energies. This approximately estimated molecular energy of a given cluster was found to be in excellent agreement with the actual calculated values. The negligibly small difference (less than 5.6 kcal mol^-1) in these two values suggests that the estimated individual HB energies in ammonia clusters are quite reliable. Furthermore, these estimated HB energies by MTA are in excellent qualitative agreement with the other indirect measures of HB strength, such as HB bond distances and angles, N-H stretching frequency and the electron density values at the (3,-1) bond critical points.

Molecular Tailoring Approach for Estimating Individual Intermolecular Interaction Energies in Benzene Clusters

There is no general method available for the estimation of individual intermolecular interaction energies in weakly bound molecular clusters, and such studies are limited only to the dimer. Recently, we proposed a molecular tailoring approach-based method for the estimation of individual O-H···O hydrogen bond energies in water clusters. In the present work, we extend the applicability of this method for estimating the individual intermolecular interaction energies in benzene clusters, which are expected to be small. The basis set superposition error (BSSE)-corrected individual intermolecular interaction energies in linear (LN) benzene clusters, LN-(Bz)_n n = 3-7, were calculated to be in the range from -1.75 to -2.33 kcal/mol with the cooperativity contribution falling between 0.05 and 0.20 kcal/mol, calculated at the MP2.5/aug-cc-pVDZ level of theory. In the case of non-linear (NLN) benzene clusters, NLN-(Bz)_n n = 3-5, the BSSE-corrected individual intermolecular interaction energies exhibit a wider range from -1.16 to -2.55 kcal/mol with cooperativity contribution in the range from 0.02 to -0.61 kcal/mol. The accuracy of these estimated values was validated by adding the sum of interaction energies to the sum of monomer energies. These estimated molecular energies of clusters were compared with their actual calculated values. The small difference (<0.3 kcal/mol) in these two values suggests that our estimated individual intermolecular interaction energies in benzene clusters are quite reliable.

Electrostatic Potential Topology for Probing Molecular Structure, Bonding and Reactivity

Following the pioneering investigations of Bader on the topology of molecular electron density, the topology analysis of its sister field viz. molecular electrostatic potential (MESP) was taken up by the authors’ groups. Through these studies, MESP topology emerged as a powerful tool for exploring molecular bonding and reactivity patterns. The MESP topology features are mapped in terms of its critical points (CPs), such as bond critical points (BCPs), while the minima identify electron-rich locations, such as lone pairs and π-bonds. The gradient paths of MESP vividly bring out the atoms-in-molecule picture of neutral molecules and anions. The MESP-based characterization of a molecule in terms of electron-rich and -deficient regions provides a robust prediction about its interaction with other molecules. This leads to a clear picture of molecular aggregation, hydrogen bonding, lone pair–π interactions, π-conjugation, aromaticity and reaction mechanisms. This review summarizes the contributions of the authors’ groups over the last three decades and those of the other active groups towards understanding chemical bonding, molecular recognition, and reactivity through topology analysis of MESP.

Molecular Tailoring Approach for the Estimation of Intramolecular Hydrogen Bond Energy

Hydrogen bonds (HBs) play a crucial role in many physicochemical and biological processes. Theoretical methods can reliably estimate the intermolecular HB energies. However, the methods for the quantification of intramolecular HB (IHB) energy available in the literature are mostly empirical or indirect and limited only to evaluating the energy of a single HB. During the past decade, the authors have developed a direct procedure for the IHB energy estimation based on the molecular tailoring approach (MTA), a fragmentation method. This MTA-based method can yield a reliable estimate of individual IHB energy in a system containing multiple H-bonds. After explaining and illustrating the methodology of MTA, we present its use for the IHB energy estimation in molecules and clusters. We also discuss the use of this method by other researchers as a standard, state-of-the-art method for estimating IHB energy as well as those of other noncovalent interactions.

Direct and Reliable Method for Estimating the Hydrogen Bond Energies and Cooperativity in Water Clusters, Wn, n = 3 to 8

No direct method for estimating the individual O-H···O hydrogen bond (H-bond) energies in water clusters (W_n) exists in the literature. In this work, we propose such a direct method based on the molecular tailoring approach, which also enables the estimation of the cooperativity contributions. The calculated H-bond energies at MP2(full)/aug-cc-pVTZ and CCSD(T)/aug-cc-pVDZ levels for W_n, n = 3 to 8, agree well with one another and fall between 0.3 and 11.6 kcal mol^-1 with the cooperativity contributions in the range of -1.2 and 7.0 kcal mol^-1. For gauging the accuracy of our H-bond energies for a cluster, the H-bond energy sum is added to the sum of monomer energies, and the results are compared with the respective total energy. These two values agree with each other to within 8.3 mH (∼5 kcal mol^-1), testifying the accuracy of our estimated H-bond energies. Further, these H-bond strengths show a good correlation with the respective O-H stretching frequencies and the molecular electron density values at the (3, -1) O-H···O H-bond critical point.

Microhydration of a Carbonyl Group: How does the Molecular Electrostatic Potential (MESP) Impact the Formation of (H2O)n:(R2C═O)Complexes?

The presence of a carbonyl group in a molecule usually leads to the identification of a π-hole on the molecular electrostatic potential (MESP) of the species. How does this electrophilic site influence the formation of microhydrated complexes? To address this point, a panel of R₂CO solutes with various MESPs was selected, and we identified the structures and properties of several complexes containing one, two, three and six water molecules. The following solutes were considered in the present study: H₂CO, F₂CO, Cl₂CO,(NC)₂CO and H₂C═CO. Geometry optimizations and frequency calculations were carried out at the LC-ωPBE/6-311++G(d,p) level, with the GD3BJ empirical correction for dispersion. For a number of n water molecules around the R₂CO solute, the structure and the features of the most stable (H₂O)_n:(R₂CO) complexes are highly dependent on the MESP of the isolated R₂CO solute. The formation of pi-hole bondings appears to play a decisive role in the initiation of a three-dimensional organization of water molecules around the solute.

Analysis of Different Fragmentation Strategies on a Variety of Large Peptides: Implementation of a Low Level of Theory in Fragment-Based Methods Can Be a Crucial Factor

We have investigated the performance of two classes of fragmentation methods developed in our group (Molecules-in-Molecules (MIM) and Many-Overlapping-Body (MOB) expansion), to reproduce the unfragmented MP2 energies on a test set composed of 10 small to large biomolecules. They have also been assessed to recover the relative energies of different motifs of the acetyl(ala)18NH2 system. Performance of different bond-cutting environments and the use of Hartree-Fock and different density functionals (as a low level of theory) in conjunction with the fragmentation strategies have been analyzed. Our investigation shows that while a low level of theory (for recovering long-range interactions) may not be necessary for small peptides, it provides a very effective strategy to accurately reproduce the total and relative energies of larger peptides such as the different motifs of the acetyl(ala)18NH2 system. Employing M06-2X as the low level of theory, the calculated mean total energy deviation (maximum deviation) in the total MP2 energies for the 10 molecules in the test set at MIM(d=3.5Å), MIM(η=9), and MOB(d=5Å) are 1.16 (2.31), 0.72 (1.87), and 0.43 (2.02) kcal/mol, respectively. The excellent performance suggests that such fragment-based methods should be of general use for the computation of accurate energies of large biomolecular systems.

Low-energy structures of benzene clusters with a novel accurate potential surface

The benzene-benzene (Bz-Bz) interaction is present in several chemical systems and it is known to be crucial in understanding the specificity of important biological phenomena. In this work, we propose a novel Bz-Bz analytical potential energy surface which is fine-tuned on accurate ab initio calculations in order to improve its reliability. Once the Bz-Bz interaction is modeled, an analytical function for the energy of the Bzn clusters may be obtained by summing up over all pair potentials. We apply an evolutionary algorithm (EA) to discover the lowest-energy structures of Bzn clusters (for n=2-25), and the results are compared with previous global optimization studies where different potential functions were employed. Besides the global minimum, the EA also gives the structures of other low-lying isomers ranked by the corresponding energy. Additional ab initio calculations are carried out for the low-lying isomers of Bz3 and Bz4 clusters, and the global minimum is confirmed as the most stable structure for both sizes. Finally, a detailed analysis of the low-energy isomers of the n = 13 and 19 magic-number clusters is performed. The two lowest-energy Bz13 isomers show S6 and C3 symmetry, respectively, which is compatible with the experimental results available in the literature. The Bz19 structures reported here are all non-symmetric, showing two central Bz molecules surrounded by 12 nearest-neighbor monomers in the case of the five lowest-energy structures.

Hybrid QTAIM and electrostatic potential-based quantum topology phase diagrams for water clusters

The topological diversity of sets of isomers of water clusters (W = H₂O)_n, 7 ≤ n ≤ 10, is analyzed employing the scalar fields of total electronic charge density ρ(r) and the molecular electrostatic potential (MESP).

Molecular tailoring approach: a route for ab initio treatment of large clusters

Conspectus Chemistry on the scale of molecular clusters may be dramatically different from that in the macroscopic bulk. Greater understanding of chemistry in this size regime could greatly influence fields such as materials science and atmospheric and environmental chemistry. Recent advances in experimental techniques and computational resources have led to accurate investigations of the energies and spectral properties of weakly bonded molecular clusters. These have enabled researchers to learn how the physicochemical properties evolve from individual molecules to bulk materials and to understand the growth patterns of clusters. Experimental techniques such as infrared, microwave, and photoelectron spectroscopy are the most popular and powerful tools for probing molecular clusters. In general, these experimental techniques do not directly reveal the atomistic details of the clusters but provide data from which the structural details need to be unearthed. Furthermore, the resolution of the spectral properties of energetically close cluster conformers can be prohibitively difficult. Thus, these investigations of molecular aggregates require a combination of experiments and theory. On the theoretical front, researchers have been actively engaged in quantum chemical ab initio calculations as well as simulation-based studies for the last few decades. To obtain reliable results, there is a need to use correlated methods such as Møller-Plesset second order method, coupled cluster theory, or dispersion corrected density functional theory. However, due to nonlinear scaling of these methods, optimizing the geometry of large clusters still remains a formidable quantum chemistry challenge. Fragment-based methods, such as divide-and-conquer, molecular tailoring approach (MTA), fragment molecular orbitals, and generalized energy-based fragmentation approach, provide alternatives for overcoming the scaling problem for spatially extended molecular systems. Within MTA, a large system is broken down into two or more subsystems that can be readily treated computationally. Finally, the properties of the large system are obtained by patching the corresponding properties of all the subsystems. Due to these approximations, the resulting MTA-based energies carry some error in comparison with calculations based on the full system. An approach for correcting these errors has been attempted by grafting the error at a lower basis set onto a higher basis set. Furthermore, investigating the growth patterns and nucleation processes in clusters is necessary for understanding the structural transitions and the phenomena of magic numbers in cluster chemistry. Therefore, systematic building-up or the introduction of stochastics for generating molecular assemblies is the most crucial step for studying large clusters. In this Account, we discuss the working principle of MTA for probing molecular clusters at ab initio level followed by a brief summary of an automated and electrostatics-guided algorithm for building molecular assemblies. The molecular aggregates presented here as test cases are generated based on either an electrostatic criterion or the basin hopping method. At MP2 level computation, the errors in MTA-based grafted energies are typically reduced to a submillihartree level, reflecting the potential of finding accurate energies of molecular clusters much more quickly. In summary, MTA provides a platform for effectively studying large molecular clusters at ab initio level of theory using minimal computer hardware.

Molecular electrostatics for probing lone pair-π interactions

An electrostatics-based approach has been proposed for probing the weak interactions between lone pair containing molecules and π deficient molecular systems. For electron-rich molecules, the negative minima in molecular electrostatic potential (MESP) topography give the location of electron localization and the MESP value at the minimum (Vmin) quantifies the electron-rich character of that region. Interactive behavior of a lone pair bearing molecule with electron deficient π-systems, such as hexafluorobenzene, 1,3,5-trinitrobenzene, 2,4,6-trifluoro-1,3,5-triazine and 1,2,4,5-tetracyanobenzene explored within DFT brings out good correlation of the lone pair-π interaction energy (E(int)) with the Vmin value of the electron-rich system. Such interaction is found to be portrayed well with the Electrostatic Potential for Intermolecular Complexation (EPIC) model. On the basis of the precise location of MESP minimum, a prediction for the orientation of a lone pair bearing molecule with an electron deficient π-system is possible in the majority of the cases studied.

Appraisal of molecular tailoring approach for large clusters

High level ab initio investigations on molecular clusters are generally restricted to those of small size essentially due to the nonlinear scaling of corresponding computational cost. Molecular tailoring approach (MTA) is a fragmentation-based method, which offers an economical and efficient route for studying larger clusters. However, due to its approximate nature, the MTA-energies carry some errors vis-à-vis their full calculation counterparts. These errors in the MTA-energies are reduced by grafting the correction at a lower basis set (e.g., 6-31+G(d)) onto a higher basis set (e.g., aug-cc-pvdz or aug-cc-pvtz) calculation at MP2 level of theory. Further, better estimates of energies are obtained by making use of many-body interaction analysis. For this purpose, R-goodness (Rg) parameters for the three- and four-body interactions in a fragmentation scheme are proposed. The procedure employing grafting and many-body analysis has been tested out on molecular clusters of water, benzene, acetylene and carbon dioxide. It is found that for the fragmentation scheme having higher three- and four-body Rg-values, the errors in MTA-grafted energies are reduced typically to ~0.2 mH at MP2 level calculation. Coupled with the advantage in terms of computational resources and CPU time, the present method opens a possibility of accurate treatment of large molecular clusters.

Rapid topography mapping of scalar fields: large molecular clusters

An efficient and rapid algorithm for topography mapping of scalar fields, molecular electron density (MED) and molecular electrostatic potential (MESP) is presented. The highlight of the work is the use of fast function evaluation by Deformed-atoms-in-molecules (DAM) method. The DAM method provides very rapid as well as sufficiently accurate function and gradient evaluation. For mapping the topography of large systems, the molecular tailoring approach (MTA) is invoked. This new code is tested out for mapping the MED and MESP critical points (CP's) of small systems. It is further applied to large molecular clusters viz. (H(2)O)(25), (C(6)H(6))(8) and also to a unit cell of valine crystal at MP2/6-31+G(d) level of theory. The completeness of the topography is checked by extensive search as well as applying the Poincaré-Hopf relation. The results obtained show that the DAM method in combination with MTA provides a rapid and efficient route for mapping the topography of large molecular systems.

Topography of scalar fields: molecular clusters and π-conjugated systems

The pioneering works due to Bader and co-workers have generated widespread interest in the study of the topography of molecular scalar fields, the first step of which is the identification and characterization of the corresponding critical points (CPs). The topography of a molecular system becomes successively richer in going from the bare nuclear potential (BNP) to the molecular electrostatic potential (MESP) through the molecular electron density (MED). The present work clearly demonstrates, through the study of some π-conjugated test molecules as well as molecular clusters, that the CPs could be economically located by following this path within ab initio level theory. Further, the topography mapping of large molecules, especially at a higher level of theory, is known to be a demanding task. However, it is rendered possible by following the above sequential mapping assisted by a divide-and-conquer-type method termed as the molecular tailoring approach (MTA). This is demonstrated with the topography mapping of β-carotene and benzene nonamer at MP2 and a (H(2)O)(32) cluster at the HF level of theory, which are rather challenging problems with contemporary off-the-shelf computer hardware.