LCOV - code coverage report
Current view: top level - wannier - wann_kptsreduc.f (source / functions) Hit Total Coverage
Test: combined.info Lines: 0 110 0.0 %
Date: 2019-09-08 04:53:50 Functions: 0 1 0.0 %

          Line data    Source code
       1             : !--------------------------------------------------------------------------------
       2             : ! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
       3             : ! This file is part of FLEUR and available as free software under the conditions
       4             : ! of the MIT license as expressed in the LICENSE file in more detail.
       5             : !--------------------------------------------------------------------------------
       6             : 
       7             :       module m_wann_kptsreduc
       8             :       use m_juDFT
       9             :       contains
      10           0 :       subroutine wann_kptsreduc(
      11           0 :      >               nop,mrot,bmat,tau,film,
      12             :      >               l_onedimens,l_nocosoc)
      13             : c*****************************************************************
      14             : c     Apply the symmetries to reduce the number of k-points.
      15             : c     Frank Freimuth
      16             : c*****************************************************************
      17             :       implicit none
      18             :       integer, intent(in) :: nop
      19             :       logical,intent(in)  :: film
      20             :       real, intent(in)    :: bmat(3,3)
      21             :       real,intent(in)     :: tau(3,nop)
      22             :       integer, intent(in) :: mrot(3,3,nop)
      23             :       logical,intent(in)  :: l_onedimens
      24             :       logical,intent(in)  :: l_nocosoc
      25             : 
      26           0 :       real,allocatable    :: weight(:)
      27             :       integer             :: reduznumk
      28           0 :       integer,allocatable :: mapk(:),mapkoper(:)
      29             :       integer             :: ikpt,nmat,k,j
      30             :       real                :: bbkpt(3)
      31             :       real                :: bkpt(3),wk,brot(3)
      32           0 :       real,allocatable    :: kpoints(:,:)
      33             :       integer             :: oper,iter
      34             :       logical             :: l_file
      35             :       real                :: scale
      36             :       integer             :: nkpts
      37             :       integer             :: sumweights
      38           0 :       integer,allocatable :: irreduc(:)
      39           0 :       integer,allocatable :: shiftkpt(:,:)
      40             :       real                :: bbmat(3,3)
      41             :       real,allocatable    :: kptslen(:)
      42             :       logical             :: l_onlysymor
      43             : 
      44           0 :       inquire(file='onlysymor',exist=l_onlysymor)
      45             : 
      46           0 :       bbmat=matmul(bmat,transpose(bmat))
      47           0 :       write(6,*) "Apply the symmetries to w90kpts"
      48             : c**********************************************************
      49             : c     read in kpoints from w90kpts file
      50             : c**********************************************************
      51           0 :       l_file=.false.
      52           0 :       inquire(file='w90kpts',exist=l_file)
      53           0 :       IF(.NOT.l_file) CALL juDFT_error("where is w90kpts?",calledby
      54           0 :      +     ="wann_kptsreduc")
      55           0 :       open(987,file='w90kpts',status='old',form='formatted')
      56           0 :       read(987,*)nkpts,scale
      57           0 :       print*,"nkpts=",nkpts
      58           0 :       allocate(mapk(nkpts),mapkoper(nkpts),kpoints(3,nkpts))
      59           0 :       allocate(weight(nkpts),irreduc(nkpts))
      60           0 :       allocate(shiftkpt(3,nkpts))
      61           0 :       do iter=1,nkpts
      62           0 :          read(987,*)kpoints(:,iter)
      63             :       enddo
      64           0 :       close(987)
      65           0 :       write(6,*) "kpoints read from w90kpts:"
      66           0 :       do iter=1,nkpts
      67           0 :          write(6,'(3f10.5)')kpoints(:,iter)/scale
      68             :       enddo
      69             : 
      70             : c*********************************************************
      71             : c     determine lengths of kpoints
      72             : c*********************************************************
      73             : c      allocate( kptslen(nkpts) )
      74             : c      do iter=1,nkpts
      75             : c         kptslen(iter)=dot_product( axmintx(kpoints(:,iter)),
      76             : c     &                matmul(bbmat,axmintx(kpoints(:,iter))) )
      77             : c      enddo
      78             : 
      79             : c***********************************************
      80             : c            determine irreducible part
      81             : c***********************************************
      82           0 :       weight(:)=1.0
      83           0 :       mapk(:)=0
      84           0 :       reduznumk=0
      85           0 :       mapkoper(:)=0
      86           0 :       irreduc(:)=0
      87           0 :       shiftkpt(:,:)=0
      88           0 :       do ikpt=1,nkpts
      89           0 :          if(mapk(ikpt).ne.0) cycle
      90           0 :          reduznumk=reduznumk+1
      91           0 :          irreduc(ikpt)=reduznumk
      92           0 :          bkpt(:)=kpoints(:,ikpt)         
      93           0 :          do oper=1,nop
      94           0 :            if (all(abs(tau(:,oper)).lt.1e-10).or..not.l_onlysymor)then
      95           0 :              brot(:)=0.0
      96           0 :              do k=1,3
      97           0 :                brot(:)=brot(:)+mrot(k,:,oper)*bkpt(k)
      98             :              enddo
      99           0 :              do j=ikpt+1,nkpts
     100           0 :               if(mapk(j).ne.0)cycle
     101             : c              if(abs(kptslen(j)-kptslen(ikpt)).gt.1.e-8)cycle
     102           0 :               if( all( axmintx( (brot(:)-kpoints(:,j))/scale ) 
     103           0 :      &                            .lt.1e-6) ) then
     104           0 :                     weight(ikpt)=weight(ikpt)+1.0
     105           0 :                     mapk(j)=ikpt
     106           0 :                     mapkoper(j)=oper
     107             :                     shiftkpt(:,j)=
     108           0 :      &                 nint( (brot(:)-kpoints(:,j))/scale )
     109             :               endif
     110             :              enddo
     111             :            endif
     112             :          enddo
     113           0 :          if(.not.l_nocosoc)then
     114           0 :           do oper=1,nop
     115           0 :            if (all(abs(tau(:,oper)).lt.1e-10).or..not.l_onlysymor)then
     116           0 :              brot(:)=0.0
     117           0 :              do k=1,3
     118           0 :                brot(:)=brot(:)-mrot(k,:,oper)*bkpt(k)
     119             :              enddo
     120           0 :              do j=ikpt+1,nkpts
     121           0 :               if(mapk(j).ne.0)cycle
     122             : c              if(abs(kptslen(j)-kptslen(ikpt)).gt.1.e-8)cycle
     123           0 :               if( all( axmintx( (brot(:)-kpoints(:,j))/scale) 
     124           0 :      &                            .lt.1e-6) ) then
     125           0 :                     weight(ikpt)=weight(ikpt)+1.0
     126           0 :                     mapk(j)=ikpt
     127           0 :                     mapkoper(j)=-oper
     128           0 :                     shiftkpt(:,j)=nint( (brot(:)-kpoints(:,j))/scale )
     129             :               endif
     130             :              enddo
     131             :            endif
     132             :           enddo
     133             :          endif 
     134             :       enddo   
     135             : 
     136             : c****************************************************
     137             : c         write results to files
     138             : c         w90kpts: whole Brillouin zone (for w90)
     139             : c         kpts: irreducible part (for fleur)
     140             : c         kptsmap: mapping from w90kpts to kpts
     141             : c****************************************************
     142           0 :       open(117,file='kptsmap',form='formatted',recl=1000)
     143           0 :       do ikpt=1,nkpts
     144           0 :          if(mapk(ikpt)==0)then
     145           0 :             write(117,*)ikpt,irreduc(ikpt),1,0,0,0
     146             :          else
     147           0 :             write(117,*)ikpt,irreduc(mapk(ikpt)),mapkoper(ikpt),
     148           0 :      &                  shiftkpt(:,ikpt)
     149             :          endif  
     150             :       enddo   
     151           0 :       close(117)
     152           0 :       open(119,file='kpts',form='formatted')
     153           0 :       if (film.and..not.l_onedimens)then
     154           0 :         write(119,'(i5,f20.10,3x,l1)')reduznumk,scale,.false.         
     155             :       else
     156           0 :         write(119,'(i5,f20.10)')reduznumk,scale
     157             :       endif   
     158           0 :       sumweights=0
     159           0 :       do ikpt=1,nkpts
     160           0 :          bkpt(:)=kpoints(:,ikpt)
     161           0 :          if (mapk(ikpt)==0)then
     162           0 :             sumweights=sumweights+weight(ikpt)
     163           0 :             write(6,*)"ikpt=",ikpt
     164           0 :             write(6,*)"irreducible"
     165           0 :             write(6,fmt='(a10,3f9.6)')"internal: ",bkpt(:)/scale
     166             : 
     167           0 :             if (film.and..not.l_onedimens)then
     168           0 :                write(119,'(3f10.5)')bkpt(1:2),weight(ikpt)
     169             :             else
     170           0 :                write(119,'(4f10.5)')bkpt(:),weight(ikpt)
     171             :             endif   
     172             :       
     173           0 :          elseif(mapkoper(ikpt).gt.0)then
     174           0 :             write(6,*)"ikpt=",ikpt
     175           0 :             write(6,*)"reducible"
     176           0 :             write(6,*)"map=",mapk(ikpt)
     177           0 :             brot(:)=0.0
     178           0 :             do k=1,3
     179             :                brot(:)=brot(:)+
     180           0 :      + mrot(k,:,mapkoper(ikpt))*kpoints(k,mapk(ikpt))
     181             :             enddo   
     182           0 :             write(6,'(a19,3f9.6)')"rotated internal: ",brot(:)/scale
     183           0 :          elseif(mapkoper(ikpt).lt.0)then
     184           0 :             write(6,*)"ikpt=",ikpt
     185           0 :             write(6,*)"reducible"
     186           0 :             write(6,*)"map=",mapk(ikpt)
     187           0 :             brot(:)=0.0
     188           0 :             do k=1,3
     189             :                brot(:)=brot(:)-
     190           0 :      + mrot(k,:,-mapkoper(ikpt))*kpoints(k,mapk(ikpt))
     191             :             enddo   
     192           0 :             write(6,'(a19,3f9.6)')"rotated internal: ",brot(:)/scale
     193             :          endif   
     194             :       enddo   
     195           0 :       close(119)
     196             : 
     197           0 :       write(6,*)"reduznumk=",reduznumk     
     198           0 :       write(6,*)"nkpts=",nkpts
     199           0 :       write(6,*)"sumweights=",sumweights
     200             :       
     201           0 :       IF(sumweights/=nkpts) CALL juDFT_error
     202             :      +     ("sum of weights differs from nkpts",calledby
     203           0 :      +     ="wann_kptsreduc")
     204           0 :       deallocate(kpoints,mapk,mapkoper,weight,irreduc)
     205             : 
     206             :       contains
     207             :       real elemental function axmintx(x)
     208             :       implicit none
     209             :       real,intent(in) :: x
     210           0 :       axmintx = abs(x-nint(x))
     211             :       end function     
     212             : 
     213             :       end subroutine wann_kptsreduc
     214             :       end module m_wann_kptsreduc

Generated by: LCOV version 1.13