DFT without eigenvectors

Introduction

Several methods have been proposed (Bekas 2005)

Theory

We define a set of orthonormal functions \(\left\{\varsigma(\vec{r})_j\right\}_j\) that span the same basis as the orbitals \(\left\{\varphi(\vec{r})_j\right\}_j\). The are related through a orthogonal transformation \({\overleftrightarrow{R}}\) such that \[\varphi_i(\vec{r}) = \sum_jR_{ij}\varsigma_j(\vec{r}) \ .\] \[\begin{aligned} n(\vec{r}) &= \sum_if_i\varphi_i^*(\vec{r})\varphi_i(\vec{r}) \\ &=\sum_{ijk}f_iR^*_{ij}\varsigma_j^*(\vec{r})R_{ik}\varsigma_k(\vec{r})\\ &=\sum_{jk}\varsigma_j^*(\vec{r})\left(\sum_iR^*_{ij}f_iR_{ik}\right)\varsigma_k(\vec{r})\\ &=\sum_{jk}\varsigma_j^*(\vec{r})F_{jk}\varsigma_k(\vec{r})\end{aligned}\]

\[F_{jk}=\sum_iR^*_{ij}f_iR_{ik} \label{eq:occupationmatrix}\]

\[\begin{aligned} \mathrm{Tr}({\overleftrightarrow{F}})&=\sum_jF_{jj}\\ &=\sum_{ij}R^*_{ij}f_iR_{ij}\\ &=\sum_{j}f_j\\ &=N\end{aligned}\]

\[\begin{aligned} \sum_if_i\epsilon_i &= \sum_if_i\int\mathrm{d}\vec{r}\varphi_i^*(\vec{r})\int\mathrm{d}\vec{r}'\mathrm{H}(\vec{r},\vec{r}')\varphi_i(\vec{r}') \\ &=\sum_{ijk}f_i\int\mathrm{d}\vec{r}R^*_{ij}\varsigma_j^*(\vec{r})\int\mathrm{d}\vec{r}'\mathrm{H}(\vec{r},\vec{r}')R_{ik}\varsigma_k(\vec{r}')\\ &=\sum_{jk}F_{jk}\int\mathrm{d}\vec{r}\varsigma_j^*(\vec{r})\int\mathrm{d}\vec{r}'\mathrm{H}(\vec{r},\vec{r}')\varsigma_k(\vec{r}')\\ &=\sum_{jk}F_{jk}{\cal{E}}_{jk}\\ &=\mathrm{Tr}({\overleftrightarrow{F}}^\intercal {\overleftrightarrow{{\cal{E}}}})\end{aligned}\]

\[{\cal{E}}_{jk}=\int\mathrm{d}\vec{r}\mathrm{d}\vec{r}'\varsigma_j^*(\vec{r})\mathrm{H}(\vec{r},\vec{r}')\varsigma_k(\vec{r}')\]

In temperature-dependent DFT, the occupations are obtained from the Fermi-Dirac distribution \[f_i=\frac1{1+e^{(\epsilon_i-\mu)/kT}} \label{eq:fermidirac}\] where the chemical potential \(\mu\) is obtained from imposing the condition that the occupations must sum to the total number of electrons. This procedure establishes a relation between the set of KS eigenvalues \(\epsilon_i\) and the occupations \(f_i\). For our approach, we must find a relation between the subspace Hamiltonian \({\overleftrightarrow{h}}\) and the occupation matrix \({\overleftrightarrow{F}}\). Formally, we can diagonalize \({\overleftrightarrow{h}}\) to obtain both \({\overleftrightarrow{R}}\) and \(\epsilon_i\), from there we can then obtain \({\overleftrightarrow{F}}\) from eqs. \ref{eq:fermidirac} and \ref{eq:occupationmatrix} as \[F_{jk}=\sum_iR^*_{ij}\frac1{1+e^{(\epsilon_i-\mu)/kT}}R_{ik}\ . \label{eq:fermidiracdiagonal}\] Of course this involves the diagonalization of \({\overleftrightarrow{{\cal{E}}}}\) which is what we wanted to avoid from the beginning. But since \({\overleftrightarrow{R}}\) and \(\epsilon_i\) are the eigenvectors and eigenvalues of \({\overleftrightarrow{{\cal{E}}}}\) we can rewrite eq. \ref{eq:fermidiracdiagonal} as a matrix formula \[{\overleftrightarrow{F}}=\left[I + \exp\left(\frac1{kT}\left({\overleftrightarrow{{\cal{E}}}}-\mu\right)\right)\right]^{-1} \ . \label{eq:fermidiracmatrix}\] with the condition that the value of \(\mu\) must be such that \(\mathrm{Tr}({\overleftrightarrow{F}})=N\).

We have converted the problem from the diagonalization of \({\cal{E}}\) to the evaluation of the function of a matrix for which several methods exists. So the question is now, can we find a method to evaluate \ref{eq:fermidiracmatrix} that is more efficient for large matrices than diagonalization.

While most approaches are focused on sparse matrices (Sidje 2010)

We can use a Padé approximation (Moler 2003) \[\exp\left({\overleftrightarrow{A}}\right)\sim{\overleftrightarrow{R}}_n({\overleftrightarrow{A}})=\left[{\overleftrightarrow{P}}_n({\overleftrightarrow{A}})\right]^{-1}\left[{\overleftrightarrow{Q}}_n({\overleftrightarrow{A}})\right]\] where \[P_n({\overleftrightarrow{A}})=\sum_{j=0}^n\frac{(2n-j)!\,n!}{(2n)!\,j!\,(n-j)!}{\overleftrightarrow{A}}^j\] and \[Q_n({\overleftrightarrow{A}})=\sum_{j=0}^n\left(-1\right)^j\frac{(2n-j)!\,n!}{(2n)!\,j!\,(n-j)!}{\overleftrightarrow{A}}^j\ .\]

\[\frac{\partial{\overleftrightarrow{F}}}{\partial \mu} = \frac1T\exp\left(\frac1{T}\left({\overleftrightarrow{{\cal{E}}}}-\mu\right)\right){\overleftrightarrow{F}}^{-2}\]

(CHECK ORDER, NON-COMMUTATIVE MULTIPLICATION)

References

  1. Constantine Bekas, Yousef Saad, Murilo L. Tiago, James R. Chelikowsky. Com