calc_gwdos.f90

Go to the documentation of this file.
00001 subroutine calc_gwdos(NWF,NTK,nsgm,Na1,Na2,Na3,nkb1,nkb2,nkb3,idlt,dmna,dmnr,FermiEnergy,a1,a2,a3,b1,b2,b3,sgmw,SK0,&
00002   KS_R,XC_R,SX_R,SC_R,shift_value,gwdos,gw_sigma_dos) 
00003   !
00004   implicit none 
00005   integer::NWF,NTK,nsgm,Na1,Na2,Na3,nkb1,nkb2,nkb3 
00006   real(8)::idlt,dmna,dmnr
00007   real(8)::a1(3),a2(3),a3(3)
00008   real(8)::b1(3),b2(3),b3(3)
00009   real(8)::sgmw(nsgm)  
00010   real(8)::SK0(3,NTK) 
00011   complex(8)::KS_R(NWF,NWF,-Na1:Na1,-Na2:Na2,-Na3:Na3)   
00012   complex(8)::XC_R(NWF,NWF,-Na1:Na1,-Na2:Na2,-Na3:Na3)   
00013   complex(8)::SX_R(NWF,NWF,-Na1:Na1,-Na2:Na2,-Na3:Na3)   
00014   complex(4)::SC_R(NWF,NWF,-Na1:Na1,-Na2:Na2,-Na3:Na3,nsgm)   
00015   !
00016   real(8),allocatable::WEIGHT_R(:,:,:)!WEIGHT_R(-Na1:Na1,-Na2:Na2,-Na3:Na3)
00017   real(8),allocatable::EKS(:,:)!EKS(NWF,NTK)           
00018   complex(8),allocatable::VKS(:,:,:)!VKS(NWF,NWF,NTK)   
00019   complex(8),allocatable::pf(:,:,:,:)!pf(-Na1:Na1,-Na2:Na2,-Na3:Na3,NTK)   
00020   complex(8),allocatable::EMK(:,:,:)!EMK(NWF,NTK,nsgm)           
00021   complex(4),allocatable::GW_R(:,:,:,:,:,:)!GW_R(NWF,NWF,-Na1:Na1,-Na2:Na2,-Na3:Na3,nsgm)   
00022   ! 
00023   integer::ia1min,ia2min,ia3min 
00024   integer::ia1,ia2,ia3,ik,ib,jb,ie  
00025   real(8)::PHASE,FermiEnergy 
00026   real(8)::SUM_REAL
00027   complex(8)::SUM_CMPX 
00028   !
00029   real(8),parameter::au=27.21151d0
00030   real(8),parameter::tpi=2.0d0*dacos(-1.0d0)
00031   complex(8),parameter::ci=(0.0D0,1.0D0) 
00032   !
00033   real(8),intent(out)::shift_value 
00034   real(8),intent(out)::gwdos(nsgm) 
00035   complex(8),intent(out)::gw_sigma_dos(nsgm) 
00036   !
00037   !1. make GR_R 
00038   !
00039   allocate(GW_R(NWF,NWF,-Na1:Na1,-Na2:Na2,-Na3:Na3,nsgm));GW_R=0.0d0   
00040   do ie=1,nsgm
00041    do ia1=-Na1,Na1
00042     do ia2=-Na2,Na2
00043      do ia3=-Na3,Na3
00044       do ib=1,NWF          
00045        do jb=1,NWF          
00046         GW_R(ib,jb,ia1,ia2,ia3,ie)&
00047        =KS_R(ib,jb,ia1,ia2,ia3)&
00048        -XC_R(ib,jb,ia1,ia2,ia3)&
00049        -SX_R(ib,jb,ia1,ia2,ia3)&
00050        +SC_R(ib,jb,ia1,ia2,ia3,ie) 
00051        enddo 
00052       enddo 
00053      enddo 
00054     enddo 
00055    enddo 
00056   enddo 
00057   write(6,*)'# finish make GW_R'
00058   !
00059   !2. WEIGHT_R BY Y.Nomura NOMURA 
00060   !
00061   allocate(WEIGHT_R(-Na1:Na1,-Na2:Na2,-Na3:Na3)); WEIGHT_R=1.0d0
00062   SUM_REAL=0.0d0 
00063   do ia1=-Na1,Na1
00064    do ia2=-Na2,Na2
00065     do ia3=-Na3,Na3
00066      if(abs(ia1)==Na1.and.mod(nkb1,2)==0.and.Na1/=0) WEIGHT_R(ia1,ia2,ia3)=WEIGHT_R(ia1,ia2,ia3)*0.5d0 
00067      if(abs(ia2)==Na2.and.mod(nkb2,2)==0.and.Na2/=0) WEIGHT_R(ia1,ia2,ia3)=WEIGHT_R(ia1,ia2,ia3)*0.5d0 
00068      if(abs(ia3)==Na3.and.mod(nkb3,2)==0.and.Na3/=0) WEIGHT_R(ia1,ia2,ia3)=WEIGHT_R(ia1,ia2,ia3)*0.5d0 
00069      SUM_REAL=SUM_REAL+WEIGHT_R(ia1,ia2,ia3)
00070     enddo!ia3
00071    enddo!ia2
00072   enddo!ia1
00073   write(6,'(a20,f15.8,i8)')'SUM_WEIGHT,NTK',SUM_REAL,NTK  
00074   if(abs(SUM_REAL-dble(NTK))>1.0d-6)then 
00075    stop 'SUM_WEIGHT/=NTK'
00076   endif 
00077   !
00078   allocate(pf(-Na1:Na1,-Na2:Na2,-Na3:Na3,NTK)); pf=0.0d0 
00079   do ik=1,NTK 
00080    do ia3=-Na3,Na3 
00081     do ia2=-Na2,Na2 
00082      do ia1=-Na1,Na1 
00083       !
00084       !NEAREST R SEARCH BY Y.Yoshimoto 
00085       !
00086       call search_Rmin(ia1,ia2,ia3,nkb1,nkb2,nkb3,a1(1),a2(1),a3(1),ia1min,ia2min,ia3min)
00087       PHASE=tpi*(SK0(1,ik)*DBLE(ia1min)+SK0(2,ik)*DBLE(ia2min)+SK0(3,ik)*DBLE(ia3min))           
00088       pf(ia1,ia2,ia3,ik)=EXP(ci*PHASE)*WEIGHT_R(ia1,ia2,ia3)          
00089      enddo!ia1 
00090     enddo!ia2 
00091    enddo!ia3 
00092   enddo!ik  
00093   write(6,*)'# finish make pf'
00094   !
00095   !3. HKS IN WANNIER BASIS. AND DIAGONALIZE
00096   !
00097   allocate(EKS(NWF,NTK)); EKS=0.0d0 
00098   allocate(VKS(NWF,NWF,NTK)); VKS=0.0d0 
00099   call make_eks(NTK,NWF,Na1,Na2,Na3,KS_R(1,1,-Na1,-Na2,-Na3),pf(-Na1,-Na2,-Na3,1),EKS(1,1),VKS(1,1,1))  
00100   !
00101   !4. HGW IN WANNIER BASIS. AND DIAGONALIZE
00102   !
00103   allocate(EMK(NWF,NTK,nsgm)); EMK=0.0d0 
00104   call make_emk(NTK,NWF,nsgm,Na1,Na2,Na3,GW_R(1,1,-Na1,-Na2,-Na3,1),pf(-Na1,-Na2,-Na3,1),EMK(1,1,1))  
00105   !
00106   !5. DETERMINE SHIFT 
00107   !
00108   call det_shift(NTK,NWF,nsgm,FermiEnergy,sgmw(1),EKS(1,1),EMK(1,1,1),shift_value) 
00109   !
00110   !6. GW-DOS CALC
00111   !
00112   gwdos=0.0d0  
00113   call calc_dos_GW(nkb1,nkb2,nkb3,NTK,NWF,nsgm,idlt,dmna,dmnr,shift_value,b1(1),b2(1),b3(1),sgmw(1),SK0(1,1),EMK(1,1,1),gwdos(1)) 
00114   ! 
00115   !7. GW-SIGMA-DOS CALC
00116   !
00117   gw_sigma_dos=0.0d0  
00118   do ie=1,nsgm 
00119    SUM_CMPX=0.0d0 
00120    do ik=1,NTK 
00121     do ib=1,NWF
00122      SUM_CMPX=SUM_CMPX+(EMK(ib,ik,ie)-cmplx(EKS(ib,ik)))   
00123     enddo!ib
00124    enddo!ik
00125    gw_sigma_dos(ie)=SUM_CMPX/dble(NTK) 
00126   enddo!ie 
00127   !
00128   deallocate(WEIGHT_R,EKS,VKS,pf,EMK,GW_R) 
00129   !
00130 return 
00131 end