Charge Analysis: Ver. 1.0

Taisuke Ozaki, RCIS, JAIST

Mulliken population

The Kohn-Sham (KS) Bloch functions $\psi_{\mu}$ are expanded in a form of linear combination of pseudo-atomic basis functions (LCPAO) $\phi_{i\alpha}$ centered on site $\tau_i$ by
    $\displaystyle \psi_{\sigma\mu}^{(\bf k)}({\bf r})
...ma\mu,i\alpha}^{(\bf k)}
\phi_{i\alpha}({\bf r}-{\bf\tau}_{i}-{\bf R}_{\rm n}),$ (1)

where $c$ and $\phi$ are an expansion coefficient and pseudo-atomic function, ${\bf R}_{\rm n}$ a lattice vector, $i$ a site index, $\sigma$ ($\uparrow$ or $\downarrow$) spin index, $\alpha\equiv (plm)$ an organized orbital index with a multiplicity index $p$, an angular momentum quantum number $l$, a magnetic quantum number $m$, and $N$ the number of repeated cells. The charge density operator ${\hat n}_{\sigma}$ for the spin index $\sigma$ is given by
    $\displaystyle {\hat n}_{\sigma}
\frac{1}{V_{\rm B}}
\int_{{\rm B}}dk^{3}
...vert \psi_{\sigma\mu}^{(\bf k)}\rangle
\langle \psi_{\sigma\mu}^{(\bf k)}\vert,$ (2)

where $\int_{\rm B}$ means the integration over the first Brillouin zone of which volume is $V_{\rm B}$, and $\sum^{\rm occ}$ means the summation over occupied states. The charge density $n_{\sigma}(\bf r)$ with the spin index $\sigma$ is found as
$\displaystyle n_{\sigma}(\bf r)$ $\textstyle =$ $\displaystyle \langle {\bf r} \vert {\hat n}_{\sigma} \vert {\bf r} \rangle,$  
  $\textstyle =$ $\displaystyle \frac{1}{V_{\rm B}}
\int_{{\rm B}}dk^{3}
\sum_{\mu}^{\rm occ}
...a\mu}^{(\bf k)}\rangle
\langle \psi_{\sigma\mu}^{(\bf k)}\vert {\bf r} \rangle,$  
  $\textstyle =$ $\displaystyle \frac{1}{V_{\rm B}}
\int_{{\rm B}}dk^{3}
\sum_{{\rm n}}^{N}
...}({\bf r}-{\bf\tau}_{j})
\phi_{i\alpha}({\bf r}-{\bf\tau}_{i}-{\bf R}_{\rm n}),$  
  $\textstyle =$ $\displaystyle \sum_{{\rm n}}^{N}
\sum_{i\alpha, j\beta}
\rho_{\sigma, i\alpha j...
...a}({\bf r}-{\bf\tau}_{i})
\phi_{j\beta} ({\bf r}-{\bf\tau}_{j}-{\bf R}_{\rm n})$ (3)

with a density matrix defined by
    $\displaystyle \rho_{\sigma, i\alpha j\beta}^{({\bf R}_{\rm n})}
...\cdot {\bf k}}
c^{(\bf k)*}_{\sigma\mu,i\alpha}
c^{(\bf k)}_{\sigma\mu,j\beta}.$ (4)

Then, Mulliken populations $M_{\sigma,i\alpha}$ are given by
    $\displaystyle M_{\sigma,i\alpha}
\sum_{{\rm n}}^{N}
\rho_{\sigma, i\alpha j\beta}^{({\bf R}_{\rm n})}
S_{i\alpha,j\beta}^{({\bf R}_{\rm n})},$ (5)

where $S_{i\alpha,j\beta}^{({\bf R}_{\rm n})}$ is an overlap integral. Since the Mulliken population can be obtained by integrating Eq. (3) over real space, and by decomposing it into each contribution specified with $\sigma$ and $i\alpha$, it can be confirmed that the sum of $M_{\sigma,i\alpha}$ gives the number of electron $N_{\rm ele}$ per unit cell as follows:
    $\displaystyle N_{\rm ele} = \sum_{\sigma} \sum_{i\alpha} M_{\sigma,i\alpha}.$ (6)

Voronoi charge

The fuzzy cell method decomposes real space into smeared Voronoi cells, called the fuzzy cell [2]. The fuzzy cell at the site $i$ is determined by a weighting function $w_{i}({\bf r})$:
    $\displaystyle w_{i}({\bf r}) = \frac{p_{i}({\bf r})}{\sum_j p_{j}({\bf r})}$ (7)

with $p_i$ defined by
    $\displaystyle p_{i}({\bf r}) = \prod_{j\ne i} s_{k}(\mu_{ij}),$ (8)
    $\displaystyle \mu_{ij} = \frac{r_i - r_j}{\tau_{ij}},$ (9)
    $\displaystyle r_i = \vert {\bf r}-\tau_i\vert, \quad
r_j = \vert {\bf r}-\tau_j\vert, \quad
\tau_{ij} = \vert \tau_i-\tau_j\vert,$ (10)
    $\displaystyle s_{k}(x) = \frac{1}{2}
\right\},$ (11)
    $\displaystyle f_{k}(x) = f_{0}(f_{k-1}(x)),\quad f_{0} = \frac{3}{2}x - \frac{1}{2}x^3,$ (12)

where $k=3$ is chosen in OpenMX. As $k$ increases the fuzzy cells defined by $w$ approach to Voronoi cells (Wigner-Seitz cells). From the definition Eq. (7) it is clear that
    $\displaystyle \sum_{i} w_{i}({\bf r}) = 1.$ (13)

Thus, the integration of the charge density Eq. (3) over real space can be decomposed by employing the weighting functions as follows:
$\displaystyle \int dr^3 n_{\sigma}(\bf r)$ $\textstyle =$ $\displaystyle \int dr^3 [\sum_{i} w_{i}({\bf r})] n_{\sigma}(\bf r),$  
  $\textstyle =$ $\displaystyle \sum_{i}\int dr^3 w_{i}({\bf r}) n_{\sigma}(\bf r).$ (14)

Thus, the Voronoi charge $N_{\sigma,i}$ at the site $i$ can be defined by
    $\displaystyle N_{\sigma,i} = \int dr^3 w_{i}({\bf r}) n_{\sigma}(\bf r).$ (15)

From Eq. (14), it is confirmed that
    $\displaystyle N_{\rm ele} = \sum_{\sigma} \sum_{i} N_{\sigma,i}.$ (16)

Electro-static potential fitting (ESP) charge

Let us consider to express the Hartree potential in a system by the sum of Coulomb potentials with an effective point charge $Q_i$ located on each atomic site $\tau_i$ as follows:
    $\displaystyle V^{\rm (ESP)}({\bf r}) =
\sum_{i=1}^{N_{\rm atom}}\frac{Q_i}{\vert {\bf r} - \tau_i \vert},$ (17)

where $N_{\rm atom}$ is the number of atoms in the system. The $Q_i$ can be found by a least square fitting with a constraint $\sum_{i}Q_i=Q_{\rm tot}$ [3,4,5], where $Q_{\rm tot}$ is the total charge in the system. The Lagrange multiplier method casts this to a minimization problem of the following function $F$:
    $\displaystyle F = \sum_{p=1}^{N_{p}}\left(V^{\rm (DFT)}({\bf r}_p)
-V^{\rm (ESP)}({\bf r}_p)\right)^2
- \lambda\left(Q_{\rm tot} - \sum_{i}Q_i\right),$ (18)

where $V_p^{\rm (DFT)}$ and $V_p^{\rm (ESP)}$ are the Hartree potential calculated by the DFT calculation and Eq. (17), respectively, $\{{\bf r}_p\}$ is a set of sampling points, and $N_{p}$ is the number of the sampling points. The sampling points are given by the grids in the real space between two shells of the first and second scale factors times van der Waals radii [6]. The conditions $\frac{\partial F}{\partial Q_{i}}=0$ and $\frac{\partial F}{\partial \lambda}=0$ lead to
    $\displaystyle \left(
a_{11} & a_{12} & \cdots & a_{1N_{\rm...
b_{N_{\rm atom}}\\
Q_{\rm tot}
\end{array}\right)$ (19)

    $\displaystyle a_{ij} = \sum_{p=1}^{N_p}
\frac{1}{\vert {\bf r}_p - \tau_i \vert \vert {\bf r}_p - \tau_j \vert }$ (20)

    $\displaystyle b_{i} = \sum_{p=1}^{N_p}
\frac{V^{\rm (DFT)}({\bf r}_p)}{\vert {\bf r}_p - \tau_i \vert}.$ (21)

By solving the linear equation Eq. (19), we can find the electro-static potential fitting (ESP) charges. It is noted that the ESP charge is an effective charge on each atom including the contribution of the core charge compared to the Mulliken and Voronoi charges.


R. S. Mulliken, J. Chem. Phys. 23, 1833 (1955).

A. D. Becke and R. M. Dickson, J. Chem. Phys. 89, 2993 (1988).

U. C. Singh and P. A. Kollman, J. Comp. Chem. 5, 129(1984).

L. E. Chirlian and M. M. Francl, J. Com. Chem. 8, 894(1987).

B. H. Besler, K. M. Merz Jr. and P. A. Kollman, J. Comp. Chem. 11, 431 (1990).