Wrong normalization of projector for f-channel of GTH pseudopotentials
It seems that GTH f-projector is not normalized correctly, this will result in large Etotal for systems with f-projector.
Example input file:
&CONTROL
calculation = 'scf'
restart_mode = 'from_scratch'
pseudo_dir = './'
outdir = './tmp'
verbosity = 'high'
disk_io = 'none'
/
&SYSTEM
ibrav = 0
nat = 2
ntyp = 1
ecutwfc = 60.0000000000
nbnd = 14
occupations = 'smearing'
smearing = 'fermi-dirac'
degauss = 0.0020000000
/
&ELECTRONS
electron_maxstep = 150
mixing_beta = 0.1
/
ATOMIC_SPECIES
Ba 1.0 Ba.gth
ATOMIC_POSITIONS bohr
Ba 0.0000000000 0.0000000000 0.0000000000
Ba 4.7506958904 4.7506958904 4.7506958904
K_POINTS automatic
6 6 6 0 0 0
CELL_PARAMETERS bohr
9.5013917807 0.0000000000 0.0000000000
0.0000000000 9.5013917807 0.0000000000
0.0000000000 0.0000000000 9.5013917807
Pseudopotential Ba.gth (parameters taken from CP2K, only the h's, header is taken from ABINIT and modified slightly):
Hartwigsen-Goedecker-Hutter psp for Ba, from PRB58, 3641 (1998)
56 10 010605
10 1 3 0 2001 0
0.5400000000 2 24.4786532500 -2.5008499400
4
0.5147761300 2 1.0467290300 0.3784745400
-0.9772170700
0.3751901100 2 -0.2024391200 -0.2874921700
0.6803306500
0.0000000000 0.0000000000
0.0000000000
0.6654031900 1 0.3784185000
0.0000000000
0.3049204300 1 -18.7952077700
0.0000000000
The results I obtained from qe-6.3:
the Fermi energy is ********** ev
! total energy = -57343.16939753 Ry
Harris-Foulkes estimate = -57343.16939753 Ry
estimated scf accuracy < 0.00000012 Ry
The total energy is the sum of the following terms:
one-electron contribution = -57370.28077685 Ry
hartree contribution = 142.50776655 Ry
xc contribution = -38.79187492 Ry
ewald contribution = -76.60421830 Ry
smearing contrib. (-TS) = -0.00029403 Ry
convergence has been achieved in 3 iterations
The problem seems to be solved when I added factor of 1/sqrt(105.0) for f-projector.
File Modules/gth.f90
!....
else if (ll==3) then lif ! f channel [ ONLY 1 PROJECTOR!! ]
!
do ii=1,nq
qr2=(qg(ii)*rrl)**2
vq(ii)=qg(ii)**3 * exp(-0.5_dp*qr2) / sqrt(105.0_dp)
end do
!
end if lif
!....
With the above modification, this is the result:
the Fermi energy is 5.5666 ev
! total energy = -102.36711325 Ry
Harris-Foulkes estimate = -102.36711341 Ry
estimated scf accuracy < 0.00000032 Ry
The total energy is the sum of the following terms:
one-electron contribution = -28.47048475 Ry
hartree contribution = 17.19287411 Ry
xc contribution = -14.48477506 Ry
ewald contribution = -76.60421830 Ry
smearing contrib. (-TS) = -0.00050924 Ry
convergence has been achieved in 16 iterations
Etotal = -51.183556625 Ha, which is quite the same with ABINIT's result.
Result from ABINIT (with the same cutoff = 30 Ha)
Components of total free energy (in Hartree) :
Kinetic energy = 1.32533731925961E+01
Hartree energy = 8.60112605940208E+00
XC energy = -7.24291638656629E+00
Ewald energy = -3.83021091383786E+01
PspCore energy = 2.81771207578505E+00
Loc. psp. energy= -3.01293429865171E+01
NL psp energy= -1.80570957421657E-01
>>>>> Internal E= -5.11827281411004E+01
-kT*entropy = -2.54776730344274E-04
>>>>>>>>> Etotal= -5.11829829178307E+01
Other information on the energy :
Total energy(eV)= -1.39275979527691E+03 ; Band energy (Ha)= -9.3210373654E+00
I have not tested the force yet, maybe someone with more experiences can confirm.