Developers Manual:Separation of the pseudopotential

From OctopusWiki
Jump to: navigation, search

We will separate the local part of the pseudpotential into a long-range and a short-range part. The long-range part will be the Coulomb potential associated with a radial Gaussian charge distribution:

\rho(r)=\frac{Z_{val}}{\left(\sqrt{2\pi}\sigma\right)^3}e^{-\frac{r^2}{\sqrt(2\pi)\sigma}}\ .

The potential is

\Phi(r) = \frac{Z_{val}}{r}\mbox{erf}\left(\frac{r}{\sqrt{2}\sigma}\right)\ ,

\Phi(r\rightarrow0) \longrightarrow \frac{2Z_{val}}{\sqrt{2\pi}\sigma}\ ,

\Phi(r\rightarrow\infty) \longrightarrow -\frac{Z_{val}}{r}\ .

In Fourier space:

\rho(q)=\frac{Z_{val}}{2\pi\sigma}e^{-\frac12{\sigma^2q^2}}\ .

We want the potential generated by this charge to be smooth so we consider the associated potential (for simplicity we consider the 1D Poisson equation)

\Phi(q) = \frac1{q^2}\left(-4\pi\rho(q)\right) = -\frac{2Z_{val}}{q^2\sigma^2}e^{-\frac12\sigma^2q^2}

What we want is that \left|\Phi(q_{max})\right| must be smaller than a certain threshold:

\frac{\Phi(q_{max})}{Z_{val}} = \frac{e^{-x}}x

with x = \frac12\sigma^2q_{max}^2.

Once we get x, we can get the optimum value of \sigma:

\sigma = \frac{\sqrt{2x}}{q_{max}} = \frac{\sqrt{2x}}{\pi}h\ .

For  \frac{e^{-x}}x = 0.001 we get x=5.2496, which implies

\sigma = 1.03 h\ .

This value will be used in the code. In principle there is a dependence on Z_{\rm val}, but solving a transcendental equation in the code is not very reliable.

Note: this has been taken from the code; in case of discrepancies the code is more likely to be correct.