# 11. Hybrid functionals

In this section, we introduce hybrid functional calculations in RESCU.
Hybrid functional calculations are so named because the exchange
functional is a weighted sum of a semi-local exchange functional (e.g.
PBE) and the exact exchange functional (employed in the Hartree-Fock
method). Such calculations are more computationally demanding because
the exact exchange functional gives rise to a non-local potential in the
Kohn-Sham equation, but they are also more accurate than semilocal DFT
calculations for a number of applications. Moreover, the non-locality of
the functional brings about some infinity-related issues in periodic
systems. RESCU implements the NOPAIN method
[C18, CCMG18] to deal with both problems. Simply put,
the NOPAIN method performs an *analytical* low-rank decomposition of the
Coulomb repulsion integral, which is the atomic orbital representation
of the Coulomb interaction kernel. Most importantly, basis functions
retain their locality during the decomposition and the coefficients are
obtained from the planewave representation of the Coulomb operator.

## 11.1. Example: diamond band gap

Certain hybrid functionals famously correct band gap in the Kohn-Sham
spectrum. In this example, we show how to calculate the DOS, and hence
band gap, of the diamond crystal. First, we perform a self-consistent
calculation to obtain the ground state density and properties. Please
create the input file `c_lcao_scf_hse`

as follows

```
LCAO.status = true % required for hybrid
info.calculationType = 'self-consistent'
info.savepath = 'results/c_lcao_scf_hse'
atom.element = [1 1]
atom.fracxyz = [0 0 0;0.25 0.25 0.25]
domain.latvec = 6.74*[.5 .5 0;.5 0 .5;0 .5 .5]
domain.lowres = 0.3
element(1).species = 'C'
element(1).path = './C_OV_LDA.mat'
functional.libxc = true % required for hybrid
functional.list = {'XC_HYB_GGA_XC_HSE06'}
kpoint.gridn = [3 3 3] % must be gamma centered
Exx.Gcutoff = 4.0
```

We comment on the purpose of each keyword in a hybrid calculation context.

`LCAO.status`

must be set to true. As mentioned in the section introduction, the NOPAIN method represents the Coulomb kernel using atomic orbitals, and hence the calculation ought to be done in an NAO basis;`element.path`

uses a pseudopotential drawn from the LDA-Troullier-Martins database. For certain atoms (e.g. transition metals, heavy atoms) hybrid calculations may require a pseudopotential including more electronic shells in the reference valence configuration and a more extensive atomic orbital basis. The standard pseudopotential file is sufficient for our simple example;`functional.libxc`

must be set to true. RESCU draws the semi-local part of hybrid functionals from the LibXC library;`functional.list`

includes a single functional (2006-version of the HSE functional) combining the semi-local exchange and correlation functionals in this example. At the moment, RESCU implements only the PBE0 and HSE functionals;`kpoint.gridn`

defines the grid sampling reciprocal space. This grid must be \(\Gamma\)-centered, so set`kpoint.shift = 0.5*[1 1 1]`

for even grids. In HSE calculations, the non-local exchange potential couples the various k-points such that the workload and memory requirements quickly grow with increasing k-sampling. Moreover, RESCU will ignore any symmetry when performing hybrid calculations, which means that the full Brillouin zone is resolved;`Exx.Gcutoff`

is the planewave cutoff, in Hartree, used in the Coulomb kernel decomposition. The higher the more accurate, but also the more computationally expensive. In a research project, convergence with respect to this parameter is crucial. The default value is 4, which is sufficient for several simple materials.

You may execute RESCU as follows

```
rescu -i c_lcao_scf_hse
```

The standard output will ensue with the addition of something like

```
Initializing YqG-matrix 68.365
Initializing Hartree-Fock potential 8.282
```

The first operation corresponds to the NOPAIN decomposition mention previously. In small systems, its cost may dominate the entire calculation, but this is less and less the case as the system size increases. The second operation is the calculation of the exchange exchange potential matrix, which is required to build the Hamiltonian in the first iteration.

Next, we write the DOS input file `c_lcao_dos_hse`

as follows

```
info.calculationType = 'dos'
info.savepath = 'results/c_lcao_dos_hse'
info.outfile = 'fc_lcao_dos_hse.out'
rho.in = './results/c_lcao_scf_hse'
kpoint.gridn = [10 10 10]
kpoint.sampling = 'tetrahedron'
dos.range = [-20,10]
units.energy = 'ev'
```

This is a standard DOS input file. We employ the tetrahedron method and set the range from -20 eV to +10 eV around the Fermi level. You can calculate the DOS and plot it typing

```
rescu -i c_lcao_dos_hse
rescu -p results/c_lcao_dos_hse
```

The band gap can be extracted as follows

```
load('results/c_lcao_dos_hse','energy')
nrg = energy.ksnrg;
bg = (min(nrg(5,:))-max(nrg(4,:)))*27.211
```

In the last line, we subtract the minimum of the fifth band (conduction band minimum) from the maximum of the fourth band (valence band maximum) and multiply by 27.211 to convert the units from Hartree to eV. You should get a value close to 6.2 eV, which is 0.7 eV above the experimental value of 5.5 eV. Following Ref. [C18], increasing the k-sampling to [7,7,7], the planewave cutoff to 6 Ha and using an optimized DZP basis brings the band gap in line with the experimental value. If you repeat the calculation using PBE, keeping all other parameters the save, you should get a value close to 4.7 eV.

Chen Ying-Chih, Chen Jing-Zhe, Vincent Michaud-Rioux, and Hong Guo. Efficient evaluation of nonlocal operators in density functional theory. Physical Review B, 97(7), 075139, 2018.

Chen Ying-Chih. NOPAIN: A Method for Efficient Evaluation of Quantum Nonlocal Operators with Applications to Solids. PhD thesis. McGill University, 2018.