Small basis set density-functional theory methods corrected with atom- centered potentials

Density-functional theory (DFT) is currently the most popular method for modeling non-covalent interactions and thermochemistry. The accurate calculation of non-covalent interaction energies, reaction energies, and barrier heights requires choosing an appropriate functional and, typically, a relatively large basis set. Deficiencies of the density-functional approximation and the use of a limited basis set are the leading sources of error in the calculation of non-covalent and thermochemical properties in molecular systems. In this article, we present three new DFT methods based on the BLYP, M062X and CAM-B3LYP functionals in combination with the 6-31G* basis set and corrected with atom-centered potentials (ACPs). ACPs are one-electron potentials that have the same form as effective-core potentials, except they do not replace any electrons. The ACPs developed in this work are used to generate energy corrections to the underlying DFT/basis-set method such that the errors in predicted chemical properties are minimized while maintaining the low computational cost of the parent methods. ACPs were developed for the elements H, B, C, N, O, F, Si, P, S, and Cl. The ACP parameters were determined using an extensive training set of 118,655 data points, mostly of complete basis set coupled-cluster level quality. The target molecular properties for the ACP-corrected methods include non-covalent interaction energies, molecular conformational energies, reaction energies, barrier heights, and bond separation energies. The ACPs were tested first on the training set and then on a validation set of 42,567 additional data points. We show that the ACP-corrected methods can predict the target molecular properties with accuracy close to complete basis set wavefunction theory methods, but at a computational cost of double-ζ DFT methods. This makes the new BLYP/6-31G*-ACP, M062X/6-31G*-ACP, and CAM-B3LYP/6-31G*-ACP methods uniquely suited to the calculation of non-covalent, thermochemical, and kinetic properties in large molecular


Introduction
Over the past few decades, density functional theory (DFT) has become the leading approach for the quantum mechanical (QM) modeling of various molecular properties. 1 DFT's success can be attributed to its favorable balance between computational cost and accuracy, along with the existence of efficient algorithmic implementations widely available in modern software packages. One of the main sources of error in DFT is the choice of exchange-correlation functional approximation, which determines the accuracy of a given DFT method for a particular purpose. For instance, it is known that common density functional approximations are unable to accurately describe dispersion forces, which are critical when modeling non-covalent interactions and chemical reactions involving large molecules. [2][3][4] As a result, several works have focused on improving the accuracy of common density functional approximations by developing various dispersion correction techniques. [5][6][7][8] It has been shown that dispersion-corrected DFT methods in combination with large basis sets can predict various molecular properties and in particular non-covalent interaction energies of small and medium-size systems with accuracy similar or slightly lower than that obtained with nearly complete basis set wavefunction theory methods. [9][10][11][12][13][14] However, even with recent advances in computer technology, the applicability of dispersion-corrected DFT methods with a large basis set is a challenge for systems containing more than ca. one hundred atoms. This is unfortunate because there are many interesting problems involving systems in this molecular size range: supra-molecular and (bio)chemical complexes, nanostructured materials and surfaces, enzyme active sites, and many more. [15][16][17] The reason for this limitation is the unfavorable increase in computational cost as approximately the third power in the number of basis functions for common DFT methods. Consequently, the development of computationally inexpensive DFT based methods that allow efficient and accurate modeling of large systems is an important area of research. [18][19][20][21] The reduction in computational cost of dispersion-corrected DFT with the use of small double-ζ basis sets, such as 6-31G*, allows for the modeling of large molecular systems. The primary sources of error in dispersion-corrected DFT plus 6-31G* are the exchange-correlation functional approximation and basis set incompleteness error caused by the small size of the basis set. Methods have been recently proposed that show that these shortcomings can be mitigated efficiently. [22][23][24] For instance, Grimme and co-workers have shown that the D3 25,26 dispersion correction when combined with two additional semiempirical corrections 27,28 designed to mitigate basis set incompleteness error (collectively known as the "3c" approach 28 ), alleviate the deficiencies of the PBEh, HSE, and B3LYP functionals with a double-ζ basis set (yielding PBEh-3c 29 , HSE-3c 30 , and B3LYP-3c 31 methods). The 3c approach, which was initially proposed for a minimal basis set Hartree-Fock (HF) method (HF-3c 28 ), has also recently been extended to triple-ζ basis set DFT methods (B97-3c 32 and r 2 SCAN-3c 33 ). Several other methodologies based on DFT have been proposed in the literature [34][35][36][37][38][39][40] , and this underscores the interest in developing computationally inexpensive DFT based methods for modeling large molecular systems.
We have shown in earlier works that atom-centered potentials 41 (ACPs) offer a useful way to mitigate the underlying shortcomings of HF and DFT methods. [42][43][44][45][46][47][48][49][50][51][52][53][54][55] ACPs are similar to one-electron effective-core potentials 56,57 (ECPs) in functional form, except ACPs do not replace any electrons. Sharing a form similar to ECPs allows ACPs to be used in any software package that implements the use of ECPs. ACPs can be developed to yield energy corrections that minimize the errors in predicted properties for a target method and basis set combination by parametrization against high-level reference data. In this way, ACPs can be used to efficiently mitigate the shortcomings of double-ζ basis set DFT methods, since the use of ACPs incurs only a ca. 10% increase 58 in the computational cost relative to the uncorrected method.
In our recent work 58 , we developed four sets of ACPs for ten elements in combination with small and minimal basis set HF. These ACPs were trained against a set of 73,832 non-covalent properties (interaction energies, molecular conformational energies, and molecular deformation energies). Only noncovalent properties were used in the training set because small basis set HF, and HF in general, is limited to applications that do not involve bond breaking or formation, such as fast geometry optimizations, highthroughput conformer screening, and prediction of non-covalent interaction strengths in large systems.
Small basis set DFT is a much more promising approach for modeling thermochemical, kinetic, as well as non-covalent properties since, as previously mentioned, the parent DFT methods are already reasonably successful in the calculation of these quantities. Therefore, ACPs developed for double-ζ basis set DFT methods that mitigate the errors associated with the density functional approximation and basis set incompleteness error can be used to model all the aforementioned non-covalent properties as well as reaction energies, transition state searches, and in barrier height calculations of large molecular systems.
In this work, we developed ACPs for three density functionals combined with a double-ζ basis set (6-31G* 59,60 ) and Grimme's D3 25,26 dispersion correction scheme where applicable. The functionals chosen were BLYP 61,62 (generalized-gradient approximation or GGA functional), M062X 63 (hybrid meta-GGA functional), and CAM-B3LYP 64 (range-separated hybrid functional), mainly due to their popularity and performance for the target properties. 13,14,22,[65][66][67] The main target molecular properties are non-covalent interaction energies, molecular conformational energies, reaction energies, barrier heights, and bond separation energies. ACPs were developed for ten elements commonly encountered in organic chemistry and biochemistry (H, C, N, O, F, P, S, Cl) plus boron and silicon. The ACP development was carried out using a training set composed of 118,655 data points calculated at a high level of theory. We used a regularized linear least-squares fitting procedure (the LASSO [68][69][70] regression method) to obtain the parameters of the ACPs, which greatly simplifies the use of such a large training set. The strengths and weaknesses of the developed ACPs are evaluated and discussed based on their performance on the training set and a validation set consisting of additional 42,567 data points.

Computational Details
The ACP development procedure employed in this article has been described in detail in our earlier works [42][43][44]58 . We summarize it here for convenience. The mathematical form of an ACP is: where ( ) = ( ) − ( ), represents atom, is the distance, and | ⟩⟨ | are projection operators using real spherical harmonics based on atom with angular momentum quantum numbers and magnetic quantum numbers. The α ( ) and ( ) terms in Equation 1 are represented by Gaussian-type functions: where the coefficients ( ) and exponents ( ) are adjustable parameters determined via a regularized least-squares fit to reference data during the ACP development.
In order to find the exponents and coefficients that best mitigate the errors in the properties predicted using the target method and basis set combination, the ACP operator (Equation 1) is first added as a perturbative correction to the Hamiltonian. To first order in the ACP perturbation, the energy correction is: where the index runs over occupied molecular orbitals .
The Δ α ( α ) = ⟨ |(| ⟩ exp(− 2 ) ⟨ |)| ⟩ integrals (Equations 4 and 5) are the ACP energy terms and they are equal to the difference between the energy when an ACP with exponent is applied and the energy in absence of any ACP, divided by the ACP coefficient. In matrix notation, the and ( ) are vectors of ACP coefficients and ACP energy terms, respectively. Equations 1-5 indicate that the energy corrections introduced by the ACPs decay exponentially with interatomic distances and depend on the molecular wavefunction as well as the local chemical environment of a given atom.
The ACP development process starts with the compilation of a training set of target molecular properties. The training set must be chemically diverse and composed of systems exclusively with atoms for which the ACPs are being developed. Next, a set of exponents ( ) and angular momenta ( ) on each atom are selected, and the ACP energy terms (Δ α ( α )) are computed for each exponent and each entry in the training set. Once all the ACP energy terms are calculated, the optimal ACP coefficients and the associated exponents are determined using a regularized least-squares fitting subject to a constraint on the sum of the absolute values of the coefficients. The ACP development process ends with the generation of the ACPs for the target method and basis set combination.  14,193 , ISOL24 14,194 , IDISP 14,161,164,193,195,196 ,EIE22 197  total atomization energies a) defined as the difference between the energy of the complex and the sum of the monomer energies. A negative interaction energy indicates the complex is more stable than the separated monomers. b) the reference data was recalculated using DLPNO-CCSD(T)/CBS (see Reference 58). c) comprises non-covalently bound dimer complexes where at least one of the monomers is negatively charged. d) the reference data was calculated at the CCSD(T)/CBS level using the same extrapolation method as in Reference 124. e) defined as the difference between the energy of a particular conformer and a lower-energy conformer of the same molecule. f) comprises negatively charged conformers. g) defined as the difference between the sum of energies of reactants minus that of products. h) forward barrier height is the energy difference between transition state and reactant(s) or pre-reaction complex; reverse barrier height is the energy difference between transition state and product(s) or post-reaction complex. i) the reference data was recalculated using DLPNO-CCSD(T)/CBS with the same extrapolation method as in Reference 74. j) the barrier heights are relative to the pre-or post-reaction complexes for forward and reverse barriers, respectively. k) the barrier heights are relative to the isolated reactant(s) or product(s) for forward and reverse barriers, respectively. l) defined as the difference between the energy of a molecule and its radical fragments formed by cleavage of a particular bond. m) includes molecular deformation energies, isomerization energies, and total atomization energies. Molecular deformation energy is the difference between the energy of a molecule deformed along a particular normal mode and the energy of the same molecule at equilibrium. Isomerization energy is the energy difference between a molecule and one of its isomers. Total atomization energy is the energy difference between a molecule and the sum of the energies of all its constituent atoms. n) includes deformations along the imaginary mode of transition state structures from Reference 74. The reference data was calculated at DLPNO-CCSD(T)/CBS level using the same extrapolation method as used in Reference 58. o) contains reference data which is negative in magnitude due to deformation along the imaginary normal mode, indicating that the deformed molecule is more stable than the transition state structure. a) defined as the difference between the energy of the complex and the sum of energy of the monomers. A negative interaction energy indicates the complex is more stable than the separated monomers. b) comprises non-covalently bound complexes with at least one negatively charged monomer. c) defined as the difference between the energy of a particular conformer and a lower-energy conformer of the same molecule. d) only 30 systems used; we could not find the rest of the systems mentioned in the supporting information of Reference 222. e) contains mostly conformational energies but also some molecular deformation energies. f) only a subset of the actual data used. g) only macrocyclic peptides used. h) only a subset from the PEPCONF 140 database for which reference data was recalculated at the DLPNO-CCSD(T)/CBS level of theory (see Reference 58 for more details). i) comprises negatively charged conformers. j) defined as the difference between the sum of energies of reactants minus that of products. k) forward barrier height is the energy difference between transition state and reactant(s) or pre-reaction complex; reverse barrier height is the energy difference between transition state and product(s) or post-reaction complex. l) the barrier heights are relative to the pre-or post-reaction complexes for forward and reverse barriers, respectively. m) contains unpublished data generated in an ongoing project, which will be published elsewhere.

Results and Discussion
Results obtained using the target methodology (BLYP-D3/6-31G*, M062X/6-31G*, and CAM-B3LYP-D3/6-31G*) with and without the proposed ACPs compared against the reference data are shown in Figure 1  Error distribution (in kcal/mol) associated with the uncorrected and ACP-corrected double-ζ DFT methods. The various molecular properties represented are: "NCI" or non-covalent interaction energies, "CONF" or molecular conformational energies, "RE" or reaction energies, "BH" or barrier height energies, and "BSE" or bond separation energies. Suffixes "train" and "val" are short for training and validation, respectively (see Table 1

(i) Non-covalent interaction energies
Regarding non-covalent interaction energies ("NCI"), the application of ACPs to the double-ζ basis set DFT methods decreases the overall MAE by 46-65%, indicating a substantial improvement in the description of non-covalent interactions. This is reasonable given that ACPs have been shown to be adequate at mitigating the basis set incompleteness error which is known to greatly affect the calculation of non-covalent interaction energies. 43,44,72 The best methods based on the overall performance for non- However, the remaining relatively high errors are probably due to the fact that the 6-31G* basis set lacks diffuse basis functions, which are known to be required for modeling anionic systems. 73 Given the overall good performance of ACPs for non-covalent interactions in the training set, we now examine their performance for systems outside the training set. Figure 1 shows that ACPs successfully bring down the overall MAE of the double-ζ basis set DFT methods by about 38-48% for the NCI validation subset ("NCI-val"). A more detailed look into the results ( Figure S3 of SI) shows that application of ACPs leads to MAE reductions (%∆MAE) in the range of about 27-67% for the NCI validation subsets containing complexes featuring a mix of common interactions found in large molecular systems (BlindNCI, DES15K, and NENCI-2021 data sets with a total of 17,413 data points). This range of %∆MAE for the subsets used in the validation resembles the %∆MAE obtained for data sets in training set with mixed character interactions such as S22x5, S66x8, and S66a8. Regarding hydrogen bonding interactions, the CE20, WaterOrg, and H2O20Bind10 subsets of the validation set feature these types of interactions. The MAEs for these three subsets are improved significantly on the application of ACPs by about 62-92%, depending on the method. The large reduction in error observed for mixed character and hydrogen bonding interactions in the validation are probably a consequence of the fact that the training set contains more data points of these two kinds than any other interaction type. In any case, the similarity between the error reduction in the validation and training set suggests that the proposed ACPs are fairly robust regarding these interactions, i.e., they can be applied to similar systems outside the training set.
Regarding the other interaction types, an assessment of the ACP performance in the validation stage for π-π stacking and pnicogen bonding interactions could not be carried out due to the scarcity of highlevel reference data in the literature. This scarcity compelled us to include all available systems containing these interaction types in the training set instead of reserving them for validation. The ACPs were further validated on systems containing other interaction types such as halogen bonding, chalcogen bonding, and close contact repulsions (XB45, CHAL336, R160x6, and R739x5 data sets) that were not specifically part of the training set. The MAEs for the data sets representing halogen bonding were mostly deteriorated (by >60%) on the application of ACPs, except for BLYP-D3/6-31G*, where the MAE was improved by 56%.
Application of ACPs to subsets representing chalcogen bonding led to a decrease in the MAEs of uncorrected methods by 6-42%, most likely due to the presence of O and S containing complexes in the training set that were not purely chalcogen bonded. On the other hand, the MAEs of the uncorrected methods for the subsets representing close contact repulsions were initially low (0.50-0.97 kcal/mol), and the application of ACPs led to either increase (by 9-30%) or decrease (by 9-33%). All these findings suggest an under-representation of halogen bonding, chalcogen bonding, and close contact repulsions in the training set compared to other interaction types. Therefore, future ACP development work will require more such systems to be included in the training set in order to increase the diversity and robustness of the resulting ACPs.
The application of ACPs to complexes containing at least one monomer with negative charge in the validation set led to mostly a reduction in the MAEs of the uncorrected methods by 19-75%. Nevertheless, applying the ACPs could only bring down the MAEs to values that were greater than 2 kcal/mol. This indicates that the performance of ACPs for anionic interactions in the validation set convey the same message as in the training set: even though errors decrease on the application of ACPs, the 6-31G* basis set is inappropriate for modeling anionic systems.
Finally, we tested the proposed ACPs for their performance regarding non-covalent interaction energies in some more challenging complexes that are significantly different from those in the training set.
In particular, we used the C60dimer set of non-covalent interaction energies between C60 dimeric complexes and the L7, S12L, S30L, and Ni2021 sets containing interaction energies between relatively large supramolecular systems. These data sets provide a more stringent test for ACPs than the other validation subsets because of the large size of the systems involved (the absolute reference energies range from 25 kcal/mol to 416 kcal/mol) as well as the multiplicity of and cooperativity between the non-covalent interactions present in these systems. The application of ACPs led to an overall reduction in the MAEs of the underlying methods for these sets by 11-66%, with only a few exceptions. For example, ACPs for CAM-B3LYP-D3/6-31G* reduce the MAEs of 4 out of 5 data sets by about 30-51%, indicating once again the robustness of the corresponding ACPs. The application of ACPs with CAM-B3LYP-D3/6-31G* led to an increase in the MAE of the C60dimer set from 1.78 kcal/mol to 3.97 kcal/mol. Although the error increases, this result also demonstrates that in the case of failure due to the systems studied being very different from those on which the ACP were trained, the results from the ACP-corrected methods are far from being catastrophic.

(ii) Conformational energies
Another molecular property included in the training set is molecular conformational energies ("CONF") with 11,161 data points. The purpose of these data is to inform ACPs about how molecular motion along rotatable bonds and torsional angles involving various effects (π-conjugation, steric interactions, intramolecular hydrogen-bonding, and electron repulsion) influence the molecular potential energy surfaces. The CONF systems include peptides, nucleotides, carbohydrates, alcohols, hydrocarbons, (H2O)11 clusters, and other molecules representing pharmaceuticals, catalysts, synthetic precursors, industrial chemicals, and organic compounds. Interestingly, the overall MAE of all uncorrected double-ζ basis set DFT methods for CONF in the training set is below 1 kcal/mol, likely due to error cancellation. In order to understand the performance of ACPs for CONF in the training set, we take a closer look at the results of CAM-B3LYP-D3/6-31G*-ACP (see Figure S1 of SI). These ACPs perform well for most

(iii) Reaction energies
We turn our attention now to covalent properties involving bond breaking and formation. One of such properties included in the training set are chemical reaction energies ("RE"). The training set contains 8,315 RE data points. Application of ACPs to the RE subsets in the training set shows an overall reduction in the MAEs of the double-ζ basis set DFT methods by about 29-34%. This improvement in the MAE mainly reflects the MAE decrease observed in two specific RE subsets, MN-RE (7,555 data points) and BH9-RE (449 data points). These subsets contain data points representing a variety of chemical reactions.
For instance, the BH9-RE subset contains reaction energies of the nine most common elementary reactions encountered in organic and bioorganic chemistry. Besides MN-RE and BH9-RE, the ACPs also perform generally well (MAE reductions between 11-87% for M062X/6-31G*-ACP and CAM-B3LYP-D3/6-31G*-ACP) for other RE subsets of various types of reactions with fewer data points, with only a few exceptions. The decrease in the MAEs for reaction energies is also observed for the relatively large RE validation data set ("RE-val"), composed of a subset of the W4-17-RE set with 5,205 data points. For this data set, which was not used in the ACP training, the MAEs decrease by about 31-50% on application of ACPs compared to the uncorrected double-ζ basis set DFT methods. This indicates that the ACPs are successful and robust for reaction energies.

(iv) Barrier heights
Besides reaction energies, the training set also included barrier heights of chemical reactions ("BH") in order to make the eventual ACP-corrected methods usable for kinetic studies. The total number of BH data points contributes nearly 50% of the training set data (58,197) and is the most dominant property overall. The application of ACPs reduces the overall MAE of the parent methods for barrier heights by about 27-50%. This improvement originates mostly from the good performance of ACPs (MAE reductions by 13-52%) on the four main BH data sets: Grambow2020-B97D3 (32,722 data points), Grambow2020-ωB97XD3 (23,922 data points), BH9 (898), and E2SN2 (418). Note that the BH9 data set was designed recently 74 to be used for the particular purpose of developing the ACPs in this work and contains various model reactions that increase the diversity of the BH data in the training set.
For the other BH data sets in the training set besides BH9, E2SN2, and the two Grambow2020 subsets, the performance of ACPs is also quite good. ACPs mostly bring down the MAEs for the other subsets by about 22-94%. Contrary to most other BH data points in the training set, the INV24 and BHROT27 subsets feature barriers for processes that do not involve bond breaking or formation. The application of ACPs on these two data sets do not show a significant reduction in MAE compared to the uncorrected methods.
To validate the performance of ACPs on barrier height prediction, we applied the uncorrected methods and their ACP-corrected counterparts to a data set with 88 data points (WaterOrgBH) not included in the training set. This data set contains to-be-published barrier height reference data for pericyclic-type As an additional test, we also applied the proposed ACPs to the systems in the work of Bistoni et al. 75 , involving the Baeyer-Villiger reaction catalyzed by the cyclohexanone monooxygenase enzyme.
The geometries used for ACP testing are the relevant stationary points of the reaction (reactant, intermediate, transition state, and product), obtained independently with three active site models of increasing size. The reference data for the relative energies along the reaction profile were obtained with DLPNO-CCSD(T0)/def2-TZVPP. Point charges were used to model the electrostatic potential from the surrounding protein environment. The errors in the calculated reaction energies and barrier heights with and without ACPs for the three differently sized active sites are shown in Table 3. The errors of various methods relative to DLPNO-CCSD(T0)/def2-TZVPP data in Table 3

(v) Bond separation energies and other properties
Lastly, we analyze ACP performance for bond separation energies ("BSE"). In our training set, we Other chemical properties that were included in the ACP training set were: molecular isomerization energies ("ISOM"), total atomization energies ("TAE"), and molecular deformation energies ("DEF"). A few ISOM and TAE data sets from the literature were included in the training set, mainly to ensure that ACPs do not lead to a degradation in the functional performance for these properties. The results from Figure S1 in the SI show that the application of the ACPs reduce the MAEs of ISOM and TAE with CAM- In addition, a large set of DEF data points were also included in the training set. These data points represent energy differences between a molecule at its equilibrium geometry and the same molecule deformed along its various normal modes. Our intention when we included these data in the training set was to improve the description of molecular potential energy surfaces to obtain reasonably accurate geometries and energy derivatives or, at least, prevent the ACP fitting procedure from deteriorating the performance of the uncorrected method for these properties. In this connection, we performed transition state searches and geometry optimizations on the stationary point structures (reactants, products, and transition states) of 100 representative chemical reactions taken from the BH9 data set. Application of ACPs resulted in geometries with average root-mean-square-deviations (RMSDs) ranging between 0.110-0.131 Å (Table S6 of SI) compared to the reference geometries. These average RMSDs are only slightly higher than the uncorrected double-ζ basis set DFT methods (0.068-0.123 Å). Note that the geometries yielded by DFT-D3 methods with a double-ζ basis set have been shown in the literature to be reasonably close to that obtained with quadruple-ζ basis sets. 22

(vi) Overall ACP performance
In this subsection, we briefly discuss the overall performance of ACPs. First, we consider the various biases of the base method as measured by the mean signed errors (MSE) and standard deviations (SD) shown in Figure 1. Regarding non-covalent interaction energies, the Figure  Lastly, we compare the performance of uncorrected and ACP-corrected double-ζ basis set DFT methods with the same DFT methods using large basis sets. We present the MAEs of the various methods in Table 4. Only a few representative data sets for which nearly complete basis set DFT results have been reported in the literature are shown in the table. For non-covalent properties (interaction and conformational energies), the MAEs of ACP-corrected methods are slightly higher than the same methods at complete basis set limit. However, a particular exception is seen in the case of large hydrogen-bonded water cluster sets, Water38 and H2O20Bind10 with CAM-B3LYP-D3/6-31G*-ACP where the MAEs are lower than CAM-B3LYP-D3/def2-QZVPDD by almost a factor of two. In the case of covalent properties (reaction energies and barrier heights), the MAEs of M062X/6-31G*-ACP and CAM-B3LYP-D3/6-31G*-ACP methods are lower than the same functionals with def2-QZVPDD basis set for a few cases like BH9-RE, BSR36, and PX13. The MAEs with CAM-B3LYP-D3/6-31G*-ACP is also lower than CAM-B3LYP-D3/def2-QZVPDD for two other data sets, namely HTBH38 and WCPT27. Another interesting result is observed for the ACP-corrected BLYP-D3/6-31G* method, where the MAEs of data sets containing reaction energies and barrier heights, with only a few exceptions, are generally lower than BLYP-D3/def2-QZVPDD MAEs, indicating that ACPs perform well in mitigating errors other than basis set incompleteness like those arising from the approximations in the density functionals itself. The overall comparison shown in Table 4 demonstrates that ACP-corrected methods have a performance close to nearly complete basis set DFT or sometimes even better, but at a ca. one order of magnitude lower computational time. a) details about the data sets can be found in Table S1 of the Supporting Information, b) from Reference 13, c) represents the value obtained with BLYP-XDM/def2-QZVPP in Reference 74,d) represents the value obtained with M062X/def2-QZVPP in Reference 74, e) represents the value obtained with CAM-B3LYP-XDM/def2-QZVPP in Reference 74.

Conclusions
The use of very accurate quantum mechanical methods, such as nearly complete basis set wavefunction theory methods, is not practical for applications involving large systems because of the steep scaling of their computational cost with system size. A low-cost alternative is the use of small basis set density-functional theory (DFT) methods, but this requires that the inherent shortcomings in these methods be mitigated. The focus of this article is to improve the performance of small basis set DFT methods without sacrificing their computational efficiency.
We developed and applied atom-centered potentials (ACPs) to mitigate the shortcomings of BLYP-D3, M062X, and CAM-B3LYP-D3 methods, in combination with the 6-31G* basis set. We expect these shortcomings to be primarily the error from the choice of density functional approximation and basis set incompleteness from the limited size of the basis set.
The ACPs presented in this work were developed for ten elements (H, B, C, N, O, F, Si, P, S, Cl).
The parametrization of the ACPs was carried out using a regularized linear least-squares fitting procedure (the LASSO regression method) using a training set of 118,655 data points calculated mostly using wavefunction theory methods extrapolated to the complete basis set limit. The main molecular properties in the training set were non-covalent interaction energies, molecular conformational energies, reaction energies, barrier heights, and bond separation energies. The performance of the proposed ACPs was assessed using the training set and an additional validation set containing 42,567 data points not used during the ACP training.
Our assessment of the new ACP-corrected methods suggests that the ACPs reduce the mean absolute errors (MAEs) of double-ζ basis set DFT methods for most subsets, and in general, lead to an improved description of all molecular properties in the training set, indicating that ACPs successfully mitigate the deficiencies of the parent double-ζ basis set DFT methods. The best performing ACP-based method, i.e., CAM-B3LYP-D3/6-31G*-ACP, yields mean absolute errors, relative to high-level of theory, of around 0.7 kcal/mol for non-covalent interaction energies, 0.4 kcal/mol for conformational energies, 5.6 kcal/mol for reaction energies, 3.0 kcal/mol for barrier heights, and 2.5 kcal/mol for bond separation energies.
Further analysis of the performance of ACP-corrected methods on the validation set shows that ACPs are relatively robust, i.e., they are suitable for applications in systems outside the training set. However, there is a performance penalty when the system under consideration is significantly different from our training set, although no catastrophic results were obtained at any point. This observation suggests that increasing the size of the training set and incorporating more diversity is a straightforward way to improve the ACPs. In addition, using basis sets containing diffuse functions could be a way to overcome the limitations of the present methods in the description of negatively charged systems. The proposed new ACPs correct small basis set DFT methods and improve the accuracy of their parent methods, and allow carrying out quantum mechanical calculations of large systems at a reasonably low computational cost.

Supporting information
The supporting Information (SI) is available free of charge on the ACS publications website at DOI: Basis set and ACP files, detailed results, and Cartesian coordinates of all the data sets used in the study. (ZIP)