Bugs in ELF.F

[chengcheng_xiao1]

I think I might have found two bugs in the `ELF.F` file.

1. There's a typo in the docstring saying the kinetic energy density is twice the value as it should be. (not very sure about this, might want to double check)

2. The coefficient for DH term seem to be wrong, DH should be

Screen Shot 2022-01-20 at 3.31.00 PM.png

To fix the second one, the following code in `ELF.F`:

Code: Select all

DH = 0.2_q*HSQDTM/PISQ* (3*PISQ* REAL( CHDEN(N) ,KIND=q) )**FIVTHI

needs to be changed to:

Code: Select all

DH = 0.2_q*HSQDTM/PISQ* (6*PISQ* REAL( CHDEN(N) ,KIND=q) )**FIVTHI

I've summarized the problem in the attached pdf file, hope it helps!

[chengcheng_xiao1]

Okay I just realized that I wasn't very clear about the origin of the second bug I found. Let me elaborate:

For spin polarized charge density (J. Chem. Phys. 92, 5397 (1990)):
D^{\sigma}=\sum_i | \nabla \psi_i | ^2 - \frac{1}{4} \frac{|\nabla \rho|^2}{\rho}
D^{\sigma}_0 = \frac{3}{5} (6 \pi^2)^{2/3} \rho^{5/3}

For spin degenerate charge density (Nature 371, 683-686 (1994)):
D=\frac{1}{2} [ \sum_i | \nabla \psi_i | ^2 - \frac{1}{4} \frac{|\nabla \rho|^2}{\rho} ]
D_0 = \frac{3}{10} (3 \pi^2)^{2/3} \rho^{5/3}

In ELF.F, a big loop over spin component is done and the charge density of each spin component is used to calculate the ELF. So, by definition, one should use the spin polarized version of D^{sigma}_0 but the the spin degenerate version D_0 is used.

[fabien_tran1]

Hi,

Yes there is indeed a bug in elf.F. The ELF function in spin-polarized calculations is wrong. However, your fix was not complete. The fix should be

Code: Select all

diff elf.F elf.F_with_fix
48a49,50
> 
>       NSPIN=MIN(WDES%NCDIJ,2)
215c217
<       CALL ELFCAL(CHDEN(1,ISP),CR,CDWORK,CNEW,GRID_SOFT)
---
>       CALL ELFCAL(CHDEN(1,ISP),CR,CDWORK,CNEW,GRID_SOFT,NSPIN)
243c245
<       SUBROUTINE ELFCAL(CHDEN,LAPLAC,CKINE,CGRDSQ,GRID)
---
>       SUBROUTINE ELFCAL(CHDEN,LAPLAC,CKINE,CGRDSQ,GRID,NSPIN)
259c261
< ! T    = - 2 --- Psi grad Psi   T+TCORR...pos.definite kinetic energy
---
> ! T    =   - --- Psi grad Psi   T+TCORR...pos.definite kinetic energy
285c287
<           DH = 0.2_q*HSQDTM/PISQ* (3*PISQ* REAL( CHDEN(N) ,KIND=q) )**FIVTHI
---
>           DH = 0.2_q/NSPIN*HSQDTM/PISQ* (3*PISQ* REAL( NSPIN*CHDEN(N) ,KIND=q) )**FIVTHI

With this fix, for non-magnetic systems the ELF functions obtained with ISPIN=1 and ISPIN=2 are the same, as expected. Thank you for having reported the bug.

[chengcheng_xiao1]

Hi Fabian,

Yes you are right! One needs to treat spin degenerate and polarized calculations differently.

Also, for the uninitiated, to patch the elf.F file, simply put the code provided by Fabian in a file named `elf.patch` under `src` folder and then issue the following command:

patch elf.F elf.patch

Hope this helps! 🥳

[foyevtsova]

Hi,

I have a question about how to interpret results of ELF calculations for a ISPIN=2 case with an antiferromagnetic configuration of localized magnetic moment.

My intuitive expectation would be that the spin-up and spin-down sets of data in the ELFCAR file would be similar but shifted in real space in accordance with the ordering wave-vector like, for instance, Q=(pi,pi) for a checkerboard type of magnetic ordering.

What I see, however, is a different picture: The first set of data shows no difference between the spin-up and the spin-down sites with localized magnetic order. The second set of data does show the expected difference but the EFL values are extremely spiky, reaching almost 1 within strongly confined regions in space.

Thank you,
Kateryna Foyevtsova