LCOV - code coverage report
Current view: top level - wannier - wann_socmatvec_rs.F (source / functions) Coverage Total Hit
Test: FLEUR test coverage Lines: 0.0 % 168 0
Test Date: 2025-06-14 04:34:23 Functions: 0.0 % 1 0

            Line data    Source code
       1              :       module m_wann_socmatvec_rs
       2              : c      USE m_fleurenv
       3              :       contains 
       4            0 :       subroutine wann_socmatvec_rs(
       5            0 :      >          rvecnum,rvec,kpoints,
       6              :      >          jspins_in,nkpts,l_bzsym,film,
       7              :      >          l_soc,band_min,band_max,neigd,
       8            0 :      >          l_socmmn0,l_ndegen,ndegen,wan90version,
       9              :      >          l_unformatted)
      10              : c*************************************************
      11              : c     Calculate the matrix elements of SOC perturbation 
      12              : c     in real space from the
      13              : c     files WF1.chk (and WF1_um.dat) (produced
      14              : c     by wannier90) and WF1.hsomtx.
      15              : c
      16              : c     Frank Freimuth
      17              : c*************************************************
      18              :       use m_constants, only:pimach
      19              :       use m_wann_read_umatrix
      20              :       use m_juDFT
      21              :       implicit none
      22              :       integer, intent(in) :: rvecnum
      23              :       integer, intent(in) :: rvec(:,:)
      24              :       real,    intent(in) :: kpoints(:,:)
      25              : 
      26              :       integer, intent(in) :: jspins_in
      27              :       integer, intent(in) :: nkpts
      28              :       logical,intent (in) :: l_bzsym,l_soc
      29              :       logical,intent(in)  :: film
      30              :       integer,intent(in)  :: band_min(2),band_max(2),neigd
      31              :       logical, intent(in) :: l_socmmn0
      32              :       logical, intent(in) :: l_ndegen
      33              :       integer, intent(in) :: ndegen(:)
      34              :       integer, intent(in) :: wan90version
      35              :       logical, intent(in) :: l_unformatted
      36              : 
      37              :       integer             :: dir
      38              :       integer             :: ikpt,jspins
      39              :       integer             :: kpts
      40              :       logical             :: l_file
      41              : c      real                :: kpoints(3,nkpts)
      42              :       integer             :: num_wann,num_kpts,num_nnmax,jspin
      43              :       integer             :: kspin,kkspin
      44              :       integer             :: num_wann2
      45              :       integer             :: i,j,k,m,info,r1,r2,r3,dummy1
      46              :       integer             :: dummy2,dummy3,dummy4,dummy5,dummy6
      47              :       integer             :: counter,m1,m2
      48              :       integer             :: num_bands2
      49              :       integer,allocatable :: iwork(:)
      50              :       real,allocatable    :: energy(:,:),ei(:)
      51              :       real,allocatable    :: eigw(:,:),rwork(:)
      52              :       complex,allocatable :: work(:),vec(:,:)
      53            0 :       complex,allocatable :: u_matrix(:,:,:,:)
      54            0 :       complex,allocatable ::    hwann(:,:,:,:,:,:)
      55            0 :       complex,allocatable :: hwannmix(:,:,:,:,:,:)
      56            0 :       complex,allocatable :: hreal(:,:,:,:,:,:)
      57            0 :       complex,allocatable :: hrealunf(:,:,:,:)
      58            0 :       complex,allocatable :: hsomtx(:,:,:,:,:,:)
      59            0 :       complex,allocatable :: hsomtx2(:,:,:,:,:,:)
      60            0 :       complex,allocatable :: hsomtxmix(:,:,:,:,:,:)
      61              :       complex             :: fac,eulav,eulav1
      62              :       real                :: tmp_omi,rdotk,tpi,minenerg,maxenerg
      63              :       real, allocatable   :: minieni(:),maxieni(:)
      64              :       character           :: jobz,uplo
      65              :       integer             :: kpt,band,lee,lwork,lrwork,liwork,n,lda
      66              :       complex             :: value(4)
      67              :       logical             :: um_format
      68              :       logical             :: repro_eig
      69              :       logical             :: l_chk,l_proj
      70              :       logical             :: have_disentangled
      71            0 :       integer,allocatable :: ndimwin(:,:)
      72            0 :       logical,allocatable :: lwindow(:,:,:)
      73              :       integer             :: chk_unit,nkp,ntmp,ierr
      74              :       character(len=33)   :: header
      75              :       character(len=20)   :: checkpoint
      76              :       real                :: tmp_latt(3,3), tmp_kpt_latt(3,nkpts)
      77              :       real                :: omega_invariant
      78            0 :       complex,allocatable :: u_matrix_opt(:,:,:,:)
      79              :       integer             :: num_bands
      80              :       logical             :: l_umdat
      81              :       real,allocatable    :: eigval2(:,:)
      82              :       real,allocatable    :: eigval_opt(:,:)
      83              :       real                :: scale,a,b
      84              :       character(len=2)    :: spinspin12(0:2)
      85              :       character(len=3)    :: spin12(2)
      86              :       character(len=6)    :: filename
      87              :       integer             :: jp,mp,kk,ii,jj,rvecind
      88              :       complex,parameter   :: ci=(0.0,1.0)
      89              :       integer :: spin2_dum,spin1_dum,num_bands1_dum,num_bands2_dum,
      90              :      &               num_dims_dum, fullnkpts_dum  
      91              :       data spinspin12/'  ','.1' , '.2'/
      92              :       data spin12/'WF1','WF2'/
      93            0 :       call timestart("wann_socmatvec_rs")
      94            0 :       tpi=2*pimach()
      95              : 
      96            0 :       jspins=jspins_in
      97            0 :       if(l_soc)jspins=1
      98              : 
      99            0 :       write(6,*)"nkpts=",nkpts
     100              : 
     101              : c*****************************************************
     102              : c     get num_bands and num_wann from the proj file
     103              : c*****************************************************
     104            0 :       do j=1,0,-1
     105            0 :           inquire(file=trim('proj'//spinspin12(j)),exist=l_file)
     106            0 :           if(l_file)then
     107            0 :             filename='proj'//spinspin12(j)
     108            0 :             exit
     109              :           endif
     110              :       enddo
     111            0 :       if(l_file)then
     112            0 :           open (203,file=trim(filename),status='old')
     113            0 :           rewind (203)
     114              :       else
     115              :          CALL juDFT_error("no proj/proj.1/proj.2",calledby
     116            0 :      +        ="wann_socmatvec_rs")
     117              :       endif
     118            0 :       read (203,*) num_wann,num_bands
     119            0 :       close (203)
     120            0 :       write(6,*)'According to proj there are ',num_bands,' bands'
     121            0 :       write(6,*)"and ",num_wann," wannier functions."
     122              : 
     123              : c****************************************************************
     124              : c        read in chk
     125              : c****************************************************************
     126            0 :       num_kpts=nkpts
     127            0 :       allocate( u_matrix_opt(num_bands,num_wann,nkpts,2) )
     128            0 :       allocate( u_matrix(num_wann,num_wann,nkpts,2) )
     129            0 :       allocate( lwindow(num_bands,nkpts,2) )
     130            0 :       allocate( ndimwin(nkpts,2) )
     131              : 
     132            0 :       do jspin=1,jspins  !spin loop
     133              :          call wann_read_umatrix2(
     134              :      >       nkpts,num_wann,num_bands,
     135              :      >       um_format,jspin,wan90version,
     136              :      <       have_disentangled,
     137              :      <       lwindow(:,:,jspin),
     138              :      <       ndimwin(:,jspin),
     139              :      <       u_matrix_opt(:,:,:,jspin),
     140            0 :      <       u_matrix(:,:,:,jspin))
     141            0 :          num_bands2=num_bands
     142              :       enddo !jspin   
     143            0 :       if(jspins.eq.1)then
     144            0 :          lwindow(:,:,2)        = lwindow(:,:,1)
     145            0 :          ndimwin(:,2)          = ndimwin(:,1)
     146            0 :          u_matrix_opt(:,:,:,2) = u_matrix_opt(:,:,:,1)
     147            0 :          u_matrix(:,:,:,2)     = u_matrix(:,:,:,1)
     148              :       endif
     149              : 
     150              : c****************************************************
     151              : c     Read the file "WF1.hsomtxvec".
     152              : c**************************************************** 
     153            0 :       allocate( hsomtx(2,2,num_bands2,num_bands2,3,nkpts) )
     154            0 :       if(l_unformatted)then
     155              : 
     156            0 :         open(304,file='WF1.hsomtxvec_unf',form='unformatted')
     157            0 :         read(304)spin2_dum,spin1_dum,num_bands1_dum,num_bands2_dum,
     158            0 :      &                   num_dims_dum, fullnkpts_dum           
     159            0 :         do nkp=1,num_kpts
     160            0 :           read(304)hsomtx(:,:,:,:,:,nkp) 
     161              :         enddo
     162            0 :         close(304)
     163            0 :         hsomtx=conjg(hsomtx)
     164              :       else
     165              : 
     166            0 :        open(304,file='WF1.hsomtxvec',form='formatted')
     167              : 
     168            0 :        read(304,*)
     169            0 :        read(304,*)
     170              : 
     171            0 :        do nkp=1,num_kpts
     172            0 :         do dir=1,3  
     173            0 :          do i=1,num_bands2
     174            0 :           do j=1,num_bands2
     175            0 :            do ii=1,2
     176            0 :             do jj=1,2
     177            0 :              read(304,*)dummy1,dummy2,dummy3,dummy4,dummy5,dummy6,a,b
     178            0 :              hsomtx(jj,ii,j,i,dir,nkp)=cmplx(a,-b)
     179              :             enddo !jj
     180              :            enddo !ii 
     181              :           enddo !j
     182              :          enddo !i
     183              :         enddo !dir
     184              :        enddo !nkp
     185            0 :        close(304)
     186              : 
     187              :       endif
     188              : c****************************************************************
     189              : c        Calculate matrix elements of SOC in the basis of
     190              : c        rotated Bloch functions.
     191              : c****************************************************************
     192            0 :       allocate( hsomtx2(2,2,num_wann,num_wann,3,nkpts) )
     193              : 
     194              :       write(6,*)"calculate matrix elements of SOC operator
     195            0 :      &between wannier orbitals"
     196              : 
     197            0 :       if(have_disentangled) then       
     198            0 :        allocate( hsomtxmix(2,2,num_wann,num_bands2,3,nkpts) )
     199            0 :        hsomtx2=0.0  
     200            0 :        hsomtxmix=0.0
     201            0 :        do nkp=1,num_kpts
     202            0 :         do dir=1,3  
     203            0 :           do jp=1,num_wann  
     204            0 :            do ii=1,2
     205            0 :             do jj=1,2
     206            0 :              do m=1,ndimwin(nkp,ii)
     207            0 :               do mp=1,ndimwin(nkp,jj)  
     208              :                   hsomtxmix(jj,ii,jp,m,dir,nkp)=
     209              :      &            hsomtxmix(jj,ii,jp,m,dir,nkp)+ 
     210              :      &            conjg(u_matrix_opt(mp,jp,nkp,jj))*
     211            0 :      &                  hsomtx(jj,ii,mp,m,dir,nkp)
     212              :               enddo !mp   
     213              :              enddo !m
     214              :             enddo !jj  
     215              :            enddo !ii
     216              :           enddo !jp 
     217              :         enddo !dir
     218              :        enddo !nkp
     219              : 
     220            0 :        do nkp=1,num_kpts
     221            0 :         do dir=1,3  
     222            0 :          do j=1,num_wann
     223            0 :           do jp=1,num_wann  
     224            0 :            do ii=1,2
     225            0 :             do m=1,ndimwin(nkp,ii)
     226            0 :              do jj=1,2
     227              :               hsomtx2(jj,ii,jp,j,dir,nkp)=
     228              :      &        hsomtx2(jj,ii,jp,j,dir,nkp)+ 
     229              :      &                  hsomtxmix(jj,ii,jp,m,dir,nkp)*
     230            0 :      &                  u_matrix_opt(m,j,nkp,ii)
     231              :              enddo !jj  
     232              :             enddo !m
     233              :            enddo !ii
     234              :           enddo !jp 
     235              :          enddo !j
     236              :         enddo !dir
     237              :        enddo !nkp
     238              : 
     239            0 :        deallocate( hsomtxmix )
     240              :       else
     241            0 :        hsomtx2 = hsomtx
     242              :       end if !have_disentangled
     243              : 
     244            0 :       allocate(hwann(2,2,num_wann,num_wann,3,num_kpts))
     245            0 :       hwann=cmplx(0.0,0.0)
     246              : 
     247            0 :       allocate(hwannmix(2,2,num_wann,num_wann,3,num_kpts))
     248            0 :       hwannmix=cmplx(0.0,0.0)
     249              : 
     250            0 :       do k=1,num_kpts
     251            0 :        do dir=1,3  
     252            0 :         do i=1,num_wann
     253            0 :          do mp=1,num_wann
     254            0 :           do j=1,num_wann
     255            0 :            do ii=1,2
     256            0 :             do jj=1,2
     257              : 
     258              :           hwannmix(jj,ii,mp,i,dir,k)=hwannmix(jj,ii,mp,i,dir,k)+
     259              : 
     260              :      *                hsomtx2(jj,ii,j,i,dir,k)*
     261            0 :      *        conjg(u_matrix(j,mp,k,jj))
     262              :              enddo !j
     263              :             enddo !i     
     264              :            enddo !jj 
     265              :           enddo !ii
     266              :          enddo !mp
     267              :        enddo !dir 
     268              :       enddo !k
     269              : 
     270              : c$$$      do k=1,num_kpts
     271              : c$$$       do dir=1,3  
     272              : c$$$        do m=1,num_wann
     273              : c$$$          do ii=1,2
     274              : c$$$           do jj=1,2
     275              : c$$$            do i=1,num_wann
     276              : c$$$             do j=1,num_wann
     277              : c$$$          hwannmix(jj,ii,j,m,dir,k)=hwannmix(jj,ii,j,m,dir,k)+
     278              : c$$$     *                hsomtx2(jj,ii,i,m,dir,k)*
     279              : c$$$     *              conjg(u_matrix(i,j,k,jj))
     280              : c$$$             enddo !j
     281              : c$$$            enddo !i     
     282              : c$$$           enddo !jj 
     283              : c$$$          enddo !ii
     284              : c$$$        enddo !m
     285              : c$$$       enddo !dir 
     286              : c$$$      enddo !k
     287              : 
     288              : 
     289            0 :       do k=1,num_kpts
     290            0 :        do dir=1,3  
     291            0 :         do m=1,num_wann
     292            0 :           do ii=1,2
     293            0 :            do jj=1,2
     294            0 :             do i=1,num_wann
     295            0 :              do j=1,num_wann
     296              :           hwann(jj,ii,j,m,dir,k)=hwann(jj,ii,j,m,dir,k)+
     297              :      *                hwannmix(jj,ii,j,i,dir,k)*
     298            0 :      *              u_matrix(i,m,k,ii)
     299              :              enddo !j
     300              :             enddo !i     
     301              :            enddo !jj 
     302              :           enddo !ii
     303              :         enddo !m
     304              :        enddo !dir 
     305              :       enddo !k
     306              : 
     307            0 :       deallocate(hwannmix)
     308              : c************************************************************
     309              : c        Calculate matrix elements in real space.
     310              : c***********************************************************      
     311            0 :       write(6,*)"calculate SOC-mat in rs"
     312              : 
     313            0 :       allocate(hreal(2,2,num_wann,num_wann,3,rvecnum))
     314            0 :       hreal=cmplx(0.0,0.0)
     315            0 :       if(l_ndegen)then
     316            0 :        do rvecind=1,rvecnum
     317            0 :         do k=1,nkpts  
     318              :          rdotk=tpi*(  kpoints(1,k)*rvec(1,rvecind)+
     319              :      +               kpoints(2,k)*rvec(2,rvecind)+
     320            0 :      +               kpoints(3,k)*rvec(3,rvecind)  )
     321            0 :          fac=cmplx(cos(rdotk),-sin(rdotk))/real(ndegen(rvecind))
     322            0 :          do dir=1,3
     323            0 :           do m2=1,num_wann
     324            0 :            do m1=1,num_wann
     325            0 :             do ii=1,2
     326            0 :              do jj=1,2
     327              :             hreal(jj,ii,m1,m2,dir,rvecind)=
     328              :      &      hreal(jj,ii,m1,m2,dir,rvecind)+
     329            0 :      &  fac*hwann(jj,ii,m1,m2,dir,k)
     330              :              enddo  !jj
     331              :             enddo !ii
     332              :            enddo !m1
     333              :           enddo !m2 
     334              :          enddo !dir 
     335              :         enddo !k
     336              :        enddo !rvecind
     337              :       else
     338            0 :        do rvecind=1,rvecnum
     339            0 :         do k=1,nkpts  
     340              :          rdotk=tpi*(  kpoints(1,k)*rvec(1,rvecind)+
     341              :      +               kpoints(2,k)*rvec(2,rvecind)+
     342            0 :      +               kpoints(3,k)*rvec(3,rvecind)  )
     343            0 :          fac=cmplx(cos(rdotk),-sin(rdotk))
     344            0 :          do dir=1,3
     345            0 :           do m2=1,num_wann
     346            0 :            do m1=1,num_wann
     347            0 :             do ii=1,2
     348            0 :              do jj=1,2
     349              :             hreal(jj,ii,m1,m2,dir,rvecind)=
     350              :      &      hreal(jj,ii,m1,m2,dir,rvecind)+
     351            0 :      &  fac*hwann(jj,ii,m1,m2,dir,k)
     352              :              enddo  !jj
     353              :             enddo !ii
     354              :            enddo !m1
     355              :           enddo !m2 
     356              :          enddo !dir 
     357              :         enddo !k
     358              :        enddo !rvecind
     359              :       endif
     360            0 :       hreal=hreal/cmplx(real(nkpts),0.0)
     361              : 
     362            0 :       if(l_unformatted)then
     363            0 :        if(l_ndegen)then  
     364              :          open(321,file='socmunfndegen',form='unformatted'
     365              : #ifdef CPP_INTEL     
     366              :      &                     ,convert='BIG_ENDIAN')  
     367              : #else
     368            0 :      &    )
     369              : #endif        
     370              :        else
     371              :          open(321,file='socmunf',form='unformatted'
     372              : #ifdef CPP_INTEL           
     373              :      &                     ,convert='BIG_ENDIAN')
     374              : #else
     375            0 :      &   )     
     376              : #endif     
     377              :        endif
     378              : 
     379            0 :        allocate(hrealunf(2*num_wann,2*num_wann,3,rvecnum))
     380            0 :        do rvecind=1,rvecnum
     381            0 :         do dir=1,3
     382            0 :          do ii=1,2
     383            0 :           do jj=1,2
     384            0 :            do m1=1,num_wann
     385            0 :             do m2=1,num_wann
     386              :              hrealunf(m2+num_wann*(jj-1),m1+num_wann*(ii-1),dir,rvecind)
     387            0 :      &  =hreal(jj,ii,m2,m1,dir,rvecind)
     388              :             enddo
     389              :            enddo
     390              :           enddo
     391              :          enddo
     392              :         enddo 
     393              :        enddo
     394              :  
     395              : 
     396            0 :        do dir=1,3
     397            0 :         do rvecind=1,rvecnum
     398            0 :          write(321)hrealunf(:,:,dir,rvecind) 
     399              :         enddo
     400              :        enddo
     401            0 :        deallocate(hrealunf)
     402              :       else
     403            0 :        open(321,file='rssocmatvec.1',form='formatted')
     404              : 
     405            0 :        do rvecind=1,rvecnum
     406            0 :          r3=rvec(3,rvecind)
     407            0 :          r2=rvec(2,rvecind)
     408            0 :          r1=rvec(1,rvecind)
     409            0 :          do i=1,num_wann
     410            0 :            do j=1,num_wann
     411            0 :             do ii=1,2
     412            0 :              do jj=1,2  
     413            0 :               do dir=1,3  
     414              :               write(321,'(i3,1x,i3,1x,i3,1x,i3,1x,i3,1x,
     415              :      &                      i3,1x,i3,1x,i3,f20.8,1x,f20.8)')
     416            0 :      &     r1,r2,r3,i,j,jj,ii,dir,hreal(jj,ii,i,j,dir,rvecind) 
     417              :               enddo
     418              :              enddo
     419              :             enddo !kk  
     420              :            enddo !j
     421              :          enddo !i
     422              :        enddo !rvecind    
     423              :       endif 
     424            0 :       close(321)
     425              : 
     426            0 :       deallocate(lwindow,u_matrix_opt,ndimwin)
     427            0 :       deallocate(u_matrix,hwann,hreal)
     428            0 :       call timestop("wann_socmatvec_rs")
     429            0 :       end subroutine wann_socmatvec_rs
     430              :       end module m_wann_socmatvec_rs
     431              : 
     432              :  
        

Generated by: LCOV version 2.0-1