Variational fitting methods for electronic structure calculations

CuGBasis: High-performance CUDA/Python library for efficient computation of quantum chemistry density-based descriptors for larger systems

CuGBasis is a free and open-source CUDA®/Python library for efficient computation of scalar, vector, and matrix quantities crucial for the post-processing of electronic structure calculations. CuGBasis integrates high-performance Graphical Processing Unit (GPU) computing with the ease and flexibility of Python programming, making it compatible with a vast ecosystem of libraries. We showcase its utility as a Python library and demonstrate its seamless interoperability with existing Python software to gain chemical insight from quantum chemistry calculations. Leveraging GPU-accelerated code, cuGBasis exhibits remarkable performance, making it highly applicable to larger systems or large databases. Our benchmarks reveal a 100-fold performance gain compared to alternative software packages, including serial/multi-threaded Central Processing Unit and GPU implementations. This paper outlines various features and computational strategies that lead to cuGBasis's enhanced performance, guiding developers of GPU-accelerated code.

Exchange-correlation kernel for perturbation dependent auxiliary functions in auxiliary density perturbation theory

CONTEXT

Analytic exchange-correlation kernel formulations are of the outermost importance for density functional theory (DFT) perturbation calculations. In this paper, the working equation for the exchange-correlation kernel of the generalized gradient approximation (GGA) for perturbation dependent auxiliary functions is derived and discussed in the framework of auxiliary density functional theory (ADFT). The presented new formulation is extended to the unrestricted approach, too. A comprehensive discussion of the implementation of the GGA ADFT kernel, using either the native exchange-correlation functional implementations in deMon2k or the ones from the LibXC library, is given. Calculations with analytic exchange-correlation kernels are compared to their finite difference counterparts. The obtained results are in quantitative agreement. Nevertheless, analytic GGA ADFT kernel implementations show substantial improvement in the computational performance. Similar results are reported for analytic second derivatives of effective core potential (ECP) and model core potential (MCP) matrix elements when compared to their finite difference counterparts in molecular frequency analyses.

METHOD

All calculations are performed in the framework of ADFT as implemented in deMon2k. In the ADFT analytic frequency calculations, auxiliary density perturbation theory was used. The underlying two-center exchange-correlation kernel matrix elements are calculated by numerical integration either with analytic or finite difference kernel expressions. Validation calculations are performed with the VWN and PBE functionals employing DFT-optimized DZVP basis sets in conjunction with automatically generated GEN-A2 auxiliary density function sets. In the (Pt3Cu)n cluster benchmark calculations, the RPBE functional was used. For Pt atoms, the quasi-relativistic LANL2DZ effective core potential with the corresponding valence basis set was employed, whereas for Cu atoms, the all-electron DFT-optimized TZVP basis was applied. The auxiliary density was expanded by the automatically generated GEN-A2* auxiliary function set. We run all benchmark calculations in parallel on 24 cores.

Use of Gaussian-Type Functions for Describing Fast Ion-Matter Irradiation with Time-Dependent Density Functional Theory

The electronic stopping power is an observable property that quantifies the ability of swift ions to penetrate matter to transfer energy to the electron cloud. The recent literature has proven the value of Real-Time Time-Dependent Density Functional Theory to accurately evaluate this property from first-principles, but questions remain regarding the capability of computer codes relying on atom-centered basis functions to capture the physics at play. In this Perspective, we draw attention to the fact that irradiation by swift ions triggers electron emission into the continuum, especially at the Bragg peak. We investigate the ability of Gaussian atomic orbitals (AOC), which were fitted to mimic continuum wave functions, to improve electronic stopping power predictions. AOC are added to standard correlation-consistent basis sets or STO minimal basis sets. Our benchmarks for water irradiation by fast protons clearly advocate for the use of AOC, especially near the Bragg peak. We show that AOC only need to be placed on the molecules struck by the ion. The number of AOC that are added to the usual basis set is relatively small compared to the total number of atomic orbitals, making the use of such a basis set an excellent choice from a computational cost point of view. The optimum basis set combination is applied for the calculation of the stopping power of a proton in water with encouraging agreement with experimental data.

Automatic Generation of Auxiliary Basis Sets in Spherical Representation Using the Cholesky Decomposition

Density fitting techniques that use automatically generated auxiliary basis sets generally rely on the formation of basis function products. Recently, Lehtola [ J. Chem. Theory Comput. 2021, 17, 6886-6900] presented a procedure making use of a purely spherical representation by adding auxiliary basis functions coupled to the required angular momentum quantum numbers for the product of spherical harmonics and then removing linear dependencies by means of a Cholesky decomposition. In this work, we extend this idea by making use of the explicit equations for the product of two spherical harmonics in the angular part of the basis function product. Some of the resulting terms are not directly accessible when popular standard integral libraries are used, which could prevent the widespread use of the exact product form. For these terms, we introduce four approximations of increasing sophistication that require integrals involving only standard Gaussian-type orbitals and thus can be computed with standard libraries. We assess the accuracy of the different schemes in the context of the aCD for the reconstruction of the electron repulsion integral matrix and absolute and relative single point energies and in the framework of optimally tuned range-separated hybrid functionals.

An information-theoretic approach to basis-set fitting of electron densities and other non-negative functions

The numerical ill-conditioning associated with approximating an electron density with a convex sum of Gaussian or Slater-type functions is overcome by using the (extended) Kullback-Leibler divergence to measure the deviation between the target and approximate density. The optimized densities are non-negative and normalized, and they are accurate enough to be used in applications related to molecular similarity, the topology of the electron density, and numerical molecular integration. This robust, efficient, and general approach can be used to fit any non-negative normalized functions (e.g., the kinetic energy density and molecular electron density) to a convex sum of non-negative basis functions. We present a fixed-point iteration method for optimizing the Kullback-Leibler divergence and compare it to conventional gradient-based optimization methods. These algorithms are released through the free and open-source BFit package, which also includes a L2-norm squared optimization routine applicable to any square-integrable scalar function.

Automatic Generation of Accurate and Cost-Efficient Auxiliary Basis Sets

We have recently discussed an algorithm to automatically generate auxiliary basis sets (ABSs) of the standard form for density fitting (DF) or resolution-of-the-identity (RI) calculations in a given atomic orbital basis set (OBS) of any form, such as Gaussian-type orbitals, Slater-type orbitals, or numerical atomic orbitals [J. Chem. Theory Comput. 2021, 17, 6886]. In this work, we study two ways to reduce the cost of such automatically generated ABSs without sacrificing their accuracy. We contract the ABS with a singular value decomposition proposed by Kállay [J. Chem. Phys. 2014, 141, 244113], used here in a somewhat different setting. We also drop the high-angular momentum functions from the ABS, as they are unnecessary for global fitting methods. Studying the effect of these two types of truncations on Hartree-Fock (HF) and second-order Møller-Plesset perturbation theory (MP2) calculations on a chemically diverse set of first- and second-row molecules within the RI/DF approach, we show that accurate total and atomization energies can be achieved by a combination of the two approaches with significant reductions in the size of the ABS. While the original approach yields ABSs whose number of functions NbfABS scales with the number of functions in the OBS, NOBSbf, as NABSbf = γNOBSbf with the prefactor γ ≈ O ( 10 ) , the reduction schemes of this work afford results of essentially the same quality as the original unpruned and uncontracted ABS with γ ≈ 5-6, while an accuracy that may suffice for routine applications is achievable with a further reduced ABS with γ ≈ 3-4. The observed errors are similar at HF and MP2 levels of theory, suggesting that the generated ABSs are highly transferable and can also be applied to model challenging properties with high-level methods.

Irradiation of Plutonium Tributyl Phosphate Complexes by Ionizing Alpha Particles: A Computational Study

The PUREX solvent extraction process, widely used for recovering uranium and plutonium from spent nuclear fuel, utilizes an organic solvent composed of tributyl phosphate (TBP). The emission of ionizing particles such as alpha particles, resulting from the decay of plutonium, makes the organic solvent vulnerable to degradation. Here, we study the ultrashort time alpha irradiation of tributylphosphate (TBP) and Pu(NO3)4(TBP)2 complex formed in the PUREX process. Electron dynamics is propagated by Real-Time-Dependent Auxiliary Density Functional Theory (RT-TD-ADFT). We investigate the use of previously proposed absorption boundary conditions (ABC) in the molecular orbital space to treat secondary electron emission. Basis set and exchange correlation functional effects with ABC are reported as well as a detailed analysis of the ABC parametrization. Preliminary results on the water molecule and then on TBP show that the phenomenological nature of the ABC parameters necessitates selecting appropriate values for each system under study. Irradiation of free and complexed TBP shows an influence of the ligands on the variation of atomic charges on the femtosecond time scale. An accumulation of atomic charges in the alkyl chains of TBP is observed in the case where the nitrate groups are predominantly irradiated. In addition, we find that the Pu atom regains its electric charge very rapidly after being hit by the projectile, with the coordination sphere serving as an electron reservoir to preserve its formal redox state. This study paves the road toward a full understanding of the degradation of organic extracants employed in the nuclear industry.

Reproducibility of Hybrid Density Functional Calculations for Equation-of-State Properties and Band Gaps

Hybrid density functional (HDF) approximations usually deliver higher accuracy than local and semilocal approximations to the exchange-correlation functional, but this comes with drastically increased computational cost. Practical implementations of HDFs inevitably involve numerical approximations─even more so than their local and semilocal counterparts due to the additional numerical complexity arising from treating the exact-exchange component. This raises the question regarding the reproducibility of the HDF results yielded by different implementations. In this work, we benchmark the numerical precision of four independent implementations of the popular Heyd-Scuseria-Ernzerhof (HSE) range-separated HDF on describing key materials' properties, including both properties derived from equations of state (EOS) and band gaps of 20 crystalline solids. We find that the energy band gaps obtained by the four codes agree with each other rather satisfactorily. However, for lattice constants and bulk moduli, the deviations between the results computed by different codes are of the same order of magnitude as the deviations between the computational and experimental results. On the one hand, this means that the HSE functional is rather accurate for describing the cohesive properties of simple insulating solids. On the other hand, this also suggests that the numerical precision achieved with current major HSE implementations is not sufficiently high to unambiguously assess the physical accuracy of HDFs. It is found that the pseudopotential treatment of the core electrons is a major factor that contributes to this uncertainty.

Massively parallel linear-scaling Hartree-Fock exchange and hybrid exchange-correlation functionals with plane wave basis set accuracy

We extend our linear-scaling approach for the calculation of Hartree-Fock exchange energy using localized in situ optimized orbitals [Dziedzic et al., J. Chem. Phys. 139, 214103 (2013)] to leverage massive parallelism. Our approach has been implemented in the onetep (Order-N Electronic Total Energy Package) density functional theory framework, which employs a basis of non-orthogonal generalized Wannier functions (NGWFs) to achieve linear scaling with system size while retaining controllable near-complete-basis-set accuracy. For the calculation of Hartree-Fock exchange, we use a resolution-of-identity approach, where an auxiliary basis set of truncated spherical waves is used to fit products of NGWFs. The fact that the electrostatic potential of spherical waves (SWs) is known analytically, combined with the use of a distance-based cutoff for exchange interactions, leads to a calculation cost that scales linearly with the system size. Our new implementation, which we describe in detail, combines distributed memory parallelism (using the message passing interface) with shared memory parallelism (OpenMP threads) to efficiently utilize numbers of central processing unit cores comparable to, or exceeding, the number of atoms in the system. We show how the use of multiple time-memory trade-offs substantially increases performance, enabling our approach to achieve superlinear strong parallel scaling in many cases and excellent, although sublinear, parallel scaling otherwise. We demonstrate that in scenarios with low available memory, which preclude or limit the use of time-memory trade-offs, the performance degradation of our algorithm is graceful. We show that, crucially, linear scaling with system size is maintained in all cases. We demonstrate the practicability of our approach by performing a set of fully converged production calculations with a hybrid functional on large imogolite nanotubes up to over 1400 atoms. We finish with a brief study of how the employed approximations (exchange cutoff and the quality of the SW basis) affect the calculation walltime and the accuracy of the obtained results.

Straightforward and Accurate Automatic Auxiliary Basis Set Generation for Molecular Calculations with Atomic Orbital Basis Sets

Density fitting (DF), also known as the resolution of the identity (RI), is a widely used technique in quantum chemical calculations with various types of atomic basis sets─Gaussian-type orbitals, Slater-type orbitals, as well as numerical atomic orbitals─to speed up density functional, Hartree-Fock (HF), and post-HF calculations. Traditionally, custom auxiliary basis sets are hand-optimized for each orbital basis set; however, some automatic schemes have also been suggested. In this work, we propose a simple yet numerically stable automated scheme for forming auxiliary basis sets with the help of a pivoted Cholesky decomposition, which is applicable to any type of atomic basis function. We exemplify the scheme with proof-of-concept calculations with Gaussian basis sets and show that the proposed approach leads to negligible DF/RI errors in HF and second-order Møller-Plesset (MP2) total energies of the non-multireference part of the W4-17 test set when used with orbital basis sets of at least polarized triple-ζ quality. The results are promising for future applications employing Slater-type orbitals or numerical atomic orbitals, as well as schemes based on the use of local fitting approaches. Global fitting approaches can also be used, in which case the high-angular-momentum functions produced by the present scheme can be truncated to minimize the computational cost.

The mystery of sub-picosecond charge transfer following irradiation of hydrated uridine monophosphate

The early mechanisms by which ionizing rays damage biological structures by so-called direct effects are largely elusive. In a recent picosecond pulse radiolysis study of concentrated uridine monophosphate solutions [J. Ma, S. A. Denisov, J.-L. Marignier, P. Pernot, A. Adhikary, S. Seki and M. Mostafavi, J. Phys. Chem. Lett., 2018, 9, 5105], unexpected results were found regarding the oxidation of the nucleobase. The signature of the oxidized nucleobase could not be detected 5 ps after the electron pulse, but only the oxidized phosphate, raising intriguing questions about the identity of charge-transfer mechanisms that could explain the absence of U⁺. We address here this question by means of advanced first-principles atomistic simulations of solvated uridine monophosphate, combining Density Functional Theory (DFT) with polarizable embedding schemes. We contrast three very distinct mechanisms of charge transfer covering the atto-, femto- and pico-second timescales. We first investigate the ionization mechanism and subsequent hole/charge migrations on a timescale of attoseconds to a few femtoseconds under the frozen nuclei approximation. We then consider a nuclear-driven phosphate-to-oxidized-nucleobase electron transfer, showing that it is an uncompetitive reaction channel on the sub-picosecond timescale, despite its high exothermicity and significant electronic coupling. Finally, we show that non-adiabatic charge transfer is enabled by femtosecond nuclear relaxation after ionization. We show that electronic decoherence and the electronic coupling strength are the key parameters that determine the hopping probabilities. Our results provide important insight into the interplay between electronics and nuclear motions in the early stages of the multiscale responses of biological matter subjected to ionizing radiation.

An improved chain of spheres for exchange algorithm

In the present work, we describe a more accurate and efficient variant of the chain-of-spheres algorithm (COSX) for exchange matrix computations. Higher accuracy for the numerical integration is obtained with new grids that were developed using global optimization techniques. With our new default grids, the average absolute energy errors are much lower than 0.1 kcal/mol, which is desirable to achieve "chemical accuracy." Although the size of the new grids is increased by roughly a factor of 2.5, the excellent efficiency of the original COSX implementation is still further improved in most cases. The evaluation of the analytic electrostatic potential integrals was significantly accelerated by a new implementation of rolled-out versions of the Dupuis-Rys-King and Head-Gordon-Pople algorithms. Compared to our earlier implementation, a twofold speedup is obtained for the frequently used triple-ζ basis sets, while up to a 16-fold speedup is observed for quadruple-ζ basis sets. These large gains are a consequence of both the more efficient integral evaluation and the intermediate exchange matrix computation in a partially contracted basis when generally contracted shells occur. With our new RIJCOSX implementation, we facilitate accurate self-consistent field (SCF) binding energy calculations on a large supra-molecular complex composed of 320 atoms. The binding-energy errors with respect to the fully analytic results are well below 0.1 kcal/mol for the cc-pV(T/Q)Z basis sets and even smaller than for RIJ with fully analytic exchange. At the same time, our RIJCOSX SCF calculation even with the cc-pVQZ basis and the finest grid is 21 times faster than the fully analytic calculation.

First-Principles Simulations of Biological Molecules Subjected to Ionizing Radiation

Ionizing rays cause damage to genomes, proteins, and signaling pathways that normally regulate cell activity, with harmful consequences such as accelerated aging, tumors, and cancers but also with beneficial effects in the context of radiotherapies. While the great pace of research in the twentieth century led to the identification of the molecular mechanisms for chemical lesions on the building blocks of biomacromolecules, the last two decades have brought renewed questions, for example, regarding the formation of clustered damage or the rich chemistry involving the secondary electrons produced by radiolysis. Radiation chemistry is now meeting attosecond science, providing extraordinary opportunities to unravel the very first stages of biological matter radiolysis. This review provides an overview of the recent progress made in this direction, focusing mainly on the atto- to femto- to picosecond timescales. We review promising applications of time-dependent density functional theory in this context.

Efficient Hybrid Density Functional Calculations for Large Periodic Systems Using Numerical Atomic Orbitals

We present an efficient, linear-scaling implementation for building the (screened) Hartree-Fock exchange (HFX) matrix for periodic systems within the framework of numerical atomic orbital (NAO) basis functions. Our implementation is based on the localized resolution of the identity approximation by which two-electron Coulomb repulsion integrals can be obtained by only computing two-center quantities-a feature that is highly beneficial to NAOs. By exploiting the locality of basis functions and efficient prescreening of the intermediate three- and two-index tensors, one can achieve a linear scaling of the computational cost for building the HFX matrix with respect to the system size. Our implementation is massively parallel, thanks to a MPI/OpenMP hybrid parallelization strategy for distributing the computational load and memory storage. All these factors add together to enable highly efficient hybrid functional calculations for large-scale periodic systems. In this work, we describe the key algorithms and implementation details for the HFX build as implemented in the ABACUS code package. The performance and scalability of our implementation with respect to the system size and the number of CPU cores are demonstrated for selected benchmark systems up to 4096 atoms.

Low-Order Scaling G0W0 by Pair Atomic Density Fitting

We derive a low-scaling G0W0 algorithm for molecules using pair atomic density fitting (PADF) and an imaginary time representation of the Green's function and describe its implementation in the Slater type orbital (STO)-based Amsterdam density functional (ADF) electronic structure code. We demonstrate the scalability of our algorithm on a series of water clusters with up to 432 atoms and 7776 basis functions and observe asymptotic quadratic scaling with realistic threshold qualities controlling distance effects and basis sets of triple-ζ (TZ) plus double polarization quality. Also owing to a very small prefactor, a G0W0 calculation for the largest of these clusters takes only 240 CPU hours with these settings. We assess the accuracy of our algorithm for HOMO and LUMO energies in the GW100 database. With errors of 0.24 eV for HOMO energies on the quadruple-ζ level, our implementation is less accurate than canonical all-electron implementations using the larger def2-QZVP GTO-type basis set. Apart from basis set errors, this is related to the well-known shortcomings of the GW space-time method using analytical continuation techniques as well as to numerical issues of the PADF approach of accurately representing diffuse atomic orbital (AO) products. We speculate that these difficulties might be overcome by using optimized auxiliary fit sets with more diffuse functions of higher angular momenta. Despite these shortcomings, for subsets of medium and large molecules from the GW5000 database, the error of our approach using basis sets of TZ and augmented double-ζ (DZ) quality is decreasing with system size. On the augmented DZ level, we reproduce canonical, complete basis set limit extrapolated reference values with an accuracy of 80 meV on average for a set of 20 large organic molecules. We anticipate our algorithm, in its current form, to be very useful in the study of single-particle properties of large organic systems such as chromophores and acceptor molecules.

Double hybrid DFT calculations with Slater type orbitals

On a comprehensive database with 1,644 datapoints, covering several aspects of main-group as well as of transition metal chemistry, we assess the performance of 60 density functional approximations (DFA), among them 36 double hybrids (DH). All calculations are performed using a Slater type orbital (STO) basis set of triple-ζ (TZ) quality and the highly efficient pair atomic resolution of the identity approach for the exchange- and Coulomb-term of the KS matrix (PARI-K and PARI-J, respectively) and for the evaluation of the MP2 energy correction (PARI-MP2). Employing the quadratic scaling SOS-AO-PARI-MP2 algorithm, DHs based on the spin-opposite-scaled (SOS) MP2 approximation are benchmarked against a database of large molecules. We evaluate the accuracy of STO/PARI calculations for B3LYP as well as for the DH B2GP-PLYP and show that the combined basis set and PARI-error is comparable to the one obtained using the well-known def2-TZVPP Gaussian-type basis set in conjunction with global density fitting. While quadruple-ζ (QZ) calculations are currently not feasible for PARI-MP2 due to numerical issues, we show that, on the TZ level, Jacob's ladder for classifying DFAs is reproduced. However, while the best DHs are more accurate than the best hybrids, the improvements are less pronounced than the ones commonly found on the QZ level. For conformers of organic molecules and noncovalent interactions where very high accuracy is required for qualitatively correct results, DHs provide only small improvements over hybrids, while they still excel in thermochemistry, kinetics, transition metal chemistry and the description of strained organic systems.

Variational Density Fitting with a Krylov Subspace Method

In this work, we present the implementation of a variational density fitting methodology that uses iterative linear algebra for solving the associated system of linear equations. It is well known that most difficulties with this system arise from the fact that the coefficient matrix is in general ill-conditioned and, due to finite precision round-off errors, it may not be positive definite. The dimensionality, given by the number of auxiliary functions, also poses a challenge in terms of memory and time demand since the coefficient matrix is dense. The methodology presented is based on a preconditioned Krylov subspace method able to deal with indefinite ill-conditioned equation systems. To assess its potential, it has been combined with double asymptotic electron repulsion integral expansions as implemented in the deMon2k package. A numerical study on a set of problems with up to 130,000 auxiliary functions shows its effectiveness to alleviate the abovementioned problematic. A comparison with the default methodology used in deMon2k based on a truncated eigenvalue decomposition of the coefficient matrix indicates that the proposed method exhibits excellent robustness and scalability when implemented in a parallel setting.

Accuracy of Localized Resolution of the Identity in Periodic Hybrid Functional Calculations with Numerical Atomic Orbitals

We present an implementation of hybrid density functional approximations for periodic systems within a pseudopotential-based, numerical atomic orbital (NAO) framework. The two-electron Coulomb repulsion integrals (ERIs) are evaluated using the localized resolution-of-the-identity (LRI) approximation. The accuracy of the LRI approximation is benchmarked unambiguously against independent reference results obtained via a computational scheme whereby the ERIs are accurately evaluated by expanding the products of NAOs in terms of plane waves. An alternative strategy for constructing auxiliary basis sets is proposed, and its accuracy is assessed and compared to the previously used procedure. Finally, the reliability of our algorithm and implementation is benchmarked against other established implementations within different numerical frameworks in terms of the calculated band gap values of a set of semiconductors and insulators.

A Quadratic Pair Atomic Resolution of the Identity Based SOS-AO-MP2 Algorithm Using Slater Type Orbitals

We report a production level implementation of pair atomic resolution of the identity (PARI) based second-order Møller-Plesset perturbation theory (MP2) in the Slater type orbital (STO) based Amsterdam Density Functional (ADF) code. As demonstrated by systematic benchmarks, dimerization and isomerization energies obtained with our code using STO basis sets of triple-ζ-quality show mean absolute deviations from Gaussian type orbital, canonical, basis set limit extrapolated, global density fitting (DF)-MP2 results of less than 1 kcal/mol. Furthermore, we introduce a quadratic scaling atomic orbital based spin-opposite-scaled (SOS)-MP2 approach with a very small prefactor. Due to a worst-case scaling of [Formula: see text], our implementation is very fast already for small systems and shows an exceptionally early crossover to canonical SOS-PARI-MP2. We report computational wall time results for linear as well as for realistic three-dimensional molecules and show that triple-ζ quality calculations on molecules of several hundreds of atoms are only a matter of a few hours on a single compute node, the bottleneck of the computations being the SCF rather than the post-SCF energy correction.

The physical stage of radiolysis of solvated DNA by high-energy-transfer particles: insights from new first principles simulations

Electron dynamics simulations based on density functional theory are carried out on nanometric molecular systems to decipher the primary processes following irradiation of bio-macromolecules by high energy transfer charged particles.

Density Functional Theory for Molecular and Periodic Systems Using Density Fitting and Continuous Fast Multipole Methods

An implementation of Kohn-Sham density functional theory within the TURBOMOLE program package with Gaussian-type orbitals (GTO) as basis functions is reported that treats molecular and periodic systems of any dimensionality on an equal footing. Its key component is a combination of density fitting/resolution of identity (DF) approximation and continuous fast multipole method (CFMM) applied for the electronic Coulomb term. This DF-CFMM scheme operates entirely in the direct space and partitions Coulomb interactions into far-field part evaluated using multipole expansions and near-field contribution calculated employing density fitting. Computational efficiency and favorable scaling behavior of our implementation approaching O(N) for the formation of Kohn-Sham matrix is demonstrated for various molecular and periodic systems including three-dimensional models with unit cells containing up to 640 atoms and 19072 GTO basis functions.

Size-dependent properties of transition metal clusters: from molecules to crystals and surfaces--computational studies with the program ParaGauss

In the so-called scalable regime the size-dependent behavior of the physical and chemical properties of transition metal clusters is described by scaling relationships. For most quantities this scalable regime is reached for cluster sizes between a few tens and a few hundreds of atoms, hence for systems for which an accurate treatment by density functional theory is still feasible. Thus, by invoking scaling relations one is able to obtain properties of very large nanoparticles and even the bulk limit from the results of a series of smaller cluster models. In this invited review we illustrate this strategy by exploiting results from computational studies that mostly were carried out with the density functional theory software ParaGauss. We address mainly the size-dependent behavior of the properties of transition metal clusters. To this end, we first present benchmark studies probing various approximations that are used in such density functional calculations. Subsequently we show how physical insight may be gained by exploring less understood types of systems. These applications range from bare clusters to nanoislands and nanoalloys to adsorption complexes.

Efficient calculation of nuclear spin-rotation constants from auxiliary density functional theory

The computation of the spin-rotation tensor within the framework of auxiliary density functional theory (ADFT) in combination with the gauge including atomic orbital (GIAO) scheme, to treat the gauge origin problem, is presented. For the spin-rotation tensor, the calculation of the magnetic shielding tensor represents the most demanding computational task. Employing the ADFT-GIAO methodology, the central processing unit time for the magnetic shielding tensor calculation can be dramatically reduced. In this work, the quality of spin-rotation constants obtained with the ADFT-GIAO methodology is compared with available experimental data as well as with other theoretical results at the Hartree-Fock and coupled-cluster level of theory. It is found that the agreement between the ADFT-GIAO results and the experiment is good and very similar to the ones obtained by the coupled-cluster single-doubles-perturbative triples-GIAO methodology. With the improved computational performance achieved, the computation of the spin-rotation tensors of large systems or along Born-Oppenheimer molecular dynamics trajectories becomes feasible in reasonable times. Three models of carbon fullerenes containing hundreds of atoms and thousands of basis functions are used for benchmarking the performance. Furthermore, a theoretical study of temperature effects on the structure and spin-rotation tensor of the H(12)C-(12)CH-DF complex is presented. Here, the temperature dependency of the spin-rotation tensor of the fluorine nucleus can be used to identify experimentally the so far unknown bent isomer of this complex. To the best of our knowledge this is the first time that temperature effects on the spin-rotation tensor are investigated.

Density perturbation theory

Despite the fundamental importance of electron density in density functional theory, perturbations are still usually dealt with using Hartree-Fock-like orbital equations known as coupled-perturbed Kohn-Sham (CPKS). As an alternative, we develop a perturbation theory that solves for the perturbed density directly, removing the need for CPKS. This replaces CPKS with a true Hohenberg-Kohn density perturbation theory. In CPKS, the perturbed density is found in the basis of products of occupied and virtual orbitals, which becomes ever more over-complete as the size of the orbital basis set increases. In our method, the perturbation to the density is expanded in terms of a series of density basis functions and found directly. It is possible to solve for the density in such a way that it makes the total energy stationary even if the density basis is incomplete.

Low-memory iterative density fitting

A new low-memory modification of the density fitting approximation based on a combination of a continuous fast multipole method (CFMM) and a preconditioned conjugate gradient solver is presented. Iterative conjugate gradient solver uses preconditioners formed from blocks of the Coulomb metric matrix that decrease the number of iterations needed for convergence by up to one order of magnitude. The matrix-vector products needed within the iterative algorithm are calculated using CFMM, which evaluates them with the linear scaling memory requirements only. Compared with the standard density fitting implementation, up to 15-fold reduction of the memory requirements is achieved for the most efficient preconditioner at a cost of only 25% increase in computational time. The potential of the method is demonstrated by performing density functional theory calculations for zeolite fragment with 2592 atoms and 121,248 auxiliary basis functions on a single 12-core CPU workstation.

Efficient calculation of the rotational g tensor from auxiliary density functional theory

The computation of the rotational g tensor with the recently developed auxiliary density functional theory (ADFT) gauge including atomic orbital (GIAO) methodology is presented. For the rotational g tensor, the calculation of the magnetizability tensor represents the most demanding computational task. With the ADFT-GIAO methodology, the CPU time for the magnetizability tensor calculation can be dramatically reduced. Therefore, it seems most desirable to employ the ADFT-GIAO methodology also for the computation of the rotational g tensor. In this work, the quality of rotational g tensors obtained with the ADFT-GIAO methodology is compared with available experimental data as well as with other theoretical results at the Hartree-Fock and coupled-cluster level of theory. It is found that the agreement between the ADFT-GIAO results and the experiment is good. Furthermore, we also show that the ADFT-GIAO g tensor calculation is applicable to large systems like carbon nanotube models containing hundreds of atom and thousands of basis functions.