The Density Functional Theory and Beyond: Example and Applications

Abstract

Density Functional Theory is one of the most widely used methods in quantum calculations of the electronic structure of matter in both condensed matter physics and quantum chemistry. Despite the importance of the density functional theory to find the correlation-exchange energy, but this quantity remains inaccurate. So we have to go beyond DFT to correct this quantity. In this framework, the random phase approximation has gained importance far beyond its initial field of application, condensed matter physics, materials science, and quantum chemistry. RPA is an approach to accurately calculate the electron correlation energy.

1. Introduction

The study of the microscopic properties of a physical system in the condensed matter branch requires the solution of the Schrödinger equation. When the studied system is composed of a large number of interacting atoms, the analytical solution of the Schrödinger equation becomes impossible. However, certain numerical calculation methods provide access to a solution to this fundamental equation for increasingly large systems. The calculation methods, called ab-initio like the density functional theory (DFT), propose to solve the Schrödinger equation without adjustable parameters. The density functional theory quickly established itself as a relatively fast and reliable way to simulate electronic and structural properties for all of the elements of the periodic table ranging from molecules to crystals. In this chapter, we recall the principle of this theory which considers electron density as a fundamental variable and that all physical properties can be expressed as a function of it.

2. Schrödinger equation

It is a fundamental equation to be solved to describe the electronic structure of a system with several nuclei and electrons and for a non-relativistic quantum description of a molecular or crystalline system and which is written:

$$H\mathrm{\Psi }=\left(-\sum _i^n\frac{{\mathrm{\hslash }}^2{\nabla }_i^2}{2m}-\sum _I^N\frac{{\mathrm{\hslash }}^2{\nabla }_I^2}{2M}-\sum _{i,I}\frac{Z_I{e}^2}{\mid \overleftarrow{r_i}-\overleftarrow{R_I}\mid }+\sum _{i<j}\frac{e^2}{\mid \overleftarrow{r_i}-\overleftarrow{r_j}\mid }+\sum _{I<J}\frac{Z_I{Z}_J{e}^2}{\mid \overleftarrow{R_I}-\overleftarrow{R_J}\mid }\right)\mathrm{\Psi },$$E1

where H is the molecular Hamiltonian and $\mathrm{\Psi }$ is the wave function. It is therefore a question of seeking the solutions of this equation. We can write the Hamiltonian in the form:

$$H=T_e+T_n+V_{e-e}+V_{n-n}+V_{n-e}.$$E2

We give the definition for each term:

$T_e=-{\sum }_i^n\frac{{\mathrm{\hslash }}^2{\nabla }_i^2}{2m}$: The kinetic energy of n electrons of mass m.

$T_n=-{\sum }_I^N\frac{{\mathrm{\hslash }}^2{\nabla }^2}{2M}$: The kinetic energy of N nuclei of mass M.

$V_{e-e}={\sum }_{i<j}\frac{e^2}{\mid \overleftarrow{r_i}-\overleftarrow{r_j}\mid }$: The electron–electron repulsive potential energy.

$V_{e-n}=-{\sum }_{i,I}\frac{Z_I{e}^2}{\mid \overleftarrow{r_i}-\overleftarrow{R_I}\mid }$: The attractive potential energy nucleus-electron.

$V_{n-n}={\sum }_{I<J}\frac{Z_I{Z}_J{e}^2}{\mid \overleftarrow{R_I}-\overleftarrow{R_J}\mid }$: The nucleus-nucleus repulsive potential energy.

For a system of N nuclei and n electrons, Schrödinger equation is too complex to be able to be solved analytically. The exact solution of this equation is only possible for the hydrogen atom and hydrogenoid systems. In order to simplify the solution of this equation, Max Born and Robert Oppenheimer [1] have proposed an approximation aiming to simplify it.

3. The Born-Oppenheimer approximation

We consider that we can decouple the movement of electrons from that of nuclei, by considering that their movement of nuclei is much slower than that of electrons: we consider them as fixed in the study of the movement of the electrons of the molecule. The inter-nuclear distances are then treated as parameters. It has an immediate computational consequence, called an adiabatic hypothesis. It is in fact the same approximation and since the Oppenheimer approximation is still used in quantum chemistry, during chemical reactions or molecular vibrations, we can consider according to the classical Born-Oppenheimer approximation that the distribution of electrons (adapts) almost instantaneously, when from the relative motions of nuclei to the resulting Hamiltonian variation. This is due to the lower inertia of the electrons $M=1800m_e$ then the electron wave function can therefore be calculated when we consider that the nuclei are immobile, from where

$$T_n=0;V_{n-n}=\mathit{\text{constant}},$$E3

and so the Hamiltonian becomes

$$H=T_e+V_{e-e}+V_{e-n}+V_{n-n}.$$E4

$$H=H_{\mathit{ele}}+V_{n-n},$$E5

with $H_{\mathit{ele}}$: electronic Hamiltonian which is equal to:

$$H_{\mathit{ele}}=T_e+V_{e-e}+V_{e-n}.$$E6

Therefore the Born Oppenheimer approximation gives us:

$$H=T_e+V_{e-e}+V_{e-n}.$$E7

We use another notation to simplify the calculations

$$H=T+V_{\mathit{ext}}+U.$$E8

$$T=T_e;U=V_{e-e}=V_H;V_{\mathit{ext}}=V_{e-n}.$$E9

The Born-Oppenheimer approximation results in the Eq. (7) which keeps a very complex form: it always involves a wave function with several electrons. This approximation significantly reduces the degree of complexity but also the new wave function of the system depends on N bodies while other additional approximations are required to be able to effectively solve this equation. The remainder of this chapter will deal with approximations allowing to arrive at a solution of this equation within the framework of the density functional theory (DFT) and the random phase approximation (RPA).

4. Density Functional theory (DFT)

Density Functional Theory is one of the most widely used methods for calculating the electronic structure of matter in both condensed matter physics and quantum chemistry. The DFT has become, over the last decades, a theoretical tool which has taken a very important place among the methods used for the description and the analysis of the physical and chemical properties for the complex systems, particularly for the systems containing a large number electrons. DFT is a reformulation of the N-body quantum problem and as the name suggests, it is a theory that only uses electron density as the fundamental function instead of the wave function as is the case in the method by Hartree and Hartree-Fock. The principle within the framework of the DFT is to replace the function of the multielectronic wave with the electronic density as a base quantity for the calculations. The formalism of the DFT is based on the two theorems of P. Hohenberg and W. Kohn [2].

4.1 Hohenberg and Kohn theorems

Hohenberg-Kohn (HK) reformulated the Schrödinger equation no longer in terms of wave functions but employing electron density, which can be defined for an N-electron system by:

$$n=2N\int {\mathit{dr}}_1\int {\mathit{dr}}_2\dots \int {\mathit{dr}}_{n-1}{\mathrm{\Psi }}^{\ast }\left(r_1,r_2,.\dots r_{n-1},r\right)\mathrm{\Psi }\left(r_1,r_2,.\dots r_{n-1},r\right),$$E10

this equation depends only on the three position parameters r = (x, y, z), position vector of a given point in space. This approach is based on two theorems demonstrated by Hohenberg and Kohn.

Theorem 1: For any system of interacting particles in an external potential $V_{\mathit{ext}}\left(\overleftarrow{r}\right)$, the potential $V_{\mathit{ext}}\left(\overleftarrow{r}\right)$ is only determined, except for an additive constant, by the electron density $n_0\left(\overleftarrow{r}\right)$ in its ground state.

The first HK Theorem can be demonstrated very simply by using reasoning by the absurd. Suppose there can be two different external potentials $V_{\mathit{ext}}^{\left(1\right)}$ and $V_{\mathit{ext}}^{\left(2\right)}$ associated with the ground state density $n\left(\overleftarrow{r}\right)$. These two potentials will lead to two different Hamiltonians $H^{\left(1\right)}$ and $H^{\left(2\right)}$ whose wave functions ${\psi }^{\left(1\right)}$ and ${\psi }^{\left(2\right)}$ describing the ground state are different. As described by the ground state of $H^{\left(1\right)}$ we can therefore write that:

$$E^{\left(1\right)}=⟨{\psi }^{\left(1\right)}|H^{\left(1\right)}|{\psi }^{\left(1\right)}⟩<⟨{\psi }^{\left(2\right)}|H^{\left(1\right)}|{\psi }^{\left(2\right)}⟩.$$E11

This strict inequality is valid if the ground state is not degenerate which is supposed in the case of the approach of HK. The last term of the preceding expression can be written:

$$⟨{\psi }^{\left(2\right)}|H^{\left(1\right)}|{\psi }^{\left(2\right)}⟩=⟨{\psi }^{\left(2\right)}|H^{\left(2\right)}|{\psi }^{\left(2\right)}⟩+⟨{\psi }^{\left(2\right)}|H^{\left(1\right)}-H^{\left(2\right)}|{\psi }^{\left(2\right)}⟩,$$E12

$$⟨{\psi }^{\left(2\right)}|H^{\left(1\right)}|{\psi }^{\left(2\right)}⟩=E^{\left(2\right)}+\int \left[V_{\mathit{ext}}^{\left(1\right)}\left(\overleftarrow{r}\right)-V_{\mathit{ext}}^{\left(2\right)}\left(\overleftarrow{r}\right)\right]n_0\left(\overleftarrow{r}\right)d^{3}r,$$E13

$$E^{\left(1\right)}<E^{\left(2\right)}+\int \left[V_{\mathit{ext}}^{\left(1\right)}\left(\overleftarrow{r}\right)-V_{\mathit{ext}}^{\left(2\right)}\left(\overleftarrow{r}\right)\right]n_0{d}^{3}r.$$E14

It will also the same reasoning can be achieved by considering $E^{\left(2\right)}$ instead of $E^{\left(1\right)}$. We then obtain the same equation as before, the symbols $\left(1\right)$ and $\left(2\right)$ being inverted:

$$E^{\left(2\right)}=⟨{\psi }^{\left(2\right)}|H^{\left(2\right)}|{\psi }^{\left(2\right)}⟩<⟨{\psi }^{\left(1\right)}|H^{\left(2\right)}|{\psi }^{\left(1\right)}⟩,$$E15

$$⟨{\psi }^{\left(1\right)}|H^{\left(2\right)}|{\psi }^{\left(1\right)}⟩=⟨{\psi }^{\left(1\right)}|H^{\left(1\right)}|{\psi }^{\left(1\right)}⟩+⟨{\psi }^{\left(1\right)}|H^{\left(2\right)}-H^{\left(1\right)}|{\psi }^{\left(1\right)}⟩,$$E16

$$⟨{\psi }^{\left(1\right)}|H^{\left(2\right)}|{\psi }^{\left(1\right)}⟩=E^{\left(1\right)}+\int \left[V_{\mathit{ext}}^{\left(2\right)}\left(\overleftarrow{r}\right)-V_{\mathit{ext}}^{\left(1\right)}\left(\overleftarrow{r}\right)\right]n_0\left(\overleftarrow{r}\right)d^{3}r.$$E17

$$E^{\left(2\right)}<E^{\left(1\right)}+\int \left[V_{\mathit{ext}}^{\left(2\right)}-V_{\mathit{ext}}^{\left(1\right)}\right]n_0{d}^{3}r,$$E18

we obtain the following contradictory equality:

$$E^{\left(1\right)}+E^{\left(2\right)}<E^{\left(1\right)}+E^{\left(2\right)}.$$E19

The initial hypothesis is therefore false; there cannot exist two external potentials differing by more than one constant leading at the same density of a non-degenerate ground state. This completes the demonstration.

$\Rightarrow$ the external potential of the ground state is a density functional.

Theorem 2: The previous theorem only exposes the possibility of studying the system via density. It only allows knowledge of the density associated with the studied system. The Hohenberg-Kohn variational principle partially answers this problem:

a universal functional for the energy $E\left[n\right]$ can be defined in terms of the density. The exact ground state is the overall minimum value of this functional.

Since the fundamental energy of the system is uniquely determined by its density, then energy can be written as a density functional. By following reasoning similar to that of the first part we show that the minimum of the functional corresponds to the energy of the ground state, indeed, the total energy can be written:

$$E_{\mathit{HK}}\left[n\right]=\int n\left(\overleftarrow{r}\right)V_{\mathit{ext}}\left(\overleftarrow{r}\right)d^{3}r+F_{\mathit{HK}}\left[n\right],$$E20

$F\left[n\right]$ is a universal functional of n(r):

$$F_{\mathit{HK}}\left[n\right]=T\left[n\right]+U\left[n\right].$$E21

And the number of particles:

$$N=\int n\left(\overleftarrow{r}\right)\mathit{dr}.$$E22

Thus, we see that by minimizing the energy of the system with respect to the density we will obtain the energy and the density of the ground state. Despite all the efforts made to evaluate this functional E[n], it is important to note that no exact functional is yet known.

4.2 Ansatz of Kohn-Sham

Since the kinetic energy of a gas of interacting electrons being unknown, in this sense, Walter Kohn and Lu Sham [3] (KS) proposed in 1965 an ansatz which consists in replacing the system of electrons in interaction, impossible to solve analytically, by a problem of independent electrons evolving in an external potential. In the case of a system without interaction, the functional E[n] is reduced to kinetic energy and the interest of the reformulation introduced by Kohn and Sham is that we can now define a monoelectronic Hamiltonian and write the equations monoelectronic Kohn-Sham. According to KS the energy is written in the following form:

$$E_{\mathit{HK}}\left[n\right]=T_s\left[n\right]+\int V_{\mathit{ext}}\left(r\right)n\left(r\right)d^{3}r+E_{\mathit{\text{hartree}}}\left[n\right]+E_{\mathit{xc}}\left[n\right],$$E23

with the functional:

$$F_{\mathit{HK}}\left[n\right]=E_{\mathit{HK}}\left[n\right]-\int V_{\mathit{ext}}\left(\overleftarrow{r}\right)n\left(\overleftarrow{r}\right)d^{3}r.$$E24

$$F_{\mathit{HK}}\left[n\right]=T_s\left[n\right]+E_c\left[n\right]+E_{\mathit{\text{hartree}}}\left[n\right]+E_x\left[n\right]=T_s\left[n\right]+E_{\mathit{\text{hartree}}}\left[n\right]+E_{\mathit{xc}}\left[n\right].$$E25

$T_s\left[n\right]:$ representing the kinetic energy of a fictitious gas of non-interacting electrons but of the same density is given by:

$$T_s\left[n\right]=\sum _i\int \mathit{dr}{\mathrm{\Psi }}_i^{\ast }\left(\overleftarrow{r}\right)\frac{-{\nabla }^2}{2}{\mathrm{\Psi }}_i\left(\overleftarrow{r}\right).$$E26

$$E_{\mathit{\text{hartree}}}\left[n\right]=\frac{e^2}{8{\mathit{\pi \epsilon }}_0}\int \int \frac{n\left(\overleftarrow{r}\right)n\left(\overleftarrow{r^{\prime }}\right)}{\mid \overleftarrow{r}-\overleftarrow{r^{\prime }}\mid }{d}^3{\mathit{rd}}^3{r}^{\prime }.$$E27

$T_s\left[n\right]:$ the kinetic energy without interaction.

$E_{\mathit{xc}}\left[n\right]:$ the exchange-correlation energy.

$E_{\mathit{\text{hartree}}}\left[n\right]:$ the electron–electron potential energy.

$$E_{\mathit{xc}}\left[n\right]=E_{\mathit{HK}}\left[n\right]-\int V_{\mathit{ext}}\left(\overleftarrow{r}\right)n\left(\overleftarrow{r}\right)d^{3}r-T_s\left[n\right]-E_{\mathit{\text{hartree}}}\left[n\right].$$E28

$$E_{\mathit{xc}}\left[n\right]=F\left[n\right]-T_s\left[n\right]-E_{\mathit{\text{hartree}}}\left[n\right].$$E29

Based on the second Hohenberg-Kohn theorem, which shows that the electron density of the ground state corresponds to the minimum of the total energy and on the condition of conservation of the number of particles

$$\mathit{\delta N}\left[n\left(\overleftarrow{r}\right)\right]=\int \mathit{\delta n}\left(\overleftarrow{r}\right)\mathit{dr}=0,$$E30

So we have:

$$\delta \left\{E_{\mathit{HK}}\left[n\right]-\mu \left(\int n\left(\overleftarrow{r}\right)d^{3}r-N\right)\right\}=0.$$E31

$$\frac{\delta E_{\mathit{HK}}\left[n\right]}{\mathit{\delta n}\left(\overleftarrow{r}\right)}=\mu ,$$E32

$$\frac{\delta T_s\left[n\right]}{\mathit{\delta n}\left(\overleftarrow{r}\right)}+v^{\mathit{eff}}\left(\overleftarrow{r}\right)=\mu ,$$E33

$$v^{\mathit{eff}}\left(\overleftarrow{r}\right)=V_{\mathit{ext}}\left(\overleftarrow{r}\right)+\frac{e^2}{8{\mathit{\pi \epsilon }}_0}\int \frac{n\left(\overleftarrow{r^{\prime }}\right)}{\mid \overleftarrow{r}-\overleftarrow{r^{\prime }}\mid }{d}^3{r}^{\prime }+\frac{\delta E_{\mathit{xc}}\left[n\right]}{\mathit{\delta n}\left(\overleftarrow{r}\right)},$$E34

therefore, the kinetic energy without interaction $T_s\left[n\right]$ is determined by:

$$\frac{\delta T_s\left[n\right]}{\mathit{\delta n}\left(\overleftarrow{r}\right)}={\left(3{\pi }^{2}n\right)}^{\frac{5}{3}}\frac{{\mathrm{\hslash }}^2}{2m}=\frac{{\mathrm{\hslash }}^2}{2m}{k}_F^2.$$E35

Finally, the mono-electronic Hamiltonian of Kohn-sham in atomic unit is put in the form:

$$H=\frac{-{\nabla }^2}{2}+v^{\mathit{eff}}.$$E36

The Hamiltonian is iteratively computed, the self-consistency of a loop is reached when the variation of the calculated quantity is lower than the fixed convergence criterion. The wave functions are calculated by a conjugate gradient method (or equivalent). The density is built from the wave functions, convergence is reached when the density is sufficiently close to the density of the previous step. When seeking to optimize the atomic structure of the system, an additional loop is added. With each iteration of this loop, the atomic positions are changed. It is said that the system is minimized when the forces are lower than the convergence criterion on the amplitude of the forces.

4.3 Expression of the exchange and correlation term

As described above, DFT is at the stage of Kohn-Sham equations, a perfectly correct theory insofar as the electron density which minimizes the total energy is exactly the density of the system of N interacting electrons. However, DFT remains inapplicable because the exchange-correlation potential remains unknown. It is therefore necessary to approximate this exchange-correlation potential. Two types of approximations exist the local density approximation or LDA and the generalized gradient approximation or GGA as well as the derived methods which are based on a non-local approach.

4.3.1 Local density approximation (LDA)

In only one model case, that of the uniform gas of electrons (corresponding quite well to the electrons of the conduction band of a metal), we know the exact expressions or with an excellent approximation of the terms of exchange and correlation respectively. In this LDA (Local Density Approximation), the electron density is assumed to be locally uniform and the exchange-correlation functional is of the form:

$$E_{\mathit{xc}}^{\mathit{LDA}}=\int n\left(\overleftarrow{r}\right){\epsilon }_{\mathit{xc}}^{\mathit{\text{homo}}}\left(n\left(\overleftarrow{r}\right)\right)\mathit{dr}.$$E37

$${\epsilon }_{\mathit{xc}}^{\mathit{\text{homo}}}\left(n\left(\overleftarrow{r}\right)\right)={\epsilon }_x\left[n\right]+{\epsilon }_c\left[n\right].$$E38

The function of ${\epsilon }_{\mathit{xc}}^{\mathit{\text{homo}}}\left(n\left(\overleftarrow{r}\right)\right)$ is determined from a quantum computation parameterization for a constant electron density $n\left(\overleftarrow{r}\right)=n$;

$${\epsilon }_{\mathit{xc}}^{\mathit{\text{homo}}}\left(n\left(\overleftarrow{r}\right)\right)=\frac{3}{4}{\left(\frac{3n\left(r\right)}{\pi }\right)}^{\frac{1}{3}}.$$E39

Ceperley-Alder [4] numerically determined the contribution of the correlations. The search for analytical functions that come as close as possible to these results leads to the development of various functionalities with varying degrees of success. In general, the LDA approximation gives good results in describing the structural properties, i.e. it allows to determine the energy variations with the crystalline structure although it overestimates the cohesion energy, also concerning the mesh parameter for the majority of solids and good values of elastic constants like the isotropic modulus of compressibility. But this model remains insufficient in in-homogeneous systems.

4.3.2 Generalized gradient approximation (GGA)

To overcome the shortcomings of the LDA method, the generalized gradient approximation considers exchange-correlation functions depending not only on the density at each point, but also on its gradient [5], of the general form.

$$E_{\mathit{xc}}^{\mathit{GGA}}\left[n_{\alpha },n_{\beta }\right]=\int n\left(\overleftarrow{r}\right){\epsilon }_{\mathit{xc}}\left[n_{\alpha },n_{\beta },{\nabla }_{n_{\alpha }\left(r\right)}{\nabla }_{n_{\beta }\left(\overleftarrow{r}\right)}\right]\mathit{dr},$$E40

$\alpha$ and $\beta$ are spins, in this case again, a large number of expressions have been proposed for this factor ${\epsilon }_{\mathit{xc}}$ leading to so many functionals. In general, the GGA improves compared to the LDA a certain number of properties such as the total energy or the energy of cohesion, but does not lead to a precise description of all the properties of a semiconductor material namely its electronic properties.

4.3.3 Functional hybrid HSE

The functions of DFT have been proved to be quite useful in explaining a wide range of molecular characteristics. The long-term nature of the exchange interaction, and the resulting huge processing needs, are a key disadvantage for periodic systems. This is especially true for metallic systems that necessitate BZ sampling. A new hybrid functionality, recently proposed by Heyd et al. [6], addresses this problem by separating the description of the exchange and the interaction into a short and long part. The expression of the exchange-correlation energy in HSE03 is given by:

$$E_{\mathit{xc}}^{\mathit{HSE}03}=\frac{1}{4}{E}_x^{\mathit{sr},\mu }+\frac{3}{4}{E}_x^{\mathit{PBE},\mathit{sr},\mu }+E_x^{\mathit{PBE},\mathit{lr},\mu }+E_c^{\mathit{PBE}}.$$E41

As can be seen from the Eq. (41) only the exchange component of the electron–electron interaction is split into a short and long (lr) range (sr) part. The full electron correlation is represented by the standard correlation portion of the density of the GGA functional. Note that the term hybrid refers to the combined use of the exact exchange energy of the Hartree-Fock model and the exchange-correlation energy at the DFT level. The construction of hybrid functionals has been a good advancement in the field of exchange-correlation energy processing by allowing an explicit incorporation of the nonlocal character through the use of the exact term of exchange energy.

5. Random phase approximation

Despite the DFT is relevance in determining the exchange-correlation energy, it is still insufficient to characterize elastic characteristics. To rectify this amount, we must go beyond DFT. This is due to the random phase approximation as a method for calculating the electronic correlation energy accurately. RPA appears in the 1950s’ [7, 8, 9, 10] as a method of solving the N-body problem and arises from the desire to describe better (i.e. better than in a mean-field approximation) the physics of uniform electron gas, where the correlation between the positions of long-range electrons is important. In fact, collective oscillations (called plasma oscillations) are observed in an electron gas, which is the direct consequence of the long-range correlation between the electrons. Bohm and Pines, who introduced RPA, propose to place these collective oscillations at the center of solving the N-body problem, hoping that a good description of one will provide a good understanding of the other. RPA has been used with some success in the literature to describe systems containing Van der Waals interactions and in particular involving [11] scattering forces, which are known to be difficult to process. RPA introduced within the framework of DFT via the fluctuation-dissipation theorem with adiabatic connection (ACFDT).

6. Adiabatic-connection fluctuation-dissipation theory

The adiabatic connection fluctuation-dissipation (AC-FDT) [12] technique will be explained in order to discover the exact exchange-correlation energy in RPA. It will serve as the starting point for introducing the random phase approximation because it provides a general formulation for the exact correlation energy.

6.1 Adiabatic-connection (AC)

The adiabatic connection (AC) is a way to express the exact exchange-correlation energy function. The central idea in this approach is to build an interpolation Hamiltonian, which connects a Hamiltonian of an independent particle (reference Hamiltonian) ${\widehat{H}}_0=\widehat{H}\left(\lambda =0\right)$ and the physical Hamiltonians (multi-body Hamiltonian) $\widehat{H}=\widehat{H}\left(\lambda =1\right)$, with $\lambda$ being a connection parameter. The AC technique can be used to derive the total energy of the ground state of a Hamiltonian of multiple interacting bodies, in which a continuous set of Hamiltonians dependent on the coupling force ($\lambda$) is introduced by:

$$\widehat{H}\left(\lambda \right)={\widehat{H}}_0+\lambda {\widehat{H}}_1\left(\lambda \right)=\widehat{T}+V_{\mathit{ne}}+\widehat{V}\left(\lambda \right)+\lambda W_{\mathit{ee}}=\sum _{i=1}^N\left[\frac{-1}{2}{\nabla }_i^2+v_{\lambda }^{\mathit{ext}\left(i\right)}\right]+\sum _{i>j=1}^N\frac{\lambda }{\mid r_i-r_j\mid }.$$E42

With N being the number of electrons, $v_{\lambda }^{\mathit{ext}}$ is an external potential with $v_{\lambda =1}^{\mathit{ext}}\left(r\right)=v^{\mathit{ext}}\left(r\right)$, being the external physical potential of the fully interactive system. Additionally, $v_{\lambda }^{\mathit{ext}}$ can be spatially non-local for $\lambda \ne 1$. Following that, the reference Hamiltonian, or the Hamiltonian for an independent particle specified by the Eq. (42) for $\lambda =0$, is of the mean field type, or is known in English as (Mean-field (MF)), i.e., a simple synthesis on a single-particle Hamiltonian:

$${\widehat{H}}_0=\sum _{i=1}^N\left[\frac{-1}{2}{\nabla }_i^2+v_{\lambda =0}^{\mathit{ext}}\left(r_i\right)\right]=\sum _{i=1}^N\left[\frac{-1}{2}{\nabla }_i^2+v^{\mathit{ext}}\left(r_i\right)+v^{\mathit{MF}}\left(r_i\right)\right].$$E43

With $v^{\mathit{MF}}$ is an average field potential resulting from the electron–electron interaction. It can be the Hartree-Fock (HF) potential ($v^{\mathit{HF}}$) or the Hartree-plus correlation-exchange potential ($v^{\mathit{Hxc}}$) in the DFT. According to the two Eqs. (42) and (43), the perturbative Hamiltonian becomes:

$$\begin{array}{ll}{\widehat{H}}_1\left(\lambda \right)=& \sum _{i>j=1}^N\frac{1}{\mid r_i-r_j\mid }+\frac{1}{\lambda }\sum _{i=1}^N\left[v_{\lambda }^{\mathit{ext}}\left(r_i\right)-v_{\lambda =0}^{\mathit{ext}}\left(r_i\right)\right]\\ & =\sum _{i>j=1}^N\frac{1}{\mid r_i-r_j\mid }+\frac{1}{\lambda }\sum _{i=1}^N\left[v_{\lambda }^{\mathit{ext}}\left(r_i\right)-v^{\mathit{ext}}\left(r_i\right)-v^{\mathit{MF}}\left(r_i\right)\right].\end{array}$$E44

In the construction of the total energy, the ground state wave function $\mid {\mathrm{\Psi }}_{\lambda }⟩$ is introduced for the system $\lambda$, such that

$$H\left(\lambda \right)\mid {\mathrm{\Psi }}_{\lambda }⟩=E\left(\lambda \right)\mid {\mathrm{\Psi }}_{\lambda }⟩.$$E45

Adopt the normalization condition, $⟨{\mathrm{\Psi }}_{\lambda }|{\mathrm{\Psi }}_{\lambda }⟩=1$, the total interacting energy of the ground state can then be obtained using the theorem of Hellmann-Feynman [13]

$$E\left(\lambda =1\right)=E_0+{\int }_0^1\mathit{d\lambda }\times ⟨{\mathrm{\Psi }}_{\lambda }|\left({\widehat{H}}_1\left(\lambda \right)+\lambda \frac{d{\widehat{H}}_1\left(\lambda \right)}{\mathit{d\lambda }}\right)|{\mathrm{\Psi }}_{\lambda }⟩,$$E46

The energy of order zero is E${}_0$. It should be noted that the adiabatic connecting path chosen in Eq. (46) is not unique. In DFT, the path is chosen so that the electron density remains constant throughout the journey. This suggests a $\lambda$ -dependency that ${\widehat{H}}_1\left(\lambda \right)$ is not aware of.

6.2 The random phase approximation in the framework of adiabatic-connection fluctuation-dissipation theory

We will quickly discuss the concept of RPA in the context of DFT, which has served as the foundation for current RPA computations. The total ground state energy for an interacting N electron system is a (implicit) function of the electron density n(r) in the Kohn-Sham approximation (KS-DFT) and can be divided into four terms:

$$E\left[n\left(r\right)\right]=T_s\left[\psi \left(r\right)\right]+E_{\mathit{ext}}\left[n\left(r\right)\right]+E_H\left[n\left(r\right)\right]+E_{\mathit{xc}}\left[{\psi }_i\left(r\right)\right].$$E47

In the KS framework, the electron density is obtained from the single particle ${\psi }_i\left(r\right)$ orbitals via $n\left(r\right)={\sum }_i^{\mathit{occ}}{\left|{\psi }_i\left(r\right)\right|}^2$. Among the four terms of the Eq. (47) only $E_{\mathit{ext}}\left[n\left(r\right)\right]$ and $E_H\left[n\left(r\right)\right]$ are explicit functions of n(r). $T_s$ is treated exactly in KS-DFT in terms of single particle ${\psi }_i\left(r\right)$ orbitals which are themselves functional of n(r). The unknown correlation-exchange (XC) energy term, which is approximated as an explicit functionality of n(r) (and its local gradients) in conventional functional functions (LDA and GGA) and as a function of ${\psi }_i\left(r\right)$ in more advanced functions, contains the complete complexity of many bodies (hybrid density functions, RPA, etc.). In DFT, several existing approximations of $E_{\mathit{xc}}$ can be categorized using a hierarchical approach called Jacob’s scale [14]. But what if we wish to improve the accuracy of $E_{x}c$ in a larger number of systems? To that purpose, starting with the technically accurate manner of generating $E_{x}c$ using the AC technique mentioned above is instructive. As previously stated, the AC path is used in KS-DFT in order to maintain the correct electron density. Reducing the Eq. (46) for the total energy of the exact ground state E = E ($\lambda$ = 1) to:

$$\begin{array}{ll}E& =E_0+{\int }_0^1\mathit{d\lambda }⟨{\mathrm{\Psi }}_{\lambda }|\frac{1}{2}\sum _{i\ne j=1}^N\frac{1}{\mid r_i-r_j\mid }|{\mathrm{\Psi }}_{\lambda }⟩+{\int }_0^1\mathit{d\lambda }⟨{\mathrm{\Psi }}_{\lambda }|\sum _{i=1}^N\frac{d}{\mathit{d\lambda }}{v}_{\lambda }^{\mathit{ext}}\left(r_i\right)|{\mathrm{\Psi }}_{\lambda }⟩\\ & =E_0+\frac{1}{2}{\int }_0^1\mathit{d\lambda }\int \int \mathit{\text{drdr}}\text{'}\times ⟨{\mathrm{\Psi }}_{\lambda }|\frac{\widehat{n}\left(r\right)[\widehat{n}\left(r\text{'}\right)-\delta \left(r-r\text{'}\right)}{\mid r-r\text{'}\mid }]|{\mathrm{\Psi }}_{\lambda }⟩+\int \mathit{drn}\left(r\right)\left[v_{\lambda }^{\mathit{ext}}\left(r\right)-v_{\lambda =0}^{\mathit{ext}}\left(r\right)\right].\end{array}$$E48

$$\widehat{n}\left(r\right)=\sum _{i=1}^N\delta \left(r-r_i\right),$$E49

$\widehat{n}\left(r\right)$ is the electron density operator and $n\left(r\right)=⟨{\mathrm{\Psi }}_{\lambda }|\widehat{n}\left(r\right)|{\mathrm{\Psi }}_{\lambda }⟩$, for any $0⩽\lambda ⩽1$. For the reference state $\mid {\mathrm{\Psi }}_0⟩$ of KS (given by the slater determinant of orbitals ${\mathrm{\Psi }}_i\left(r\right)$ occupied by a single particle, we get

$$E_0=⟨{\mathrm{\Psi }}_0|\sum _{i=1}^N\left[-\frac{1}{2}{\nabla }^2+v_{\lambda =0}^{\mathit{ext}}\left(r_i\right)\right]|{\mathrm{\Psi }}_0⟩=T_s\left[{\mathrm{\Psi }}_i\left(r\right)\right]+\int \mathit{drn}\left(r\right)v_{\lambda =0}^{\mathit{ext}}\left(r\right).$$E50

$$E=T_s\left[{\mathrm{\Psi }}_i\left(r\right)\right]+\int \mathit{drn}\left(r\right)v_{\lambda =1}^{\mathit{ext}}+\frac{1}{2}{\int }_0^1\mathit{d\lambda }\int \int \mathit{\text{drdr}}\text{'}⟨{\mathrm{\Psi }}_{\lambda }\mid \frac{\widehat{n}\left(r\right)\left[\widehat{n}\left(r\text{'}\right)-\delta \left(r-r\text{'}\right)\right]}{\mid r-r\text{'}\mid }\mid {\mathrm{\Psi }}_{\lambda }.$$E51

From the Eqs. (47) and (51), we obtained:

$$E_H\left[n\left(r\right)\right]=\frac{1}{2}\int \mathit{\text{drdr}}\text{'}\frac{n\left(r\right)n\left(r\text{'}\right)}{\mid r-r\text{'}\mid }.$$E52

$$E_{\mathit{ext}}\left[n\left(r\right)\right]=\int \mathit{drn}\left(r\right)v_{\lambda =1}^{\mathit{ext}}\left(r\right).$$E53

We get the formally exact correlation-exchange energy expression XC;

$$E_{\mathit{xc}}=\frac{1}{2}\int \mathit{d\lambda }\int \int \mathit{\text{drdr}}\text{'}\frac{n_{\mathit{xc}}^{\lambda }\left(r,r\text{'}\right)n\left(r\right)}{\mid r-r\text{'}\mid },$$E54

with $n_{\mathit{xc}}^{\lambda }\left(r,r\text{'}\right)$ is defined by

$$n_{\mathit{xc}}^{\lambda }\left(r,r\text{'}\right)=\frac{⟨{\mathrm{\Psi }}_{\lambda }|\delta \widehat{n}\left(r\right)\delta \widehat{n}\left(r\text{'}\right)|{\mathrm{\Psi }}_{\lambda }⟩}{n\left(r\right)}-\delta \left(r-r\text{'}\right).$$E55

The mathematical expression for the so-called XC-hole is 55, with $\delta \widehat{n}=\widehat{n}\left(r\right)-n\left(r\right)$ denoting the fluctuation of the density operator $\widehat{n}\left(r\right)$ around its expectation value n(r). The hole (XC) is also related to the density-density correlation function, as shown by the Eq. (55). It illustrates how the presence of an electron at point r reduces the density of all other electrons at point $r\text{'}$ in physical terms. The temperature fluctuation-dissipation (FDT) theorem is used to relate the density-density correlations (fluctuations) in the Eq. (55) to the response (dissipation) features of the system in the second step. In statistical physics, FDT is a powerful approach. It shows that the reaction of a system in thermodynamic equilibrium to a tiny external disturbance is the same as the response to spontaneous internal fluctuations in the absence of disturbance [15]. FDT is applicable to both thermal and quantum mechanical fluctuations and shows itself in a variety of physical phenomena. A good example of the latter is the dielectric formulation of the many-body problem by Nozières and Pines [16]. The FDT at zero temperature performed at [16] is relevant in this situation.

$$⟨{\mathrm{\Psi }}_{\lambda }|\delta \widehat{n}\left(r\right)\delta \widehat{n}\left(r\text{'}\right)|{\mathrm{\Psi }}_{\lambda }⟩=-\frac{1}{\pi }{\int }_0^{\infty }\mathit{\text{dωIm}}{\chi }^{\lambda }\left(r,r\text{'},\omega \right),$$E56

with ${\chi }^{\lambda }\left(r,r\text{'},\omega \right)$, is the linear density-response function of the system. Using the Eqs. (54) and (55) and $v\left(r,r\text{'}\right)=\frac{1}{\mid \mathit{rr}\text{'}\mid }$, we arrive at the renamed ACFD expression for XC energy in DFT

$$\begin{array}{l}{E}_{\mathit{xc}}=\frac{1}{2}{\int }_0^1\mathit{d\lambda }\int \int \mathit{\text{drdr}}\text{'}v\left(r,r\text{'}\right)\times \left[-\frac{1}{\pi }{\int }_0^{\infty }\mathit{\text{dωIm}}{\chi }^{\lambda }\left(r,r\text{'},\omega \right)-\delta \left(r-r\text{'}\right)n\left(r\right)\right]=\frac{1}{2\pi }{\int }_0^1\mathit{d\lambda }\int \int \mathit{\text{drdr}}\text{'}\\ \times v\left(r,r\text{'}\right)\times \left[-\frac{1}{\pi }{\int }_0^{\infty }d{\omega }^{\lambda }\left(r,r\text{'},\mathit{i\omega }\right)-\delta \left(r-r\text{'}\right)n\left(r\right)\right].\end{array}$$E57

The analytical structure of ${\chi }^{\lambda }\left(r,r\text{'},\omega \right)$ and the fact that it becomes real on the imaginary axis are the reasons why the above frequency integration can be conducted along the imaginary axis. The problem of computing the energy XC on one of the response functions of a succession of fictional systems along the path AC is transformed by the expression ACFD in the Eq. (57), which must also be tackled in practice. RPA is a particularly basic approximation of the response function in this context:

$${\chi }_{\mathit{RPA}}^{\lambda }\left(r,r\text{'},\mathit{i\omega }\right)={\chi }^0\left(r,r\text{'},\mathit{i\omega }\right)+\int {\mathit{dr}}_1{\mathit{dr}}_2{\chi }^0\left(r,r_1,\mathit{i\omega }\right)\times \mathit{\lambda v}\left(r_1-r_2\right){\chi }_{\mathit{RPA}}^{\lambda }\left(r_2,r\text{'},\omega \right),$$E58

${\chi }^0\left(r,r_1,\mathit{i\omega }\right)$, is the response function of independent particles of KS of the reference system $\lambda =0$ and is known explicitly in terms of orbitals ${\psi }_i\left(r\right)$ single particle (KS), orbital energies ${\epsilon }_i$ and occupancy factors $f_i$:

$${\chi }^0\left(r,r\text{'},\mathit{i\omega }\right)=\sum _{\mathit{ij}}\frac{\left(f_i-f_j\right){\psi }_i^{\ast }\left(r\right){\psi }_j\left(r\right){\psi }_j^{\ast }\left(r\text{'}\right){\psi }_i\left(r\text{'}\right)}{{\epsilon }_i-{\epsilon }_j-\mathit{i\omega }}.$$E59

From the Eqs. (57) and (58), the energy XC in RPA can be split into an exchange-exact (EX) and the correlation term RPA:

$$E_{\mathit{xc}}^{\mathit{RPA}}=E_x^{\mathit{EX}}+E_{\mathit{RPA}}^c.$$E60

$$E_x^{\mathit{EX}}=-\sum _{\mathit{ij}}{f}_i{f}_j\int \int \mathit{\text{drdr}}\text{'}{\psi }_i^{\ast }\left(r\right){\psi }_j\left(r\right)v\left(r,r\text{'}\right){\psi }_j^{\ast }\left(r\text{'}\right){\psi }_i\left(r\text{'}\right).$$E61

$$E_c^{\mathit{RPA}}=\frac{1}{2\pi }{\int }_0^{\infty }\mathit{\text{dωTr}}\left[\mathit{ln}\left(1-{\chi }^0\left(\mathit{i\omega }\right)v\right)+{\chi }^0\left(\mathit{i\omega }\right)v\right].$$E62

7. Approximation of pseudo-potentials

The goal is to study the ground state of a system made up of nuclei, core electrons and valence electrons. The heart electrons are often closely linked to nuclei, they are considered (frozen). This approximation makes it possible to develop the valence wave functions on a reduced number of plane waves having a kinetic energy lower than the energy of the cut-off (E ${}_{\mathit{cut}}$${\mathrm{E}}_{\mathit{cut}}⩾\frac{{\mathrm{\hslash }}^2}{2m}{\left|K+G\right|}^2$), which allows correct treatment of the problem depends on the pseudo-potential used and the system studied. It consists in replacing the ionic potential $V_{\mathit{el},\mathit{nu}}$ by a pseudopotential $V^{\mathit{ps}}$ (see Figure 1) which acts on a set of wave pseudo-functions instead and places true wave functions and having the same eigenstates in the atomic Schrdinger equation. This idea has been developing since the end of the 1950s. This potential is constructed so as to reproduce the scattering properties for the true valence wave functions, while ensuring that the pseudo-wave function does not have a node in the core region defined by a cutoff radius $r_c$ which is optimized for each orbital. Beyond the core region, the pseudopotential is reduced to the ionic potential so that the pseudo-wave function is equal to the true wave function. The use of a pseudo-potential reduces on the one hand the number of electrons considered in the problem by taking into account only the valence electrons and on the other hand makes it possible to restrict the base of plane waves for the electrons of valence by eliminating most of the oscillations of the wave functions in the heart region.

Figure 1.

Schematic illustration of all-electron potential (solid lines) and pseudo-potential (broken lines) and their corresponding wave functions.

8. Projection-augmented plane wave method (PAW)

We are always looking for the precision of the computation and thus to minimize the reasonable computation time as much as possible, we advise to use the least hard potentials possible, that is to say requiring the fewest plane waves to ensure convergence. The method of plane waves augmented by projection (PAW) (Projector Augmented Waves) [17] best meets this requirement, which explains their use in this thesis. The projection augmented plane wave method is an approach developed by P. Blochl, which models heart states from wave functions for an isolated atom; it assumes that these states are not modified when the atom is placed in a compound (approximation of frozen hearts). The PAW approach consists of using two kinds of basic functions: one uses partial atomic wave functions inside cores and the other uses a set of functions. Wave planes in the interstitial region. There is then a separation of space into two parts. The cause of this separation is the behavior of the effective potential seen by the valence shell electrons: in the interstitial zone, the total electron density is low. The wave functions oscillate rapidly which makes it possible to give strong variations in the amplitude of the total electron density. Because of these variations, the Fourier decomposition of the potential created by this charge distribution has components large wave vectors, which disadvantages the treatment in plane waves. We are therefore led to write in a different way the wave functions of the valence electrons in these two regions of space:

• In the interstitial region, plane waves are used.

• Inside the hearts, a partial wave decomposition solutions of the Schrödinger equation for the isolated atom which are much better adapted to the potential which prevails there.

9. Application

In spite of the significance of diamond and silicon, there’s still a need for solid hypothetical and exploratory information on the elastic constants of these materials, in specific on the versatile constants of the third-order. This lack of experimental data limits the capacity of researchers to create modern materials with a focus on mechanical reactions. Besides, this lack triggered interest in other theoretical calculations. M. Barhoumi et al. [18] have proposed to calculate the elastic properties with different approximations of DFT and beyond with ACFDT in RPA, since the RPA has significant advantages, especially for those interested in functional density theory. It correctly describes the dispersion and van der Waals interactions. In this direction, they have found that the results obtained with RPA are in good agreement with the previous published [19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33]. Also, it should be noted that from the calculated elastic constants, other structural properties such as elastic modulus, shear modulus, Young’s modulus and Poisson’s ratio can be derived.

10. Conclusion

In this chapter, we have introduced the general method of calculating the ground state energy of a crystalline solid by application of DFT. We have just described how it is possible to determine the energy of the ground state of a solid by studying a fictitious system of independent particles giving rise to the same density as the real electronic system. On the other hand, we have highlighted the approximations necessary to be able to apply this theory. Despite the importance of the DFT to find the exchange-correlation energy, but this quantity remains inaccurate to describe the elastic properties. So we have to go beyond DFT with RPA to correct this handicap. In this direction, we have shown that RPA is a good description of electronic correlation energy.