Real-space Density Functional Theory Research Articles

Efficient real space formalism for hybrid density functionals.

We present an efficient real space formalism for hybrid exchange-correlation functionals in generalized Kohn-Sham density functional theory (DFT). In particular, we develop an efficient representation for any function of the real space finite-difference Laplacian matrix by leveraging its Kronecker product structure, thereby enabling the time to solution of associated linear systems to be highly competitive with the fast Fourier transform scheme while not imposing any restrictions on the boundary conditions. We implement this formalism for both the unscreened and range-separated variants of hybrid functionals. We verify its accuracy and efficiency through comparisons with established planewave codes for isolated as well as bulk systems. In particular, we demonstrate up to an order-of-magnitude speedup in time to solution for the real space method. We also apply the framework to study the structure of liquid water using abinitio molecular dynamics, where we find good agreement with the literature. Overall, the current formalism provides an avenue for efficient real-space DFT calculations with hybrid density functionals.

Improving the precision of forces in real-space pseudopotential density functional theory.

The high-order finite difference real-space pseudopotential density functional theory (DFT) approach is a valuable method for large-scale, massively parallel DFT calculations. A significant challenge in the approach is the oscillating "egg-box" error introduced by aliasing associated with a coarse grid spacing. To address this issue while minimizing computational cost, we developed a finite difference interpolation (FDI) scheme [Roller et al., J. Chem. Theory Comput. 19, 3889 (2023)] as a means of exploiting the high resolution of the pseudopotential to reduce egg-box effects systematically. Here, we show an implementation of this method in the PARSEC code and examine the practical utility of the combination of FDI with additional methods for improving force precision and/or reducing its computational cost, including orbital-based forces, compensating charges (namely, adding and subtracting a judiciously chosen charge density such that the total density is unaltered), and a modified spatial domain in which the real-space grid is defined. Using selected small molecules, as well as metallic Li, as test cases, we show that a combination of all four aspects leads to a significant reduction in computational cost while retaining a high level of precision that supports accurate structures and vibrational spectra, as well as stable and accurate molecular dynamics runs.

What can lattice DFT teach us about real-space DFT?

In this paper we establish a connection between density functional theory (DFT) for lattice models and common real-space DFT. We consider the lattice DFT description of a two-level model subject to generic interactions in Mermin's DFT formulation in the grand canonical ensemble at finite temperature. The case of only density-density and Hund's rule interaction studied in earlier work is shown to be equivalent to an exact-exchange description of DFT in the real-space picture. In addition, we also include the so-called pair-hopping interaction which can be treated analytically and, crucially, leads to non-integer occupations of the Kohn-Sham (KS) levels even in the limit of zero temperature. Treating the hydrogen molecule in a minimal basis is shown to be equivalent to our two-level lattice DFT model. By means of the fractional occupations of the KS orbitals (which, in this case, are identical to the many-body ones) we reproduce the results of full configuration interaction, even in the dissociation limit and without breaking the spin symmetry. Beyond the minimal basis, we embed our HOMO-LUMO model into a standard DFT calculation and, again, obtain results in overall good agreement with exact ones without the need of breaking the spin symmetry.

Calculation of phonons in real-space density functional theory.

We present an accurate and efficient formulation for the calculation of phonons in real-space Kohn-Sham density functional theory. Specifically, employing a local exchange-correlation functional, norm-conserving pseudopotential in the Kleinman-Bylander representation, and local form for the electrostatics, we derive expressions for the dynamical matrix and associated Sternheimer equationthat are particularly amenable to the real-space finite-difference method, within the framework of density functional perturbation theory. In particular, the formulation is applicable to insulating as well as metallic systems of any dimensionality, enabling the efficient and accurate treatment of semi-infinite and bulk systems alike, for both orthogonal and nonorthogonal cells. We also develop an implementation of the proposed formulation within the high-order finite-difference method. Through representative examples, we verify the accuracy of the computed phonon dispersion curves and density of states, demonstrating excellent agreement with established plane-wave results.

Finite Difference Interpolation for Reduction of Grid-Related Errors in Real-Space Pseudopotential Density Functional Theory.

The real-space pseudopotential approach is a well-known method for large-scale density functional theory (DFT) calculations. One of its main limitations, however, is the introduction of errors associated with the positioning of the underlying real-space grid, a phenomenon usually known as the "egg-box" effect. The effect can be controlled by using a finer grid, but this raises the cost of the calculations or even undermines their feasibility altogether. Therefore, there is ongoing interest in the reduction of the effect per a given real-space grid. Here, we present a finite difference interpolation of electron orbitals as a means of exploiting the high resolution of the pseudopotential to reduce egg-box effects systematically. We implement the method in PARSEC, a finite difference real-space pseudopotential DFT code, and demonstrate error mitigation and improved convergence at a low additional computational cost.

Chemical Bonding in Large Systems Using Projected Population Analysis from Real-Space Density Functional Theory Calculations.

We present an efficient and scalable computational approach for conducting projected population analysis from real-space finite-element (FE)-based Kohn-Sham density functional theory calculations (DFT-FE). This work provides an important direction toward extracting chemical bonding information from large-scale DFT calculations on materials systems involving thousands of atoms while accommodating periodic, semiperiodic, or fully nonperiodic boundary conditions. Toward this, we derive the relevant mathematical expressions and develop efficient numerical implementation procedures that are scalable on multinode CPU architectures to compute the projected overlap and Hamilton populations. The population analysis is accomplished by projecting either the self-consistently converged FE discretized Kohn-Sham orbitals or the FE discretized Hamiltonian onto a subspace spanned by a localized atom-centered basis set. The proposed methods are implemented in a unified framework within the DFT-FE code where the ground-state DFT calculations and the population analysis are performed on the same FE grid. We further benchmark the accuracy and performance of this approach on representative material systems involving periodic and nonperiodic DFT calculations with LOBSTER, a widely used projected population analysis code. Finally, we discuss a case study demonstrating the advantages of our scalable approach to extract the quantitative chemical bonding information of hydrogen chemisorbed in large silicon nanoparticles alloyed with carbon, a candidate material for hydrogen storage.

Efficient single-grid and multi-grid solvers for real-space orbital-free density functional theory

To improve the computational efficiencies of the real-space orbital-free density functional theory, this work develops a new single-grid solver by directly providing the closed-form solution to the inner iteration and using an improved bisection method to accelerate the line search process in the outer iteration, and extended the single-grid solver to a multi-grid solver. Numerical examples show that while maintaining high level of accuracy, the proposed single-grid solver can improve the computational efficiencies by two orders of magnitude comparing with the methods in the literature (For example, the ground state energy of a 45×45×45 FCC Aluminum cluster (376786 atoms) can be computed within 4h27min on a desktop computer.) and the multi-grid solver can improve the computational efficiencies even once for the cases where high-resolution electron densities are needed.

Gaussian-Approximated Poisson Preconditioner for Iterative Diagonalization in Real-Space Density Functional Theory

Various real-space methods optimized on massive parallel computers have been developed for efficient large-scale density functional theory (DFT) calculations of materials and biomolecules. The iterative diagonalization of the Hamiltonian matrix is a computational bottleneck in real-space DFT calculations. Despite the development of various iterative eigensolvers, the absence of efficient real-space preconditioners has hindered their overall efficiency. An efficient preconditioner must satisfy two conditions: appropriate acceleration of the convergence of the iterative process and inexpensive computation. This study proposed a Gaussian-approximated Poisson preconditioner (GAPP) that satisfied both conditions and was suitable for real-space methods. A low computational cost was realized through the Gaussian approximation of a Poisson Green's function. Fast convergence was achieved through the proper determination of Gaussian coefficients to fit the Coulomb energies. The performance of GAPP was evaluated for several molecular and extended systems, and it showed the highest efficiency among the existing preconditioners adopted in real-space codes.

Accelerating self-consistent field iterations in Kohn-Sham density functional theory using a low-rank approximation of the dielectric matrix

We present an efficient preconditioning technique for accelerating the fixed point iteration in real-space Kohn-Sham density functional theory (DFT) calculations. The preconditioner uses a low rank approximation of the dielectric matrix (LRDM) based on G\^ateaux derivatives of the residual of fixed point iteration along appropriately chosen direction functions. We develop a computationally efficient method to evaluate these G\^ateaux derivatives in conjunction with the Chebyshev filtered subspace iteration procedure, an approach widely used in large-scale Kohn-Sham DFT calculations. Further, we propose a variant of LRDM preconditioner based on adaptive accumulation of low-rank approximations from previous SCF iterations, and also extend the LRDM preconditioner to spin-polarized Kohn-Sham DFT calculations. We demonstrate the robustness and efficiency of the LRDM preconditioner against other widely used preconditioners on a range of benchmark systems with sizes ranging from $\sim$ 100-1100 atoms ($\sim$ 500--20,000 electrons). The benchmark systems include various combinations of metal-insulating-semiconducting heterogeneous material systems, nanoparticles with localized $d$ orbitals near the Fermi energy, nanofilm with metal dopants, and magnetic systems. In all benchmark systems, the LRDM preconditioner converges robustly within 20--30 iterations. In contrast, other widely used preconditioners show slow convergence in many cases, as well as divergence of the fixed point iteration in some cases. Finally, we demonstrate the computational efficiency afforded by the LRDM method, with up to 3.4$\times$ reduction in computational cost for the total ground-state calculation compared to other preconditioners.

TTDFT: A GPU accelerated Tucker tensor DFT code for large-scale Kohn-Sham DFT calculations

We present the Tucker tensor DFT (TTDFT) code which uses a tensor-structured algorithm with graphic processing unit (GPU) acceleration for conducting ground-state DFT calculations on large-scale systems. The Tucker tensor DFT algorithm uses a localized Tucker tensor basis computed from an additive separable approximation to the Kohn-Sham Hamiltonian. The discrete Kohn-Sham problem is solved using Chebyshev filtered subspace iteration method that relies on matrix-matrix multiplications of a sparse symmetric Hamiltonian matrix and a dense wavefunction matrix, expressed in the localized Tucker tensor basis. These matrix-matrix multiplication operations, which constitute the most computationally intensive step of the solution procedure, are GPU accelerated providing ∼8-fold GPU-CPU speedup for these operations on the largest systems studied. The computational performance of the TTDFT code is presented using benchmark studies on aluminum nano-particles and silicon quantum dots with system sizes ranging up to ∼7,000 atoms. Program summaryProgram Title: TTDFT: Tucker tensor density functional theory codeCPC Library link to program files:https://doi.org/10.17632/8dgmcs8ys2.1Licensing provisions: LGPLProgramming language: C/C++External routines/libraries: TuckerMPI (https://gitlab.com/tensors/TuckerMPI), cuBLAS (https://docs.nvidia.com/cuda/cublas/index.html), cuSparse (https://docs.nvidia.com/cuda/cusparse/index.html), ALGLIB (http://www.alglib.net/), Boost (https://www.boost.org/), BLAS (http://www.netlib.org/blas/), LAPACK (http://www.netlib.org/lapack/), PETSc (https://www.mcs.anl.gov/petsc), SLEPc (http://slepc.upv.es)Nature of problem: Real-space Kohn-Sham density functional theory calculations using localized Tucker tensor basis.Solution method: We present a real-space Kohn-Sham density functional theory code based on tensor-structured techniques with GPU acceleration. Tensor-structured techniques are adopted for computing a Tucker tensor basis, representing the eigenfunctions of an additive separable approximation to the Kohn-Sham Hamiltonian. The Tucker tensor basis is further localized using L1 regularization to improve the sparsity of the Kohn-Sham Hamiltonian matrix, and improve the computational efficiency and parallel scalability of the proposed algorithm. The solution to the Kohn-Sham problem in the localized Tucker tensor basis is computed using the Chebyshev filtered subspace iteration (ChFSI) method.Additional comments including restrictions and unusual features: The code works with Troullier-Martin (TM) pseudopotentials in Kleinman-Bylander form. The current release supports only non-periodic DFT calculations with the local density approximation (LDA) for exchange-correlation functional.This TTDFT project uses GitHub via Git, a free distributed version control software. The archived version at the time of submission of this work can be found on the CPC program library through program files DOI provided above. The GitHub repository of this project can be found on https://github.com/ttdftdev/ttdft_public.

Ab-Initio Investigation of Finite Size Effects in Rutile Titania Nanoparticles with Semilocal and Nonlocal Density Functionals

In this work, we employ hybrid and generalized gradient approximation (GGA) level density functional theory (DFT) calculations to investigate the convergence of surface properties and the electronic gap of rutile titania (TiO2) nanoparticles with particle size. The surface energies and electronic gaps are calculated for cuboidal particles with minimum dimension ranging from 3.7 Å (24 atoms) to 10.3 Å (384 atoms) using a highly parallel real-space DFT code to enable hybrid-level DFT calculations of larger nanoparticles than are typically practical. We deconvolute the geometric and electronic finite size effects in surface energy and evaluate the influence of defects on the electronic gap and density of states (DOS). The electronic finite size effects in surface energy vanish when the minimum length scale of the nanoparticles becomes greater than 10 Å. We show that this length scale is consistent with a computationally efficient numerical analysis of the characteristic length scale of electronic interactions. The surface energy of nanoparticles having a minimum dimension beyond this characteristic length can be approximated using slab calculations that account for the geometric defects. In contrast, the finite size effects on the electronic gap and DOS are highly dependent on the shape and size of these particles. The DOS for cuboidal particles and more realistic particles constructed using the Wulff algorithm reveal that defect states within the electronic gap play a key role in determining the eigenvalue distribution of nanoparticles and the electronic gap does not converge to the bulk limit for the particle sizes investigated.

Microscopic mechanism of adatom diffusion on stepped SiC surfaces revealed by first-principles calculations

We report first-principles total-energy calculations based on real-space density-functional theory that unveil the atom-scale mechanisms of surface diffusion of adatoms on the Si-faced 3C-SiC(111) stepped surface. The quality of the epitaxial layer of SiC affects the device performance. Therefore, fundamental knowledge of the microscopic mechanisms of epitaxial growth is crucial for the improvement of the quality of SiC power devices. However, adatom diffusion on the growing stepped SiC surfaces is still unknown. We identify diffusion pathways for three important adatom species of SiC chemical vapor deposition (CVD), Si, C, and H, and obtain the corresponding energy profiles for both on the surface terraces and near the surface steps, providing a complete picture of the adatom diffusion. We find that the Si adatom is most mobile on the terrace and shows prominent Ehrlich-Schwoebel (ES) effect in the inter-terrace diffusion, whereas the C and H adatoms show less ES effect and in an ascending diffusion even the inverse ES effect appears for the C adatom. The results obtained are fundamentals to explore the microscopic mechanism of the epitaxial growth on 3C-SiC(111) surfaces and equivalently hexagonal SiC(0001) surfaces.

Band engineering of large scale graphene/hexagonal boron nitride in-plane heterostructure: Role of the connecting angle

The modulation of the interface topography in two-dimensional (2D) lateral heterostructure (LHS) is significant for adjusting their electronic properties. However, previous theoretical studies were usually based on dozens of atoms and ignore the irregular interface configurations. In this work, we constructed a series of graphene/hexagonal boron nitride (G/h-BN) LHS of various misorientation angles, with thousands of atoms and uniform atomic proportion. The influence of connecting angles on the electronic properties has been investigated in detail by using real-space density functional theory. Results show that the bandgap can be modified efficiently by changing the connecting angle between graphene and h-BN in LHS, on account of the mixture of diverse interfacial atomic configurations. The wave function characteristics of different interfaces can be interpreted by the strong quantum confinement and gauge field for the edge states. These findings provide a theoretical basis for elucidating the relationship between the atomic construction and electronic properties of planar G/h-BN heterostructures, which could pave a way for further controllable and tunable 2D electronic devices.

Real-space density functional theory adapted to cyclic and helical symmetry: Application to torsional deformation of carbon nanotubes

We present a cyclic and helical symmetry-adapted formulation and large-scale parallel implementation of real-space Kohn-Sham density functional theory for one-dimensional (1D) nanostructures, with application to the mechanical and electronic response of carbon nanotubes subject to torsional deformations. Specifically, employing a semilocal exchange correlation and a local formulation of the electrostatics, we derive symmetry-adapted variants for the energy functional, variational problem governing the electronic ground state, Kohn-Sham equations, atomic forces, and axial stress, all posed on the fundamental domain. In addition, we develop a representation for twisted nanotubes of arbitrary chirality within this framework. We also develop a high-order finite-difference parallel implementation capable of performing accurate cyclic and helical symmetry-adapted Kohn-Sham calculations in both the static and dynamic settings, and verify it through numerical tests and comparisons with established codes. We use this implementation to perform twist-controlled simulations for a representative set of achiral and chiral carbon nanotubes, in both the small and large deformation regimes. In the linear regime, we find that the torsional moduli are proportional to the cube of the diameter; metallic nanotubes undergo metal-insulator transitions; and both the band gap as well as effective mass of charge carriers are proportional to the shear strain and sine of three times the chiral angle. In the nonlinear regime, we find that there is significant Poynting effect, particularly at the ultimate strain, the value of which is determined by the chiral angle; torsional deformations provide a possible mechanism for the irreversible phase transformation from armchair to zigzag nanotubes; and both the band gap as well as effective mass have an oscillatory behavior, with the period for metal-insulator transitions being inversely proportional to the square of the diameter and sine of three times the chiral angle. Wherever available, the results are in good agreement with experimental observations and measurements. Overall, this opens an avenue for the highly accurate and efficient first-principles study of 1D nanostructures that have cyclic and/or helical symmetry, as well as their response to torsional deformations.

Negative resistance for colloids driven over two barriers in a microchannel.

When considering the flow of currents through obstacles, one core expectation is that the total resistance of sequential single resistors is additive. While this rule is most commonly applied to electronic circuits, it also applies to other transport phenomena such as the flow of colloids or nanoparticles through channels containing multiple obstacles, as long as these obstacles are sufficiently far apart. Here we explore the breakdown of this additivity for fluids of repulsive colloids driven over two energetic barriers in a microchannel, using real-space microscopy experiments, particle-resolved simulations, and dynamical density functional theory. If the barrier separation is comparable to the particle correlation length, the resistance is highly non-additive, such that the resistance added by the second barrier can be significantly higher or lower than that of the first. Surprisingly, in some cases the second barrier can even add a negative resistance, such that two identical barriers are easier to cross than a single one. We explain this counterintuitive observation in terms of the structuring of particles trapped between the barriers.

Real-space formulation of the stress tensor for O(N) density functional theory: Application to high temperature calculations.

We present an accurate and efficient real-space formulation of the Hellmann-Feynman stress tensor for O(N) Kohn-Sham density functional theory (DFT). While applicable at any temperature, the formulation is most efficient at high temperature where the Fermi-Dirac distribution becomes smoother and the density matrix becomes correspondingly more localized. We first rewrite the orbital-dependent stress tensor for real-space DFT in terms of the density matrix, thereby making it amenable to O(N) methods. We then describe its evaluation within the O(N) infinite-cell Clenshaw-Curtis Spectral Quadrature (SQ) method, a technique that is applicable to metallic and insulating systems, is highly parallelizable, becomes increasingly efficient with increasing temperature, and provides results corresponding to the infinite crystal without the need of Brillouin zone integration. We demonstrate systematic convergence of the resulting formulation with respect to SQ parameters to exact diagonalization results and show convergence with respect to mesh size to the established plane wave results. We employ the new formulation to compute the viscosity of hydrogen at 106 K from Kohn-Sham quantum molecular dynamics, where we find agreement with previous more approximate orbital-free density functional methods.

On preconditioning the self-consistent field iteration in real-space Density Functional Theory

We present a real-space formulation for isotropic Fourier-space preconditioners used to accelerate the self-consistent field iteration in Density Functional Theory calculations. Specifically, after approximating the preconditioner in Fourier space using a rational function, we express its real-space application in terms of the solution of sparse Helmholtz-type systems. Using the truncated-Kerker and Resta preconditioners as representative examples, we show that the proposed real-space method is both accurate and efficient, requiring the solution of a single linear system, while accelerating self-consistency to the same extent as its exact Fourier-space counterpart.

Symmetry-adapted real-space density functional theory for cylindrical geometries: Application to large group-IV nanotubes

We present a symmetry-adapted real-space formulation of Kohn-Sham density functional theory for cylindrical geometries and apply it to the study of large X (X=C, Si, Ge, Sn) nanotubes. Specifically, starting from the Kohn-Sham equations posed on all of space, we reduce the problem to the fundamental domain by incorporating cyclic and periodic symmetries present in the angular and axial directions of the cylinder, respectively. We develop a high-order finite-difference parallel implementation of this formulation, and verify its accuracy against established planewave and real-space codes. Using this implementation, we study the band structure and bending properties of X nanotubes and Xene sheets, respectively. Specifically, we first show that zigzag and armchair X nanotubes with radii in the range 1 to 5 nm are semiconducting. In particular, we find an inverse linear dependence of the bandgap with respect to the radius for all nanotubes, other than the armchair and zigzag type III carbon variants, for which we find an inverse quadratic dependence. Next, we exploit the connection between cyclic symmetry and uniform bending deformations to calculate the bending moduli of Xene sheets in both zigzag and armchair directions. We find Kirchhoff-Love type bending behavior for all sheets, with graphene and stanene possessing the largest and smallest moduli, respectively. In addition, other than graphene, the sheets demonstrate significant anisotropy, with larger bending moduli along the armchair direction. Finally, we demonstrate that the proposed approach has very good parallel scaling and is highly efficient, enabling ab initio simulations of unprecedented size for systems with a high degree of cyclic symmetry. In particular, we show that even micron-sized nanotubes can be simulated with modest computational effort.

On the calculation of the stress tensor in real-space Kohn-Sham density functional theory.

We present an accurate and efficient formulation of the stress tensor for real-space Kohn-Sham density functional theory calculations. Specifically, while employing a local formulation of the electrostatics, we derive a linear-scaling expression for the stress tensor that is applicable to simulations with unit cells of arbitrary symmetry, semilocal exchange-correlation functionals, and Brillouin zone integration. In particular, we rewrite the contributions arising from the self-energy and the nonlocal pseudopotential energy to make them amenable to the real-space finite-difference discretization, achieving up to three orders of magnitude improvement in the accuracy of the computed stresses. Using examples representative of static and dynamic calculations, we verify the accuracy and efficiency of the proposed formulation. In particular, we demonstrate high rates of convergence with spatial discretization, consistency between the computed energy and the stress tensor, and very good agreement with reference planewave results.

Structural determination of phosphosilicate glass based on first-principles molecular dynamics calculation

We report static and dynamic first-principles calculations determining the structure and electronic properties of phosphosilicate glass (PSG). A stable structure of PSG with composition of (SiO2)x(P2O5)y was calculated by simulated annealing based on Car–Parrinello molecular dynamics (CPMD) with a real-space density functional theory (RSDFT) code. As a result, phosphorus (P) atoms in PSG were four-fold coordinated with oxygen (O) atoms and one of the four PO bonds was dangling toward the interstitial region (an –O3PO configuration). We also compared the structures having a different number of O atoms to discuss the stable configuration of PSG in the realistic situation of phosphoryl chloride (POCl3) annealing. We found that –O3PO is most stable for O2-rich conditions, whereas three-coordinate P (an –O3P configuration) becomes stable in O2-poor conditions. The candidates of the states which create levels inside the bandgap of PSG were also calculated.