<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN" "http://www.w3.org/TR/html4/loose.dtd">

<html><head>

    <meta http-equiv="Content-Type" content="text/html; charset=utf-8">

    <meta name="generator" content="Osso Notes">

    <title></title></head>

<body>

<p>

<br>ciao stefano,

<br>

<br>indeed, Oliviero is retesting in these days everything -

<br>when we checked it with Ismaila (fig 5 of the 2008 prb)

<br>it seemed correct, but we need to recheck units, and the

role of the quadrupole.

<br>

Soon we'll have an update.

<br>

<br>Nicola

<br>

<br>----- Original message -----

<br>> dear Nicola, dear all,

<br>>     admittedly I know nothing about your 2008 PRB with Ismaila and there 

> may be some problem of units or other glitch in th eimplementation ..

<br>>     but on the simple tests that are contained in the cluster_example 

<br>> directory (H20, NH4+ and atomic N)   the new MP correction is clearly 

<br>> inferior to the old one in the only case where it is not negligible, the 

<br>> charged NH4+ system.

<br>> 

<br>>   stefano

<br>> 

<br>> comparison follows between total energies of a few molecule in simple 

<br>> cubic boxes with 16,20,24 a.u. sides

<br>> The reference results are those that can be obtained   by running   

<br>> run_example in cluster_example and correspond to Martyna Tuckeman 

<br>> correction.

<br>> 

<br>> NEW MP correction

<br>> results-mp_new/h2o.out-16:!       Total+Makov-Payne energy   =         

<br>> -44.00269835 Ry

<br>> results-mp_new/h2o.out-20:!       Total+Makov-Payne energy   =         

<br>> -44.00238516 Ry

<br>> results-mp_new/h2o.out-24:!       Total+Makov-Payne energy   =         

<br>> -44.00256923 Ry

<br>> results-mp_new/nh4+.out-16:!       Total+Makov-Payne energy   =         

<br>> -32.45635049 Ry

<br>> results-mp_new/nh4+.out-20:!       Total+Makov-Payne energy   =         

<br>> -32.45333443 Ry

<br>> results-mp_new/nh4+.out-24:!       Total+Makov-Payne energy   =         

<br>> -32.45292630 Ry

<br>> results-mp_new/n.out-16:!       Total+Makov-Payne energy   =       

<br>> -27.82758316 Ry results-mp_new/n.out-20:!       Total+Makov-Payne energy   =

<br>>         -27.82739183 Ry results-mp_new/n.out-24:!       Total+Makov-Payne

<br>> energy   =         -27.82770125 Ry

<br>> 

<br>> OLD MP correction

<br>> results-mp_old/h2o.out-16:!       Total+Makov-Payne energy   =         

<br>> -44.00242684 Ry

<br>> results-mp_old/h2o.out-20:!       Total+Makov-Payne energy   =         

<br>> -44.00224753 Ry

<br>> results-mp_old/h2o.out-24:!       Total+Makov-Payne energy   =         

<br>> -44.00249007 Ry

<br>> results-mp_old/nh4+.out-16:!       Total+Makov-Payne energy   =         

<br>> -32.45145525 Ry

<br>> results-mp_old/nh4+.out-20:!       Total+Makov-Payne energy   =         

<br>> -32.45083466 Ry

<br>> results-mp_old/nh4+.out-24:!       Total+Makov-Payne energy   =         

<br>> -32.45148127 Ry

<br>> results-mp_old/n.out-16:!       Total+Makov-Payne energy   =       

<br>> -27.82758316 Ry results-mp_old/n.out-20:!       Total+Makov-Payne energy   =

<br>>         -27.82739183 Ry results-mp_old/n.out-24:!       Total+Makov-Payne

<br>> energy   =         -27.82770125 Ry

<br>> 

<br>> Martyna Tuckerman approach

<br>> reference/h2o.out-16:!       total energy                           =         -44.00239103 Ry

<br>> reference/h2o.out-20:!       total energy                           =         -44.00224062 Ry

<br>> reference/h2o.out-24:!       total energy                           =         -44.00248682 Ry

<br>> reference/nh4+.out-16:!       total energy                           =         -32.45134332

<br>> Ry reference/nh4+.out-20:!       total energy                           =       

<br>> -32.45079847 Ry reference/nh4+.out-24:!       total energy                           =   

<br>>     -32.45146630 Ry reference/n.out-16:!       total energy                           =     

<br>>   -27.82757465 Ry reference/n.out-20:!       total energy                           =       

<br>> -27.82739133 Ry reference/n.out-24:!       total energy                           =       

<br>> -27.82770111 Ry

<br>> 

<br>> 

<br>> 

<br>> 

<br>> Nicola Marzari wrote:

<br>> > Ciao Paolo,

<br>> > 

<br>> > 

<br>> > I spoke today with Oliviero Andreussi (that wrote the email

> > below) and he was also going to do some tests to try out.

<br>> > 

<br>> > Fig. 5 in the 2008 PRB Dabo at al was showing the discrepancy between 

> > the original M-P formula and our corrected one.

<br>> > 

<br>> > 

<br>> >             nic

<br>> > 

<br>> > 

<br>> > 

<br>> > On 12/3/10 7:04 PM, Paolo Giannozzi wrote:

<br>> > 

<br>> > > It would be useful if anybody with previous Makov-Payne

<br>> > > experience could make some tests

<br>> > > 

<br>> > > P.

<br>> > > 

<br>> > > 

<br>> > > > Dear all,

<br>> > > > 

<br>> > > > I checked the equations in the subroutine PW/makov-payne.f90 and it

<br>> > > > seems to me that the equation implemented is a modified version of

<br>> > > > the original Makov-Payne correction of 1995. The modification is

<br>> > > > due to the fact that in the original derivation the second

<br>> > > > correction was computed assuming that pointlike charge that

<br>> > > > compensate the molecular density was placed such as to give a zero

<br>> > > > net dipole moment, while the implemented version is more general

<br>> > > > and thus include a dipolar term. Moreover the correction has the

> > > > correct sign (the opposite of the one in the original equation).

<br>> > > > This to say that the equation does not contain the factor one half

<br>> > > > suggested in the paper by Ismaila and Nicola.

<br>> > > > I would suggest the following changes to some of the lines of the

<br>> > > > subroutine:

<br>> > > > 

<br>> > > > from

<br>> > > > 

<br>> > > > 

<br>> > > > !

<br>> > > > ! ... Makov-Payne correction, PRB 51, 4014 (1995)

<br>> > > > ! ... Note that Eq. 15 has the wrong sign for the quadrupole term

<br>> > > > !

<br>> > > > corr1 = - madelung(ibrav) / alat * qq**2

<br>> > > > !

<br>> > > > aa = quadrupole

<br>> > > > bb = dipole(1)**2 + dipole(2)**2 + dipole(3)**2

<br>> > > > !

<br>> > > > corr2 = ( 4.D0 / 3.D0 * pi )*( qq*aa - bb ) / alat**3

<br>> > > > 

<br>> > > > 

<br>> > > > to (I also made explicit the factor e2/2.0 to convert from Hartree

<br>> > > > to Rydberg, that was implicitly in the original implementation, so

<br>> > > > e2 should also be included in the declaration of the constants)

<br>> > > > 

<br>> > > > !

<br>> > > > ! ... Makov-Payne correction, PRB 51, 4014 (1995)

<br>> > > > ! ... Note that Eq. 15 has the wrong sign for the quadrupole term

<br>> > > > ! ... The prefactor of the quadrupole term is half the one in Eq.

<br>> > > > 15 ! ... according to the results derived by Dabo et al., PRB 77,

<br>> > > > 115139 (2008)

<br>> > > > !

<br>> > > > corr1 = - madelung(ibrav) / alat * qq**2 / 2.0D0 * e2

<br>> > > > !

<br>> > > > aa = quadrupole

<br>> > > > bb = dipole(1)**2 + dipole(2)**2 + dipole(3)**2

<br>> > > > !

<br>> > > > corr2 = ( 1.D0 / 3.D0 * pi )*( qq*aa - bb ) / alat**3 * e2

<br>> > > > 

<br>> > > > Similarly, for the routine CPV/makov_payne.f90

<br>> > > > 

<br>> > > > !

<br>> > > > ! ... Makov-Payne correction, PRB 51, 43014 (1995)

<br>> > > > ! ... Note that Eq. 15 has the wrong sign for the quadrupole term

<br>> > > > !

<br>> > > > ! 1 / 2 Ry ->   a.u.

<br>> > > > corr1 = - 2.8373D0 / alat * charge**2 / 2.0D0

<br>> > > > !

<br>> > > > aa = quadrupole

<br>> > > > bb = dipole(1)**2 + dipole(2)**2 + dipole(3)**2

<br>> > > > !

<br>> > > > ! 4 / 3 ->   2 / 3 Ry ->   a.u.

<br>> > > > corr2 = ( 2.D0 / 3.D0 * pi )*( charge*aa - bb ) / alat**3

<br>> > > > 

<br>> > > > to (in this case I assumed that e2 is defined in CP and with a

<br>> > > > value equal to 1.0,

<br>> > > > otherwise, probably is just better to forget about it)

<br>> > > > 

<br>> > > > !

<br>> > > > ! ... Makov-Payne correction, PRB 51, 43014 (1995)

<br>> > > > ! ... Note that Eq. 15 has the wrong sign for the quadrupole term

<br>> > > > ! ... The prefactor of the quadrupole term is half the one in Eq.

<br>> > > > 15 ! ... according to the results derived by Dabo et al., PRB 77,

<br>> > > > 115139 (2008)

<br>> > > > !

<br>> > > > ! 1 / 2 Ry ->   a.u.

<br>> > > > corr1 = - 2.8373D0 / alat * charge**2 / 2.0D0 * e2

<br>> > > > !

<br>> > > > aa = quadrupole

<br>> > > > bb = dipole(1)**2 + dipole(2)**2 + dipole(3)**2

<br>> > > > !

<br>> > > > ! 2/ 3 ->   1 / 3 Ry ->   a.u.

<br>> > > > corr2 = ( 1.D0 / 3.D0 * pi )*( charge*aa - bb ) / alat**3 * e2

<br>> > > > 

<br>> > > > Please, if someone else can also check these corrections before

<br>> > > > implementing them it would be better.

<br>> > > > 

<br>> > > ---

<br>> > > Paolo Giannozzi, Democritos and University of Udine, Italy

<br>> > > 

<br>> > > 

<br>> > > _______________________________________________

<br>> > > Q-e-developers mailing list

<br>> > > <a href="mailto:Q-e-developers@qe-forge.org">Q-e-developers@qe-forge.org</a>

<br>> > > <a href="http://qe-forge.org/mailman/listinfo/q-e-developers">http://qe-forge.org/mailman/listinfo/q-e-developers</a>

<br>> > > 

<br>> > 

<br>> > 

<br>> > 

<br>> 

<br><br></p>

</body>

</html>