m_rd_dat_zvo.f90

Go to the documentation of this file.
00001 module m_rd_dat_zvo 
00002   implicit none
00003   private 
00004   public::rd_dat_hr
00005   public::rd_dat_geom 
00006   public::rd_dat_bandkpts
00007   public::rd_dat_mkkpts 
00008   public::rd_dat_ef 
00009   !h_mat_r(300) 
00010   integer,public::NWF 
00011   complex(8),public,allocatable::HR(:,:,:,:,:)!HR(NWF,NWF,-Na1:Na1,-Na2:Na2,-Na3:Na3) 
00012   !geom(303)
00013   integer,public::N_wannier_center 
00014   real(8),public::a1(3),a2(3),a3(3) 
00015   real(8),public::b1(3),b2(3),b3(3)
00016   real(8),public::VOLUME 
00017   real(8),public,allocatable::wcenter_lat(:,:)!wcenter_lat(3,N_wannier_center) 
00018   !bandkpts(304)
00019   integer,public::Ndiv,N_sym_points,NSK_BAND_DISP
00020   real(8),public,allocatable::SK_BAND_DISP(:,:)!SK_BAND_DISP(3,NSK_BAND_DISP) 
00021   !mkkpts(305)  
00022   integer,public::NTK 
00023   integer,public::nkb1,nkb2,nkb3 
00024   integer,public::Na1,Na2,Na3 
00025   real(8),public,allocatable::SK0(:,:)!SK0(3,NTK) 
00026   !ef(306) 
00027   real(8),public::FermiEnergy_bandcalc 
00028   contains
00029   !--
00030   subroutine rd_dat_hr 
00031     implicit none 
00032     integer::dum1,dum2,dum3,dumi,dumj
00033     integer::N_element 
00034     character(len=80)::dum_ch
00035     integer,allocatable::unit_vec(:)!unit_vec(N_element) 
00036     integer::i,ia1,ia2,ia3,ib,jb 
00037     real(8),parameter::au=27.21151d0
00038     !
00039     !OPEN(300,R,FILE='./dir-model/zvo_hr.dat') 
00040     !
00041     OPEN(300,FILE='./dir-model/zvo_hr.dat') 
00042     rewind(300) 
00043     read(300,'(a)') dum_ch 
00044     read(300,'(i10)') NWF 
00045     read(300,'(i10)') N_element 
00046     allocate(unit_vec(N_element)); unit_vec=1
00047     read(300,'(15i5)')(unit_vec(i),i=1,N_element) 
00048     deallocate(unit_vec) 
00049     allocate(HR(NWF,NWF,-Na1:Na1,-Na2:Na2,-Na3:Na3));HR(:,:,:,:,:)=0.0d0  
00050     do ia1=-Na1,Na1
00051      do ia2=-Na2,Na2
00052       do ia3=-Na3,Na3 
00053        do ib=1,NWF!n_occ
00054         do jb=1,NWF!n_occ
00055          read(300,'(i5,i5,i5,i5,i5,2f15.10)') dum1,dum2,dum3,dumi,dumj,HR(ib,jb,ia1,ia2,ia3) 
00056         enddo!jb
00057        enddo!ib
00058       enddo!ia3
00059      enddo!ia2
00060     enddo!ia1 
00061     HR=HR/au !au<-eV
00062   end subroutine
00063   ! 
00064   subroutine rd_dat_geom 
00065     implicit none 
00066     real(8),parameter::tpi=2.0d0*dacos(-1.0d0)
00067     integer::ib,i 
00068     !
00069     !OPEN(303,R,FILE='./dir-model/zvo_geom.dat') 
00070     !
00071     OPEN(303,FILE='./dir-model/zvo_geom.dat') 
00072     read(303,'(3f15.10)')(a1(i),i=1,3) 
00073     read(303,'(3f15.10)')(a2(i),i=1,3) 
00074     read(303,'(3f15.10)')(a3(i),i=1,3) 
00075     read(303,'(i10)') N_wannier_center 
00076     allocate(wcenter_lat(3,N_wannier_center)) 
00077     do ib=1,N_wannier_center 
00078      read(303,*)(wcenter_lat(i,ib),i=1,3) 
00079     enddo 
00080     !
00081     call OUTER_PRODUCT(a2(1),a3(1),b1(1))
00082     VOLUME=a1(1)*b1(1)+a1(2)*b1(2)+a1(3)*b1(3)
00083     b1(:)=b1(:)*tpi/VOLUME 
00084     call OUTER_PRODUCT(a3(1),a1(1),b2(1))
00085     b2(:)=b2(:)*tpi/VOLUME 
00086     call OUTER_PRODUCT(a1(1),a2(1),b3(1))
00087     b3(:)=b3(:)*tpi/VOLUME 
00088   end subroutine
00089   !
00090   subroutine rd_dat_bandkpts 
00091     implicit none 
00092     integer::ik,i 
00093     !
00094     !OPEN(304,R,FILE='zvo_bandkpts.dat') 
00095     !
00096     OPEN(304,FILE='./dir-model/zvo_bandkpts.dat') 
00097     read(304,'(3i10)') NSK_BAND_DISP,Ndiv,N_sym_points 
00098     allocate(SK_BAND_DISP(3,NSK_BAND_DISP));SK_BAND_DISP=0.0d0 
00099     do ik=1,NSK_BAND_DISP
00100      read(304,*)(SK_BAND_DISP(i,ik),i=1,3) 
00101     enddo 
00102   end subroutine
00103   !
00104   subroutine rd_dat_mkkpts 
00105     implicit none 
00106     integer::ik,i 
00107     !
00108     !OPEN(305,R,FILE='zvo_mkkpts.dat') 
00109     !
00110     OPEN(305,FILE='./dir-model/zvo_mkkpts.dat') 
00111     read(305,'(i10)') NTK 
00112     allocate(SK0(3,NTK));SK0=0.0d0 
00113     do ik=1,NTK 
00114      read(305,*)(SK0(i,ik),i=1,3) 
00115     enddo 
00116     call est_nkbi(NTK,SK0,nkb1,nkb2,nkb3)  
00117     Na1=nkb1/2; Na2=nkb2/2; Na3=nkb3/2
00118   end subroutine
00119   !
00120   subroutine rd_dat_ef 
00121     implicit none 
00122     !
00123     !OPEN(306,R,FILE='zvo_ef.dat') 
00124     !
00125     OPEN(306,FILE='./dir-model/zvo_ef.dat') 
00126     read(306,'(f15.10)') FermiEnergy_bandcalc 
00127   end subroutine rd_dat_ef 
00128   !
00129   subroutine OUTER_PRODUCT(vec_x,vec_y,vec_z)
00130     implicit none 
00131     real(8)::vec_x(3),vec_y(3),vec_z(3) 
00132     !
00133     vec_z(1)=vec_x(2)*vec_y(3)-vec_x(3)*vec_y(2)
00134     vec_z(2)=vec_x(3)*vec_y(1)-vec_x(1)*vec_y(3) 
00135     vec_z(3)=vec_x(1)*vec_y(2)-vec_x(2)*vec_y(1)
00136     !
00137     return
00138   end subroutine 
00139   !
00140   subroutine est_nkbi(N,SK,nkb1,nkb2,nkb3)  
00141     implicit none 
00142     integer,intent(in)::N
00143     real(8),intent(in)::SK(3,N) 
00144     integer,intent(out)::nkb1,nkb2,nkb3
00145     integer::NTK  
00146     integer::i 
00147     real(8)::x 
00148     real(8),parameter::dlt_BZ=1.0d-6 
00149     x=1.0d0 
00150     do i=1,N
00151      if(abs(SK(1,i))<dlt_BZ) cycle 
00152      if(abs(SK(1,i))<x) then 
00153       x=abs(SK(1,i))  
00154      endif 
00155     enddo    
00156     nkb1=nint(1.0d0/x)  
00157     x=1.0d0 
00158     do i=1,N
00159      if(abs(SK(2,i))<dlt_BZ) cycle 
00160      if(abs(SK(2,i))<x) then 
00161       x=abs(SK(2,i))  
00162      endif 
00163     enddo    
00164     nkb2=nint(1.0d0/x)  
00165     x=1.0d0 
00166     do i=1,N
00167      if(abs(SK(3,i))<dlt_BZ) cycle 
00168      if(abs(SK(3,i))<x) then 
00169       x=abs(SK(3,i))  
00170      endif 
00171     enddo    
00172     nkb3=nint(1.0d0/x)  
00173     NTK=nkb1*nkb2*nkb3 
00174     write(6,'(a50,4i10)')'nkb1, nkb2, nkb3, NTK:',nkb1,nkb2,nkb3,NTK  
00175     return 
00176   end subroutine 
00177   !
00178 end module 

Generated on 17 Nov 2020 for respack by  doxygen 1.6.1