Abstract
Objectives
According to a recent expert review on Density Functional Theory (DFT) [[i]], one of the current challenges is "the need to improve the description of ... dispersion/van der Waals interactions". In recent years a number of approaches addressing this issue have appeared. One of promising methods (BH-DFT-D) have been developed by my host, Dr. Alisa Krishtal and coworkers [[ii]]. This approach is based on distributed atomic polarizabilities obtained from perturbed densities. First, molecular polarizabilities are calculated within the finite-field approach. Subsequently, atomic components of the polarizability tensors are obtained using atomic weighting functions defined according to Hirshfeld [[iii]] or iterative Hirshfeld partitioning [[iv]]. From the atomic polarizabilities, intermolecular interaction energies are obtained.

There are indications in the literature [[v]] that most common density functionals overestimate molecular polarizabilities and this deficiency can be improved by applying a density functional with a correct asymptotic behavior. This suggests that the atomic polarizabilities and, eventually, intermolecular interaction energies can be also improved if an asymptotically corrected functional is used. Therefore, the central idea of the current project is to implement and to test asymptotically corrected functionals within the BH-DFT-D method.

[i]. Cohen A. J., Mori-Sánchez P., Yang W., Chem. Rev. 2012, 112, 289–320

[ii]. Krishtal A., Vanommeslaeghe K., Olasz A., Veszprémi T., van Alsenoy C., Geerlings P., J. Chem. Phys. 2009, 130, 174101

[iii]. Hirshfeld F. L. Theoret. Chim. Acta, 1977, 44, 129--139.

[iv]. Bultinck P., van Alsenoy C., Ayers P. W., Carbó-Dorca R., J. Chem. Phys. 2007, 126, 144111.

[v]. Vasiliev I., Chelikowski J. R., Phys. Rev. A, 2010, 82, 012502.

Achievements
1. Implementation
We implemented two different asymptotically corrected functionals. The first one is the van Leeuwen–Baerends correlation functional (LB94)

[[i]]. However, it has a problematic performance in the high-density region. Therefore, Casida and Salahub proposed a method that combines a shifted LDA (or GGA) potential in the core region with the LB94 potential in the low-density region [[ii], [iii]].

The implementation of LB94 in the BRABO code was relatively straightforward. BRABO already has subroutines necessary to calculate the electron density and gradient. We wrote a subroutine implementing the LB94 functional in a close analogy to existing subroutines in BRABO implementing other exchange-correlation functionals. The subroutine was tested by comparison with MOLPRO results.

The implementation of the Casida–Salahub functional is a little more tricky, since it depends, in addition, on the ionization potential IP. In our implementation, the ionization potential was obtained from a Hartree–Fock calculation according to the Koopmans theorem. The values of IP and are input for a Casida–Salahub calculation. Therefore, a Casida–Salahub calculation consists essentially of three different ones: a Hartree–Fock run (to get the IP), an LDA run (to obtain ), and the proper Casida–Salahub run. We implemented the above equation in the BRABO program with an energy shift as an input. The preceding Hartree–Fock and LDA calculation are performed in separate BRABO runs, but the IP and were extracted automatically from the output by a shell script. The same script also creates an input for the final Casida–Salahub calculation.

In order to obtain all the necessary molecular polarizabilities within the finite-field approach, a large series of BRABO runs is necessary, since about 40 various perturbations (±x, ±y, ±z, ±xx, ±xy, ±xz, ..., ±xxx, ±xxy, ±xyz, ..., ±zzz) must be run in order to obtain the perturbed density matrices. These calculations are also governed by a script. We refrained from parallelization of BRABO, since about 40 independent calculations, which can be run simultaneously, are needed anyway.

The atomic polarizabilities are compute from the perturbed density matrices by the STOCK program. An MPI-parallelized version developed by Dr. Krishtal already existed. Upon a slight modification it was successfully compiled on the Nehalem computer.

2. Numerical calculations
For a number of reasons, calculations took much more time and effort, and I eventually managed to calculate only a part of what had been planned. The calculation were performed for a selection of monomers from the S22 set [[iv]] with the 6-311+G*. The BRABO jobs (perturbed-density evaluation) were carried out using Hartree–Fock, B3LYP, and Casida–Salahub functionals, while the final STOCK jobs (polarizability evaluation ) have been so far computed only for B3LYP.

By now, the polarizabilities were computed for the following molecules: C2H2 (several geometries), C2H4 (several geometries), C6H6 (several geometries), (CH3)2C(H)(C2H5), N(CH2)4. For other monomers from the S22 set, a part of calculations is finished.

[i]. van Leeuwen R., Baerends E. J., Phys. Rev. A, 1994, 49, 2421--2431

[ii]. Casida M. E., Salahub D. R., J. Chem. Phys., 2000, 113, 8918.

[iii]. Casida M. E., Casida K. C., Salahub D. R., Int. J. Quantum Chem., 1998, 70, 933

[iv]. Jurečka P., Šponer J., Černý J., Hobza P. Phys.Chem.Chem.Phys 2006, 8, 1985