Introduction
The reason for applying pseudopotential(PP) method arose from the complexity of “wiggles” in core wavefunction. This can be interpreted qualitatively by a simple 1d infinity square well model. Core electron will behave more like in a infinity square well model with shorter distance between electron and core. Thus rapid oscillation occurs in the core wavefunction region. The valence wavefunction must be orthogonal to core wavefunction, thus the dimensionality of the Hamiltonian matrix will be tremendous, resulting much unnecessary computing time.
The PP method is created to deal with such problem. The word ”pseudo” means “fake” in Greek. We use a “fake” potential to approximate the “real” potential without the ”wiggles” in the core wavefunction region. Since it’s a approximation to the real wavefunction, we need to have the same results as the real wavefunction. This is rule (1): real and pesudo eigenvalues must be the same. Also, since we are trying to cancel the “wiggles” in the core region, beyond this region, the pesudo wavefunction should be the same as real wavefunction, which leads to rule (2): beyond cutoff radius $r_c$ , real and pesudo wavefunction should be the same.
Another import issue is the transferability of the PP, i.e. the ability for this PP to be applied in other chemical environment. For PP, it may be pure empirical so that one PP can only be used in some limited cases, which is not desirable. In order to achieve optimum transferability, we should make the scattering properties of the pseudo wavefunction agree with the real one. So we have rule (3): in core region($r < r_c$ ), the pesudo probability of the wavefunction should agree with the real one. Beyond the core region($r > r_c$), we have rule (4): the logarithmic derivatives of the real and pesudo wave function and their first energy derivatives agree.(See Derivation section) Thus we can summarize the four criteria:
- $\varepsilon_p s = \varepsilon e$
- $\psi_p$ is nodeless
- $\psi_p s = \psi_r e$ for $r > r_c $
- $\int_{r < r_c} \lvert\psi_ps\rvert^2r^2dr = \int_{r < r_c} \lvert\psi_re\rvert^2r^2dr$
Deriviation
Relation Between Logarithmic Derivative and Scattering Properties
Recall a radial Schördinger equation:
If V is spherical, we can rewrite $\psi $ in spherical coordinate
where $\mathcal{Y}_{lm}(\theta,\phi)$ is the spherical harmonic.
Then we can further simplify the equation by putting $\textbf{p}$ parallel to $z$-axis. $\theta$ is the angle between $\textbf{p}$ and the axis, then $\phi$ becomes independent and can be reduced.
Inserting into Eq.$\ref{radialSch}:$
Here we make $\epsilon_p=p^2/2m$ so that $u_{l,p}(0)=0$ at the origin.
Then we consider $r\rightarrow\infty$ in Eq.$\ref{radialSch2}$. Assume the potential vanishes for $r\rightarrow\infty$, then Eq.$\ref{radialSch2}$ can be rewritten as:
Eq.$\ref{radialSchInfty}$ has the following solution:
where $j_l(pr)$ and $n_l(pr)$ are the spherical Bessel and Neumann functions of order $l$. For $pr\rightarrow0$, the Neumann function is eliminated, and for $pr\rightarrow\infty$, the Bessel function is eliminated so that the wavefunctions can be “regular”.
Thus the radial components has the following asymptotic behavior:
which is equivalent to
where $\delta_l(p)$ is the \textit{phase shift} due to the potential, describing a shift with respect to the free wave(no potential). This phase shift depend not only on potential but also the scattering energy $E_p$. So we could rewrite the Eq.$\ref{radialsolution}$ as:
Compare Eq.$\ref{solution1}$ with Eq.$\ref{solution2}$, we can get
As $A_{l,p}$ is a free cte from the normalization of the wavefunction, we can write the radial solution as:
From Eq.$\ref{shifteqn}$, we can get the phase shift $\delta_l$
Practically, the phase shift $\tan\delta_l$ is obtained by solving the radial Schödinger when the potential $V(r)$ vanished(just like Eq.$\ref{radialSchInfty}$). Then $\tan\delta_l$ can be obtained from Eq.$\ref{deltader}$ together with its derivative:
where $\gamma$ is the logarithmic derivative of the solution $u_{l,p}(r)/r$
This is the relation between the logarithmic derivative of the wavefunction and the scattering property. The logarithmic derivative must agree so that the transferability is guaranteed.
Identity Relating Rule(3) and (4)
The identity relating rule(3) and (4) is given as Eq.(1) in the NCPP paper.
No explanation is given in the paper, but I manage to find a detailed derivation in the appendix of W. C. Topp and J. J. Hopfield, Phys. Rev. B 7, 1295 (1974).(DOI:10.1103/PhysRevB.7.1295).
This equation is acquired by integrate the Schrödinger equation in the spherical coordinate. Use Green’s theorem, and put all energy term in the LHS, we can get
So that if we can find correct ground and excited state energy for a certain symmetry configuration, we can get the correct normalization of the wavefunction. And for $r$ beyond $r_c$, the wavefunction is determined by the Coulomb tail to the potential, thus we can have correct charge distribution in the tail region, which all together gives the generalized charge density in the unit cell of the solid.
