Introduction

KSSOLV2.0 :cite:`KSSOLV2` is a MATLAB toolbox for performing first-principles density functional theory (DFT) ground- and excited- state calculations for molecules and solids within plane-wave basis sets under periodic boundary conditions.

Ground state calculation

The ground state calculation is based on the Kohn-Sham density functional theory. The main equations to be solved are the Kohn-Sham equations of the form

(1)\[H(\rho) \psi_i({\bf r}) = \varepsilon_i \psi_i({\bf r}),\]

where \(\psi_i\), \(i=1,2,...n_e\), with \(n_e\) being the number of electrons, are orthonormal quasi-electron orbitals, and \(\varepsilon_i\) are the corresponding quasi-electron energies. The function \(\rho({\bf r})\) is the electron density defined to be

\[\rho({\bf r}) = \sum_{i=1}^{n_e} \left| \psi_i({\bf r}) \right|^2.\]

In the so called local density approximation (LDA) or general gradient approximation (GGA), the Kohn-Sham Hamiltonian \(H(\rho)\) has the form

(2)\[H(\rho) = -\nabla^2 + \int \frac{\rho({\bf r}')}{\left| {\bf r}-{\bf r}' \right|} d{\bf r}' + v_{xc}(\rho({\bf r})) + v_{ext}({\bf r}),\]

where the first term in the above equation is associated with the kinectic energy of the electrons, the second term corresponds to the electro-static repulsion among electrons, the third terms is the so called exchange-correlation potential that accounts for other many-body properties of the system, and the last term represents the ionic potential contributed by nuclei. There are a number of expressions for \(v_{xc}\) since the exact form of the exchange-correlation potential is unknown. These expressions also depend on the way \(v_{ext}\) is approximated. By treating inner electrons as part of an ionic core, an approach known as the pseudopotential method, we obtain \(v_{ext}\) that is easier to compute.

The Kohn-Sham equations form a nonlinear eigenvalue problem in which the Kohn-Sham Hamiltonian to be diagonalized is a function of the electron density \(\rho({\bf r})\), which is in turn a function of the eigenfunction \(\psi_i\)’s to be computed.

These equations are the first order necessary condition of the contrained minimization problem:

(3)\[\min_{\langle\psi_i,\psi_j\rangle = \delta_{i,j}} E_{total}(\{\psi_j\}) \equiv \frac{1}{2}\sum_{i=1}^{n_e} \int \left|\nabla \psi_i ({\bf r}) \right|^2 dr + E_{ion} + E_{coul}(\rho) + E_{xc}(\rho),\]

where

\[E_{ion} = \int \rho({\bf r}) v_{ext} ({\bf r}) dr,\]

\[E_{coul} = \frac{1}{2}\int \int \frac{\rho(r)\rho(r')}{\left|r - r'\right|} dr dr'.\]

Therefore, the Kohn-Sham nonlinear eigenvalue problem can also be solved as a constrained optimization problem.