## Construction of the exchange matrix

LDA and GGA exchange-correlation functionals typically underestimate the bandgap of a material and may even predict a metal even though the material is an insulator. One way to overcome this issue is the exployment of the LDA+U method. But this approach needs material-specific parameters. Another approach is to mix LDA and GGA with the Hartree-Fock approach which typically overestimates the bandgap. The resulting hybrid functionals require the incorporation of nonlocal (NL), spin $\sigma$ dependent Hartree-Fock-like exchange

in the setup of the Hamiltonian. The sums hereby go over all occupied states $\psi_{\nu,\vec{q}}^\sigma$ within the Brillouin zone (BZ).

In hybrid functionals calculations with Fleur this is finally added together with a scaling factor in terms of a nonlocal exchange matrix $V_{\text{x},\vec{G},\vec{G'}}^{\text{NL},\sigma}(\vec{k})$ to the Hamiltonian matrix. However, for reasons of efficiency the computation of this matrix is first performed in terms of an auxilliary basis of Kohn-Sham orbitals determined in a previous iteration of the SCF loop. It is then brought into the final form by a basis set transformation. In detail, the matrix in terms of the auxilliary basis is

Of course, the size of the auxilliary basis is a cutoff parameter that has to be considered when ensuring parameter convergence.

In the input file the number of Kohn-Sham orbitals used for the auxilliary basis is set in calculationSetup/prodBasis/@bands. Setting it to a certain value also requires setting calculationSetup/cutoffs/@numbands to a slightly larger value.

To construct this matrix another basis set is introduced: The mixed product basis (MPB) which is capable of efficiently representing wavefunction products. It consists of plane waves projected onto the interstitial region

and additional functions with Bloch function character to represent the products in each MT sphere $M_{l,m,p}^{\alpha,\vec{k}}$. These functions are directly constructed by considering the products of all possible radial functions contributing in a product to the respective $(l,m)$ channel. Since the resulting set of functions is typically numerically linearly dependent, they are subsequently used to construct a basis (indexed by $p$) spanning approximately the same space. Finally these functions are used to construct the respective Bloch functions.

Indexed by $I$ overall the mixed product basis is given as $\left\lbrace M_I^\vec{k}(\vec{r}) \right\rbrace = \left\lbrace M_{l,m,p}^{\alpha,\vec{k}}, M_\vec{G}^\vec{k} \right\rbrace$.

With this basis the expression for the exchange matrix becomes

where

is the Coulomb matrix expressed in terms of the mixed product basis.

The mixed product basis comes with several cutoff parameters determining its size and composition. These have to be considered when optimizing the parameter convergence.

The reciprocal plane-wave cutoff for the interstitial basis functions within the MPB is set in /calculationSetup/prodBasis/@gcutm. For the MT part atomSpecies/species/prodBasis/@lcutm defines the angular momentum cutoff for the MPB. (TBC)

A detailed discussion on the construction of the exchange matrix in the context of the PBE0 hybrid functional is available in M. Betzinger, C. Friedrich, and S. Blügel, Hybrid functionals within the all-electron FLAPW method: implementation and applications of PBE0, Phys. Rev. B 81, 195117 (2010).