surfstat_jdos compatibility with Quantum Espresso

[dengtq]

Dear Dr Wu

I encountered strange results when computing the spin component of surface spectrum. I found the issue may originate from the surfstat_jdos subroutine where the spin projectors were only defined in VASP fashion between line 705-717 in surfstat.f90:

!> spin operator matrix
!> Note, the basis here should be |↑↑↓↓>
nw_half = Num_wann/2
DO i=1, Np
    DO j=1, nw_half
        sigma_x( Num_wann*(i-1)+j        , Num_wann*(i-1)+j+nw_half ) =  1.0d0
        sigma_x( Num_wann*(i-1)+j+nw_half, Num_wann*(i-1)+j         ) =  1.0d0
        sigma_y( Num_wann*(i-1)+j        , Num_wann*(i-1)+j+nw_half ) = -zi
        sigma_y( Num_wann*(i-1)+j+nw_half, Num_wann*(i-1)+j         ) =  zi
        sigma_z( Num_wann*(i-1)+j        , Num_wann*(i-1)+j         ) =  1.0d0
        sigma_z( Num_wann*(i-1)+j+nw_half, Num_wann*(i-1)+j+nw_half ) = -1.0d0
    ENDDO
ENDDO

After replacing these lines (following fermiarc.f90) with the following, the result then looks reasonable:

nw_half= Num_wann/2
if (index( Package, 'VASP')/=0.or. index( Package, 'Wien2k')/=0 &
   .or. index( Package, 'Abinit')/=0.or. index( Package, 'openmx')/=0) then
   do i=1, Np
      do j=1, nw_half
         sigma_x(Num_wann*(i-1)+j,         Num_wann*(i-1)+nw_half+j) = 1.0d0
         sigma_x(Num_wann*(i-1)+j+nw_half, Num_wann*(i-1)+j        ) = 1.0d0
         sigma_y(Num_wann*(i-1)+j,         Num_wann*(i-1)+nw_half+j) =-zi
         sigma_y(Num_wann*(i-1)+j+nw_half, Num_wann*(i-1)+j        ) = zi
         sigma_z(Num_wann*(i-1)+j,         Num_wann*(i-1)+j        ) = 1.d0
         sigma_z(Num_wann*(i-1)+j+nw_half, Num_wann*(i-1)+j+nw_half) =-1.d0
      enddo
   enddo
elseif (index( Package, 'QE')/=0.or.index( Package, 'quantumespresso')/=0 &
     .or.index( Package, 'quantum-espresso')/=0.or.index( Package, 'pwscf')/=0) then
   do i=1, Np
      do j=1, nw_half
         sigma_x(Num_wann*(i-1)+(2*j-1), Num_wann*(i-1)+2*j    ) = 1.0d0
         sigma_x(Num_wann*(i-1)+2*j,     Num_wann*(i-1)+(2*j-1)) = 1.0d0
         sigma_y(Num_wann*(i-1)+(2*j-1), Num_wann*(i-1)+2*j    ) =-zi
         sigma_y(Num_wann*(i-1)+2*j,     Num_wann*(i-1)+(2*j-1)) = zi
         sigma_z(Num_wann*(i-1)+(2*j-1), Num_wann*(i-1)+(2*j-1)) = 1.0d0
         sigma_z(Num_wann*(i-1)+2*j,     Num_wann*(i-1)+2*j    ) =-1.0d0
      enddo
   enddo
else
   if (cpuid.eq.0) write(stdout, *)'Error: please report your software generating tight binding and wannier90.wout to me'
   if (cpuid.eq.0) write(stdout, *)'wuquansheng@gmail.com'
   stop 'Error: please report your software and wannier90.wout to wuquansheng@gmail.com'
endif

Is this a correct fix for Quantum Espresso?

Best Regards,

Deng Tianqi, Institute of High Performance Computing, A*STAR, Singapore