I use vasp5.2 to calculate bcc Fe, and I want to obtain optic information with spin-orbit couping. How I can output optic matrix element, i.e. < phi_i | nabla | phi_j >?
I have tried to add some code in "optics.F" aroud the lines 156 ~ 191,
Code: Select all
!-----------------------------------------------------------------------
! calculate < phi_i | nabla | phi_j >
!-----------------------------------------------------------------------
!ad---------------by fwx
      io_begin
      open(91,file="nabij.out",status="unknown")
        open(92,file="Px.dat",status="unknown")
        open(93,file="Py.dat",status="unknown")
        open(94,file="Pz.dat",status="unknown")
      write(92,*) KPOINTS%NKPTS, WDES%NB_TOT
      write(93,*) KPOINTS%NKPTS, WDES%NB_TOT
      write(94,*) KPOINTS%NKPTS, WDES%NB_TOT
      io_end
!ad------------------ by fwx
! loop over all (current) k-points, spin component and cartesian directions
      DO NK=1,KPOINTS%NKPTS
      DO ISP=1,INFO%ISPIN
      !  Do isp = 1,2
        DO IDIR=1,3
! get the matrix element of the nabla-operator in eV/Angstroem units
         CALL NABIJ_SOFT(NABIJ,IDIR,NK,ISP,W,WDES,GRID_SOFT,LATT_CUR)
         CALL NABIJ_AUG_ADD(NABIJ,IDIR,NK,ISP,W,WDES,P,T_INFO)
         NB_TOT=WDES%NB_TOT
         CALLMPI( M_sum_g(WDES%COMM,NABIJ(1,1),NB_TOT*NB_TOT))
!ad-----by fwx-------------
        io_begin
        if(idir==1)then
           write(92,"(I5,4E15.7)")Nk,(WDES%VKPT(I,NK),I=1,3),WDES%WTKPT(NK)
           do i=1,NB_TOT
           do j=1,NB_TOT
             write(92,"(2I5,2E15.7)") i,j,nabij(i,j)
           enddo
           enddo
           write(92,*)
        else if(idir==2)then
           write(93,"(I5,4E15.7)")Nk,(WDES%VKPT(I,NK),I=1,3),WDES%WTKPT(NK)
           do i=1,NB_TOT
           do j=1,NB_TOT
             write(93,"(2I5,2E15.7)") i,j,nabij(i,j)
           enddo
           enddo
           write(93,*)
         else
           write(94,"(I5,4E15.7)")Nk,(WDES%VKPT(I,NK),I=1,3),WDES%WTKPT(NK)
           do i=1,NB_TOT
           do j=1,NB_TOT
             write(94,"(2I5,2E15.7)") i,j,nabij(i,j)
           enddo
           enddo
           write(94,*)
        endif
        io_end
!ad-----by fwx-------------
         SCALE_NABIJ_TO_VIJ=AUTOA    ! hbar/m_e * 1/length in whatever unit
         NABIJ=NABIJ*SCALE_NABIJ_TO_VIJ   ! rescaling to bohr radii/Hartree
! write OPTIC
         io_begin
         IF (IU>=0) THEN
           CALL OUTOPT(NABIJ,IDIR,NK,ISP,W%CELTOT,W%FERTOT,NB_TOT, &
                       INFO%ISPIN,KPOINTS%NKPTS,WDES%VKPT, &
                       WDES%WTKPT,IRECLO,NBVAL,NBCON,NKOPT,NKOFF,IU)
         ENDIF
         io_end
        ENDDO    ! IDIR
       ENDDO     ! ISP
      ENDDO      ! NK
      io_begin
      IF (IU>=0) CLOSE(IU)
!ad-----by fwx-------------
      CLOSE(91)
      CLOSE(92)
      CLOSE(93)
      CLOSE(94)
!d----------by fwx
      io_end
But in my output file "Px.dat, Py.dat, Pz.dat", the matrix elements < phi_i | nabla | phi_j > (i.e. the array nabij(:,:) in above code) has not the same valus with two equivalent k point. My KPOINTS is:
Code: Select all
fixed line  mesh
6
Reciprocal lattice
   0.00000000  0.00000000  0.00000000       0.037
   0.00000000  0.00000000  0.00000000       0.037
   0.00000000  0.00000000  0.00000000       0.037
   0.50000000  0.50000000 -0.50000000       0.037
   0.50000000  0.50000000 -0.50000000       0.037
   0.50000000  0.50000000 -0.50000000       0.037
for example, the first and second k point are equivalent, the gamma point, but their "nabij(:,:)" has not smae values.
Another problem, the "nabij(:,:)" has different values in two same calcuation with the same file INCAR. My INCAR is:
Code: Select all
ENCUT=350
#IBRIONÂ =Â Â -1
NSW=0
#ISPIN=2
#MAGMON=0Â 16*5Â 6*0
LSORBIT=T
SAXIS=Â 0Â 0Â 1
ALGO=Fast
ICHARG=11
NSW=0
LWAVE=F
LMAXMIX=4
VOSKOWNÂ Â =Â Â 1
NBANDS=Â Â 30
NELMDL=6
#ISMEAR=0
AMIXÂ =Â 0.2
BMIXÂ =Â 0.0001
AMIX_MAGÂ =Â 0.8
BMIX_MAGÂ =Â 0.0001
#LNABLAÂ =Â .TRUE.
LOPTICSÂ =Â .TRUE.
NPARÂ =Â 1
#MAXMIX=50
Any suggestions will be gratefull!
jz zhao
<span class='smallblacktext'>[ Edited ]</span>