[Q-e-developers] Fwd: Makov-Payne correction

[Paolo Giannozzi]

It would be useful if anybody with previous Makov-Payne
experience could make some tests

P.

> Dear all,
>
> I checked the equations in the subroutine PW/makov-payne.f90 and it  
> seems to me that the equation implemented is a modified version of  
> the original Makov-Payne correction of 1995. The modification is  
> due to the fact that in the original derivation the second  
> correction was computed assuming that pointlike charge that  
> compensate the molecular density was placed such as to give a zero  
> net dipole moment, while the implemented version is more general  
> and thus include a dipolar term. Moreover the correction has the  
> correct sign (the opposite of the one in the original equation).
> This to say that the equation does not contain the factor one half  
> suggested in the paper by Ismaila and Nicola.
> I would suggest the following changes to some of the lines of the  
> subroutine:
>
> from
>
>
>   !
>   ! ... Makov-Payne correction, PRB 51, 4014 (1995)
>   ! ... Note that Eq. 15 has the wrong sign for the quadrupole term
>   !
>   corr1 = - madelung(ibrav) / alat * qq**2
>   !
>   aa = quadrupole
>   bb = dipole(1)**2 + dipole(2)**2 + dipole(3)**2
>   !
>   corr2 = ( 4.D0 / 3.D0 * pi )*( qq*aa - bb ) / alat**3
>
>
> to (I also made explicit the factor e2/2.0 to convert from Hartree  
> to Rydberg, that was implicitly in the original implementation, so  
> e2 should also be included in the declaration of the constants)
>
>   !
>   ! ... Makov-Payne correction, PRB 51, 4014 (1995)
>   ! ... Note that Eq. 15 has the wrong sign for the quadrupole term
>   ! ... The prefactor of the quadrupole term is half the one in Eq. 15
>   ! ... according to the results derived by Dabo et al., PRB 77,  
> 115139 (2008)
>   !
>   corr1 = - madelung(ibrav) / alat * qq**2 / 2.0D0 * e2
>   !
>   aa = quadrupole
>   bb = dipole(1)**2 + dipole(2)**2 + dipole(3)**2
>   !
>   corr2 = ( 1.D0 / 3.D0 * pi )*( qq*aa - bb ) / alat**3 * e2
>
> Similarly, for the routine CPV/makov_payne.f90
>
>    !
>    ! ... Makov-Payne correction, PRB 51, 43014 (1995)
>    ! ... Note that Eq. 15 has the wrong sign for the quadrupole term
>    !
>    ! 1 / 2 Ry -> a.u.
>    corr1 = - 2.8373D0 / alat * charge**2 / 2.0D0
>    !
>    aa = quadrupole
>    bb = dipole(1)**2 + dipole(2)**2 + dipole(3)**2
>    !
>    ! 4 / 3 -> 2 / 3 Ry -> a.u.
>    corr2 = ( 2.D0 / 3.D0 * pi )*( charge*aa - bb ) / alat**3
>
> to (in this case I assumed that e2 is defined in CP and with a  
> value equal to 1.0,
> otherwise, probably is just better to forget about it)
>
>    !
>    ! ... Makov-Payne correction, PRB 51, 43014 (1995)
>    ! ... Note that Eq. 15 has the wrong sign for the quadrupole term
>   ! ... The prefactor of the quadrupole term is half the one in Eq. 15
>   ! ... according to the results derived by Dabo et al., PRB 77,  
> 115139 (2008)
>    !
>    ! 1 / 2 Ry -> a.u.
>    corr1 = - 2.8373D0 / alat * charge**2 / 2.0D0 * e2
>    !
>    aa = quadrupole
>    bb = dipole(1)**2 + dipole(2)**2 + dipole(3)**2
>    !
>    ! 2/ 3 -> 1 / 3 Ry -> a.u.
>    corr2 = ( 1.D0 / 3.D0 * pi )*( charge*aa - bb ) / alat**3 * e2
>
> Please, if someone else can also check these corrections before  
> implementing them it would be better.

---
Paolo Giannozzi, Democritos and University of Udine, Italy