Electron-nuclear potential (-Z/r)

Previously, you solved the Kohn-Sham equations in the very simple potential coming from a simple harmonic oscillator. For most electronic structure calculations, however, the potential comes from the interaction between electrons and nuclei with the form

\begin{displaymath}
V(\vec r) = - \sum_I \frac{Z_I}{\vert\vec r-\vec r_I\vert},
\end{displaymath} (1)

where $Z_I$ and $\vec r_I$ are the charge and location of each nucleus.

Our approach to funding the potential Vdual from the nuclei for our periodic calculations is based on our standard approach to solving Poisson's equation,

phinuc=cI(Linv((-4*pi)*O(cJ(rhonuc))));
and our the definition of Vdual from Problem Set #8,
Vdual=cJdag(O(cJ(-phinuc)));
where the extra `-' in front of phinuc corresponds to the minus sign in Eq. (1). Combining the above two results gives
Vdual=cJdag(O(Linv((4*pi)*O(cJ(rhonuc)))));
where we have used the fact that ${\cal J}({\cal I}(\vec v))=\vec v$.

The only difficulty with this approach is that rhonuc represents a point charge, whose density is zero everywhere but at a single point where the density is infinite. We can, however, find the expansion coefficients $[$cJ(rhonuc) $]_{\vec G} \equiv \hat \rho_{\vec G}$ mathematically. To do this, we note that

\begin{displaymath}
\rho(\vec r) \equiv \sum_{\vec G} \hat \rho_{\vec G} e^{i\vec G \cdot
\vec r}
\end{displaymath}

so that

\begin{eqnarray*}
\int e^{-i\vec G \cdot \vec r} \rho(\vec r)\,dV & = & \int e^{...
...' \cdot
\vec r}\,dV\\
& = & \det({\bf R})\, \hat \rho_{\vec G},
\end{eqnarray*}


where $\det({\bf R})$ is the volume of the cell. Therefore,

\begin{displaymath}
\hat \rho_{\vec G} = \frac{1}{\det({\bf R})} \int e^{-i\vec G \cdot \vec r} \rho(\vec r)\,dV.
\end{displaymath}

In the case of a single nucleus, $\rho(\vec r)$ is a point charge centered at location $\vec r_I$ with total integral $Z_I$. Thus, we have

\begin{displaymath}
\hat \rho_{\vec G} = \frac{Z_I}{\det({\bf R})} e^{-i\vec G \cdot \vec r_I}.
\end{displaymath}

In octave notation, the above expression becomes
rhohat=ZI/det(R)*exp(-i*G*rI);
where ZI is the charge of the nucleus and rI is a column vector containing its location. Finally, for the case of multiple nuclei, one should sum a series of such contributions, one for each nucleus.



Subsections
Tomas Arias 2005-04-22