Density Functional Theory and Materials Modeling at Atomistic Length Scales

We discuss the basic concepts of density functional theory (DFT) as applied to materials modeling in the microscopic, mesoscopic and macroscopic length scales. The picture that emerges is that of a single unified framework for the study of both quantum and classical systems. While for quantum DFT, the central equation is a one-particle Schrodinger-like Kohn-Sham equation, the classical DFT consists of Boltzmann type distributions, both corresponding to a system of noninteracting particles in the field of a density-dependent effective potential, the exact functional form of which is unknown. One therefore approximates the exchange-correlation potential for quantum systems and the excess free energy density functional or the direct correlation functions for classical systems. Illustrative applications of quantum DFT to microscopic modeling of molecular interaction and that of classical DFT to a mesoscopic modeling of soft condensed matter systems are highlighted.


Introduction
Designing materials with tailored properties has always been a long cherised dream and even its partial realisation heavily depends on the research on modeling and simulation.At the fundamental level, the properties of materials are of course ultimately determined by the nature of the electron density reorganisation that takes place during their formation from the constituent atoms, and hence in principle are obtainable through a detailed quantum mechanical calculation.However, a complete ab-initio calculation by solving the many-electron Schrodinger equation or its variants even with the state of the art computational advances is often too ambitious for the wide range of length scales involved in materials modeling and considerations involving a separation of length scales become essential.Thus, in the shorter microscopic length scales, one can aim for a complete ab-initio electronic structure description considering interacting electrons and the nuclei using quantum mechanics within the framework of approximate theories that have been developed over the years.On the other hand, for the intermediate mesoscopic length scale, an atomistic description with the concept of particles (atoms and molecules), their movements and distributions using classical and statistical mechanics is often most suitable.The interparticle interaction which is needed as input in this description can be generated through a quantum mechanical calculation within the microscopic length scale.The properties obtained at the mesoscopic length scale can in turn serve as input for the investigation in the macroscopic length scales, where one considers matter as a continuous medium and the conventional approaches of continuum mechanics and hydrodynamics of classical physics are used as theoretical tools.The objective of this work is to present a unified theoretical framework for the study of materials at these three length scales encompassing electronic, atomistic as well as continuum descriptions.Besides modeling the same material using different length scales, the recent upsurge in nanomaterials of the mesoscopic size regime has also led to greater opportunities for materials modeling.
Even in the microscopic length scale, the solution for the wavefunction of a many-electron system is difficult and methods such as density functional theory (DFT) [1][2][3][4][5] which employs the single-particle electron density [6] as the basic variable bypassing the many-electron wavefunction have been known to provide valuable tools for the electronic structure calculations.The scope of DFT has been further broadened [7] through its application in the mesoscopic domain where the density distribution of the atomistic particles serve as the basic variable.A density based description as used, for example, in hydrodynamics is also possible in the domain of macroscopic length scale.In this work, we present a unified view of DFT using the electronic, atomistic and macroscopic density variables for the respective length scales.
In the microscopic scale, one deals with the electron density based descriptions leading to conceptual simplicity, computational economy and a framework for strengthening the foundation of widely used chemical concepts.In the mesoscopic scale, the corresponding basic variable is the single-particle density of the atoms (or molecules) and the effective interaction potentials can be obtained either by a detailed DFT calculation using the electron density or by modeling the chemical binding in terms of the atomic parameters defined within the DFT framework.The field of soft condensed matter [8] physics has found tremendous applications of this mesoscopic classical DFT where mainly empirical interaction potentials are used.An important field dealing with the DFT based ab-initio simulation proposed by Car and Parrinello [9] has also emerged where the interatomic force is obtained using DFT for the electrons during simulation instead of an empirical pair potential.In the macroscopic domain, the mass, current, energy and other property densities are the basic variables for a continuous medium description.While the Car-Parrinello method [9] bridges the microscopic and mesoscopic length scales, many of the areas of materials science and engineering dealing with the material properties needs a bridging of the mesoscopic and macroscopic length scales as well.
Thus, the building blocks for a density functional treatment in the short, intermediate and large length scales are respectively the electron, atom and the volume element, with the corresponding basic variables being the electron density, the single particle density and the property density.Most of the applications of DFT involves the common feature of a collection of building blocks with an inhomogeneous density distribution arising due to the field of an external potential.There also is the density inhomogeneity arising through stabilisation of thermal density fluctuations, examples being the phenomena of freezing, nucleation, crystallisation etc.A unified view of the DFT encompassing the diverse situations and length scales as mentioned above is presented here for quantum as well as classical systems.We discuss here the basic structure of DFT as applied to quantum and classical systems covering the microscopic and mesoscopic length scales, with the main concern being on the atomistic length scale of materials modeling.

Density Functional Theory of Quantum and Classical Systems
The basic variable in DFT is the single-particle density ρ(r) which is defined, for an N -particle system, by integrating the N -particle distribution function and represents locally the number of particles per unit volume, thus integrating to the total number N as ρ(r)dr = N .This concept of single-particle density is valid for both quantum and classical systems and for the former, one has P (r 1 , r 2 , . . ., r N ) = |ψ(r 1 , r 2 , . . ., r N )| 2 where ψ is the manyelectron wavefunction.The possibility of describing a quantum or classical many-particle system completely in terms of the reduced one-particle density has been rigorously proved by Hohenberg and Kohn [1] and Mermin [2] (HKM) who demonstrated for the first time a one to one mapping between the density and the potential.For a many-particle system characterised by an external potential v(r) (arising due to the nuclei in the case of electrons and due to the walls or pores in the case of confined classical fluid particles), the ground state energy E v [ρ] (for a many-electron quantum system) or, the grand potential Ω v [ρ] (for a many-particle classical system) can be expressed as the unique functionals of density, given respectively by and where µ is the chemical potential.The respective true equilibrium densities are determined by the Euler-Lagrange equation which corresponds to the minima of the energy density functionals.For many-electron systems, the quantity F [ρ] is always a universal functional of density whereas for classical particles, it is so only for a specified interparticle interaction.For many-electron systems, F [ρ] can be expressed as where the functionals T s [ρ], E coul [ρ] and E xc [ρ] represent the noninteracting kinetic energy, classical Coulomb energy and the exchange-correlation (XC) energy respectively.The energy functional E v [ρ] with the known exact expressions for T s [ρ] and E coul [ρ] as given by yields, on minimization, the effective one-particle Kohn-Sham equation [3] given by where the density ρ(r) is obtained as the sum and the effective potential and the energy functional are given respectively by and For a classical fluid [7] with inhomogeneous density distribution, the functional F [ρ] representing the intrinsic Helmholtz free energy can be expressed as where the ideal-gas free-energy functional F id [ρ] (corresponding to absence of internal interactions) is the analogue of the noninteracting kinetic energy T s [ρ] of Eqs. ( 5) and ( 6), and is given by the exact expression where β 0 (= 1/k B T ) is the inverse temperature and Λ is the thermal de-Broglie wavelength.The quantity F ex [ρ] representing the excess free-energy for the classical system is analogous to the interaction energy functional ( ) of the quantum system and is a universal functional of density for a specified interparticle interaction.The Euler-Lagrange equation ( 4) for this classical system is given by which essentially represents a Boltzmann-like distribution of an ideal gas in the field of an effective potential where the first order direct correlation function (DCF) c (1) (r) is defined as the functional derivative This extra contribution of the first order DCF c (1) (r; [ρ]) to the effective potential, arising from interparticle interactions and correlations is analogous to the exchange-correlation contribution to the Kohn-Sham potential for quantum systems.It may be noted that a spin-polarized quantum system and a two-component classical fluid mixture can be treated using an extended version of the above respective frameworks with two density components as the basic variables.
In both quantum and classical DFT, the one-particle picture has emerged through an exact mapping of the actual system of N interacting particles in the field of an external potential v(r) into another system of N noninteracting particles of the same density ρ(r) but moving in an effective potential v ef f (r; [ρ]), which itself depends on the density, thus requiring a self-consistent iterative procedure for solution of the Kohn-Sham or the Boltzmann-like density equations.
Due to lack of knowledge of the exact form of v ef f (r; [ρ]) as a density functional for inhomogeneous systems, it is however essential to approximate the XC energy functional E xc [ρ] for quantum systems, and the excess free energy functional F ex [ρ] or its derivatives for classical systems.The knowledge of the functionals for specific systems with homogeneous density are often useful in approximating the functionals for the corresponding inhomogeneous systems.For many-electron quantum systems, a standard simple scheme is the local density approximation [10] (LDA) where the expression for the XC energy functional of the homogeneous system is directly evaluated using the inhomogeneous density, viz.
where ˜ xc denotes the XC energy per electron for the homogeneous electron gas.This approach is however not suitable for approximating the excess free energy of the classical systems such as liquids for which the particles are of finite size and the peak densities might even exceed the density of close packing.For such systems, the corresponding procedure is the weighted density approximation (WDA) [11][12][13] where the actual density distribution ρ(r) is smoothed out by coarse graining with a weight function w(r, r ; ρ(r)) to obtain an effective density ρ(r) as ρ(r) = dr ρ(r )w(r, r ; ρ(r)) (17) which is then used to obtain an approximation to the excess free energy or the first order DCF as by evaluating the corresponding expressions for fex and c(1) of the homogeneous fluid.Other procedures involve the gradient corrections to the LDA in electronic structure theory [14] and functional perturbation techniques [15,16] in the theory of complex fluids.

Electron Density and Microscopic Modeling of Intra-and Inter-molecular Interaction Potential
One of the important areas of research in materials modeling involves computer simulation using interparticle potentials as input which are obtained either empirically or through a quantum mechanical calculation.While the interatomic or intermolecular potential can be obtained through a direct detailed DFT calculation based on electron density with suitable approximations for the XC energy functionals, other simplified alternatives include a perturbative approach involving the electron density reorganisation arising out of interatomic or intermolecular perturbations.Thus, considering the up-and down-spin densities ρ α (r) and ρ β (r), (which integrate to the corresponding number of electrons N α and N β ), as the basic variables, the spin-polarised version of the energy functional of Eq. ( 2) can be written as Now if the system is perturbed by a perturbing potential δv(r) with consequent density changes δρ α (r) and δρ β (r), the resulting energy change can be expressed [17,18] using the functional Taylor expansion (retaining terms upto second order) as which on functional differentiation and using the result [19] (δ 2 E/δρ ν (r)δv(r )) = δ(r − r ), yields the up-and down-spin chemical potentials given (upto first order) by where µ 0 α and µ 0 β denote the corresponding chemical potentials of the unperturbed system, viz.
also identified [20] as the negative of the two electronegativities (∂E/∂N α ) and (∂E/∂N β ) respectively.These are the key equations providing a justification for the chemical potential equalisation within a molecular species.Here the hardness kernel [21] η µν (r, r ) represents the energy functional derivative with the subscripts µ and ν denoting any of the spins α and β, and determines the local hardness [22] η µν (r) as well as the cross hardness η αβ = (∂ 2 E/∂N α ∂N β ) as where the spin-dependent Fukui function [23] f ν (r) = (∂ρ ν (r)/∂N ν ).
The expressions for the energy and chemical potential in terms of density perturbation are useful in the studies of molecule formation, response properties and also intermolecular interaction.The basic approach is the same in all the cases and hence we now specialise to the simple case of molecule formation from its constituent atoms.The objective is to develop a lattice model with the lattice sites given by the locations of the M atoms and P bonds between them at the fixed positions {R i } corresponding to a particular (say, equilibrium) configuration of the molecule and therefore the expressions of chemical potential given by Eqs. ( 22) and ( 23) are now recast to obtain the same for the i-th atomic or bond site as given by with a similar expression for µ β,i obtained by interchanging α with β.This expression corresponds to the partitioning of the density changes δρ ν (r) on molecule formation as a sum of the atomic and bond region components as δρ ν (r) = i δρ ν,i (r) = i δρ ν,i (r i ), with r i denoting the atom or bond region around the i-th site location R i .Also, it is assumed without any loss of generality that the density components δρ ν,i (r) vanish outside this region r i belonging to the i-th site.Introducing the assumption η µν (r i , r j ) ≈ η µν (R i , R j ) for the hardness kernel, the chemical potential for the i-th site as given by Eq. ( 28) can be further simplified as with a similar expression for µ β,i given by where q ν,j denotes the atomic or bond site charges.One thus arrives at the lattice model of a molecular species with point charges (zeroeth moment of the density) located at the lattice sites.An equalisation of the effective chemical potentials of each spin as given by Eqs. ( 29) and ( 30) for all the M atomic and P bond sites leads to M + P − 1 linear equations in the charge variables for each spin which along with the charge conservation yield the individual site charges.The expression for the energy change as given by Eq. ( 21) as reexpressed in terms of these charges is given by While for molecule formation, δv(i) = 0, one can calculate the polarizability [18,24] or chemical reactivity etc, by considering this contribution [18].The energy expression of Eq. ( 31) is reminiscent of semiempirical quantum chemistry and numerical results based on this procedure are found [25,26] to reproduce the binding energies as well as partial atomic charges of homonuclear as well as heteronuclear simple molecular systems quite well.These DFT based formalisms can also be used for predicting the properties of solids.One such example is the Miedema [27] model for the prediction of heat of formation of simple alloy systems.The ingredients of the present theory are the chemical potentials of each of the atomic and bond sites and also the self and cross components of the hardness matrix corresponding to different spins and sites.The atomic up-and down-spin chemical potentials µ 0 α,i and µ 0 β,i as well as the corresponding hardness parameters η 0 α,i , η 0 β,i and the cross term η 0 αβ,i which are essentially the diagonal elements of the hardness kernels, viz.η αα (i, i), η ββ (i, i) and η αβ (i, i) of the i-th atom can be obtained from spin polarised DFT calculations as has been done earlier [26].The bond site chemical potential as well as the hardness parameters can be approximated [28] by suitable averaging [29] of the corresponding values for the neighbouring atoms.For the off-diagonal (i = j) elements of the hardness kernel η µν (i, j), one can essentially employ the atom-in-molecule hardness matrix concept of Nalewajski [30] generalised for the spin-dependence and model along similar lines [31] following basically an electrostatic analogy details of which have been discussed elsewhere [18].For a nonbonded pair of sites (atom or bond), one considers a Coulomb potential η αβ (i, j) = 1/ R ij , with as the dielectric constant of the electron cloud medium.while for bonded sites, a better modeling [31] involves the Mataga-Nishimoto-Ohno formula as ).Thus, the DFT framework with the simplifications through the associated concepts is shown to lead to a simple equation (Eq.( 31)) for predicting the intra-as well as inter-molecular interaction energies.

Density Functional Theory of Soft Matter: A Mesoscopic Domain
While the quantum DFT discussed in the previous section with electron density as the basic variable provides a simplified route to the estimation of the intermolecular potential, the present prescription which considers DFT in the mesoscopic length scale uses this as input.For the soft matter systems considered here, the collective interaction and correlation effects play important role making the DFT formalism an attractive tool for study of these many-body systems.Some of the applications of DFT in this area include [7] the density distributions at interfaces consisting of liquid-gas, solid-liquid, charged electrode-electrolyte [32] interfaces, oscillating solvation force [33,34], sedimentation [35] and adsorption isotherm in porous medium [36] etc.The swelling of polymer gels on absorbing water, display using liquid crystals, field responsive electro-/magnetorheological fluids etc involve a distribution of the complex fluid particles in external fields.The structure of complex fluids and their mixtures at solid-fluid interfaces is directly relevant to adsorption, wetting, catalysis etc and is thus of both fundamental and technological importance.We therefore consider equilibrium arrangement of the fluid particles, which are the building blocks in the intermediate mesoscopic length scale, at interfaces due to confinement.
The system of interest here is a two-component fluid mixture, with the interparticle interaction (such as hard sphere, Lennard-Jones etc for neutral species and charged hard sphere, DLVO screened Coulomb (for colloids) etc for charged species) u µν (r 12 ) with µ, ν = α, β denoting the two components.The fluid mixture is confined in a cavity (a slit between two planar walls or a spherical or cylindrical pore) which provides the external potential v ν (r) and the resulting inhomogeneous density distribution is ρ ν (r) for the ν-th component.This inhomogeneous fluid mixture is considered to be in equilibrium with the corresponding bulk phases of densities ρ 0 α and ρ 0 β and chemical potentials µ α and µ β for the two components.The DFT for this system is quite analogous to the spin polarised DFT for the many-electron systems considered in the previous section.
The grand potential density functional Ω[ρ α , ρ β ] for this two-component system can be written as where F ex [ρ α , ρ β ] is the excess free energy arising from interparticle interaction and correlation effects.The grand potential, on minimisation, with respect to the component densities leads to the Euler-Lagrange equations, the final forms of which as obtained after equating the component While these approaches are either fully or partially perturbative, nonperturbative WDA based methods can also be obtained by generalising Eqs.(18) or (19) to the case of mixtures.The basic approach here is the same and terms of all orders are approximated by evaluating the excess free energy or the first order DCF at suitable effective weighted densities.In the partially nonperturbative approaches [39], WDA is used to evaluate the short range contribution to the first-order DCF and for the long range contribution, one employs a functional perturbative expansion.A typical application that has been considered involves the structure of electrode-electrolyte interface, with the ionic components modeled as charged hard spheres, and the solvent as a dielectric continuum.For this restricted primitive model [32] the long range contribution leads to an electrostatic potential analogous to that in the many-electron systems, and the short range correlation is evaluated through WDA or functional perturbative method.For walls or cavities of planar, spherical or cylindrical shapes, the density variation is along only one coordinate and hence the numerical calculation becomes much simpler.
An important consequence of the inhomogeneous density distribution at the interfaces is the solvation or structural force which can be directly measured experimentally.The interaction between colloidal particles in a suspension is also essentially an interaction between the overlapping electric double layers.DFT has been used successfully to predict the solvation force for simple liquids and the results are in good agreement with simulation results.We have also calculated the interaction energies between the walls for a neutral liquid and have been able to reproduce the experimental results for a real liquid quite well [34] by choosing a suitable value for the hard sphere diameter.The DFT approach has also been used to calculate the interaction between colloidal particles.There is also interest in the study of the effect of surface induced inhomogeneity in density and the surface forces on the dynamics in condensed phase.
Besides these cases involving density inhomogeneity as a consequence of external potential due to confinement, there are situations where one has spontaneous appearance of density inhomogeneity, for example, in phase transitions.Here thermal density fluctuations leading to density inhomogeneity get stabilised at suitable values of the external parameters such as temperature, density etc. DFT has been quite successful in its prediction of the phase transition from liquid phase to an ordered phase due to the pioneering works of Ramakrishnan and Yussouff [40] who laid the foundation of the DFT of freezing.The theory employs the free energy or the grand potential functionals to locate the physical parameters at which the ordered phase with an inhomogeneous periodic density distribution is energetically more stabilised as compared to the uniform density fluid phase.The original theory has been extended and applied extensively to a large class of soft matter systems including colloidal suspensions [41], liquid crystals and other forms of complex fluids.With a simple modeling of the interfacial density as a decaying oscillatory function in the perpendicular direction and periodic in the parallel plane, DFT has been able to predict [42] the surface energy as a function of the number of interfacial layers in good agreement with simulation results.Application to various other phase transitions, nucleation and many other phenomena has also been found to be quite successful.

Density Functional Theory of Dynamical Phenomena
The applicability of DFT is not restricted to equilibrium situations alone, but is extended to timedependent (TD) situations for both quantum and classical systems.The basic variables involve the TD density ρ(r, t) and TD current density j(r, t).The TD DFT [43] for quantum systems involve the TD one-particle Kohn-Sham like equations of the form where v ef f (r, t) is the effective TD potential dependent on the densities determined by the orbitals.The equivalent quantum hydrodynamic formulation [44] consists of the continuity equation corresponding to the time evolution of the density given by ∂ρ(r, t) ∂t = ∇.j(r,t) (40) and an Euler type equation for the time-evolution of the current density as ∂j(r, t) ∂t = −ρ(r, t)∇ v(r, t) + δF δρ − ∇. jj ρ (41) which differs from the conventional hydrodynamic equation in the additional force of quantum origin which augments the classical force.For classical TD DFT, the continuity equation is the same as Eq.( 40) and the equation for the current density is given by ∂j(r, t) ∂t = −ρ(r; t)∇ v(r, t) + δF δρ − t 0 dτ dr Γ(r, r ; t − τ ) δF δj(r , τ) (42) and involves the dissipation kernel Γ(r, r ; t).The quantum and classical TD DFT thus involve closely analogous frameworks, suitable for bridging the microscopic and mesoscopic length scales.The TD DFT equations have also been extended to spin-polarised situations, and multicomponent fluid mixtures leading to many applications to dynamical phenomena.The interconnections of TD DFT formalism with mode coupling theories [45,46] of diffusion in fluid mixtures have also been demonstrated recently.

Concluding Remarks
This work presents a brief review of DFT as applied to quantum as well as classical systems encompassing the microscopic, mesoscopic and macroscopic length scales, the emphasis being on the unified nature of the framework and its promise for future applications in diverse interdisciplinary areas of research.The present review is however only illustrative and not exhaustive.Materials modeling is an exciting and promising field of research and with the advancement of theoretical and computational techniques it is expected to be more and more challenging and further developments in the field of DFT are likely to play an important role.