LCOV - code coverage report
Current view: top level - propcalc/dos - types_eigdos.F90 (source / functions) Hit Total Coverage
Test: FLEUR test coverage Lines: 106 115 92.2 %
Date: 2024-05-15 04:28:08 Functions: 11 14 78.6 %

          Line data    Source code
       1             : !--------------------------------------------------------------------------------
       2             : ! Copyright (c) 2018 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             : MODULE m_types_eigdos
       7             :   USE m_juDFT
       8             :   use m_constants
       9             :   IMPLICIT NONE
      10             :   PRIVATE
      11             :   PUBLIC:: t_eigdos,t_eigdos_list,t_eigdos_make_dos
      12             : 
      13             :   TYPE,abstract :: t_eigdos
      14             :     CHARACTER(len=20)            :: name_of_dos="unnamed"
      15             :     !each eigenvalue might be described by weights
      16             :     REAL,ALLOCATABLE             :: eig(:,:,:)
      17             :     REAL,ALLOCATABLE             :: dos_grid(:) !This is the grid the DOS part uses internally (FOR IO use the routine below)
      18             :     REAL,ALLOCATABLE             :: dos(:,:,:) !(grid,spins,weights)
      19             :   CONTAINS
      20             :     procedure(get_weight_name),DEFERRED    :: get_weight_name
      21             :     procedure(get_num_weights),DEFERRED    :: get_num_weights
      22             :     procedure(get_weight_eig),DEFERRED     :: get_weight_eig !should be overwritten in derived type
      23             :     procedure          :: get_spins
      24             :     procedure          :: get_neig
      25             :     procedure          :: get_eig
      26             :     procedure          :: get_dos_grid
      27             :     procedure          :: make_dos=>t_eigdos_make_dos
      28             :     procedure          :: smooth=>dosdata_smooth
      29             :     procedure          :: write_raw   !should be implemented later to allow eig66 functionality
      30             :     procedure          :: write_dos
      31             :     procedure          :: write_band
      32             :     procedure          :: write_EVData
      33             :     procedure          :: sym_weights,sym_weights_eigdos
      34             :   END TYPE
      35             : 
      36             :   type::t_eigdos_list
      37             :     class(t_eigdos),POINTER:: p
      38             :   end type
      39             : 
      40             :   INTERFACE
      41             :     function get_weight_eig(this,id)
      42             :       import t_eigdos
      43             :       class(t_eigdos),intent(in):: this
      44             :       INTEGER,intent(in)         :: id
      45             :       real,allocatable:: get_weight_eig(:,:,:)
      46             :     end function
      47             :   END interface
      48             :   INTERFACE
      49             :     integer function get_num_weights(this)
      50             :       import t_eigdos
      51             :       class(t_eigdos),intent(in):: this
      52             :     end function
      53             :   end interface
      54             :   INTERFACE
      55             :     character(len=20) function get_weight_name(this,id)
      56             :       import t_eigdos
      57             :       class(t_eigdos),intent(in):: this
      58             :       INTEGER,intent(in)         :: id
      59             :     end function
      60             :   end interface
      61             : CONTAINS
      62           1 :   subroutine sym_weights(this)
      63             :     class(t_eigdos),INTENT(INOUT):: this 
      64           1 :   end subroutine
      65             : 
      66         142 :   subroutine sym_weights_eigdos(this,weight)
      67             :     !Subroutine to make all weights equal for degenerate states
      68             :     class(t_eigdos),INTENT(IN):: this 
      69             :     real,intent(inout)        :: weight(:,:,:)
      70             : 
      71             :     INTEGER:: ispin,ikpt,i,j
      72             : 
      73         304 :     DO ispin=1,size(this%eig,3)
      74        3092 :       DO ikpt=1,size(this%eig,2)   
      75             :         ! Make sure equivalent states have same weight
      76             :         i=1   
      77       53868 :         DO while(i<size(this%eig,1))
      78             :          j=1
      79       68300 :          do while (abs(this%eig(i,ikpt,ispin)-this%eig(i+j,ikpt,ispin))<1E-5)
      80       18850 :             j=j+1
      81       68300 :             if (i+j>size(this%eig,1)) exit
      82             :          ENDDO
      83       50918 :          if (j>1) THEN
      84       16864 :             j=j-1
      85             :             !Make sure all weights are equal
      86       88292 :             weight(i:i+j,ikpt,ispin)=sum(weight(i:i+j,ikpt,ispin))/(j+1)
      87             :             i=i+j
      88             :          endif
      89       50918 :          i=i+1   
      90             :         enddo 
      91             :       ENDDO
      92             :     ENDDO
      93         142 :   end subroutine
      94             : 
      95          71 :   function get_neig(this)
      96             :     CLASS(t_eigdos),INTENT(IN)::this
      97             :     integer,allocatable :: get_neig(:,:)
      98          71 :     real,allocatable::ev(:,:,:)
      99             :     integer ::k,j
     100             : 
     101       30252 :     ev=this%get_eig()
     102         284 :     allocate(get_neig(size(ev,2),size(ev,3)))
     103         144 :     DO j=1,this%get_spins()
     104        1572 :       DO k=1,size(ev,2)
     105       30181 :         get_neig(k,j)=count(ev(:,k,j)<1E99)
     106             :       ENDDO
     107             :     ENDDO
     108             :   end function
     109             : 
     110         158 :   pure integer function get_spins(this)
     111             :     CLASS(t_eigdos),INTENT(IN)::this
     112         158 :     get_spins=size(this%eig,3)
     113         158 :   END function
     114         528 :   function get_eig(this)
     115             :     CLASS(t_eigdos),INTENT(IN):: this
     116             :     real,allocatable:: get_eig(:,:,:)
     117      123219 :     get_eig=this%eig
     118             :   END function
     119             : 
     120         184 :   function get_dos_grid(this)
     121             :     CLASS(t_eigdos),INTENT(IN):: this
     122             :     real,allocatable:: get_dos_grid(:)
     123      121900 :     get_dos_grid=this%dos_grid
     124             :   END function
     125             : 
     126           5 : subroutine dosdata_smooth(eigdos,banddos)
     127             :   use m_smooth
     128             :   use m_types_banddos
     129             :   class(t_eigdos),INTENT(INOUT)  :: eigdos
     130             :   type(t_banddos),INTENT(IN)     :: banddos
     131             : 
     132             :   integer :: jspin,i
     133           5 :   real,allocatable :: dos_grid(:)
     134             : 
     135           5 :   if (.not.allocated(eigdos%dos)) return
     136          20 :   if (size(eigdos%dos)==0) return
     137          20 :   if (size(eigdos%dos)<1) return
     138           5 :   if (banddos%sig_dos==0.0) RETURN
     139        6610 :   dos_grid=eigdos%get_dos_grid()
     140             : 
     141       13220 :   IF ( NINT((maxval(dos_grid) - minval(dos_grid))/banddos%sig_dos) > size(dos_grid) ) THEN
     142           0 :     WRITE(oUnit,*) 'sig_dos too small for DOS smoothing:'
     143           0 :     WRITE(oUnit,*) 'Reduce energy window or enlarge banddos%sig_dos!'
     144           0 :     WRITE(oUnit,*) 'For now: no smoothing done'
     145           0 :     return
     146             :   ENDIF
     147             : 
     148          87 :   DO i=1,size(eigdos%dos,3)
     149         176 :     DO jspin=1,eigdos%get_spins()
     150         171 :       CALL smooth(dos_grid,eigdos%dos(:,jspin,i),banddos%sig_dos,size(eigdos%dos_grid))
     151             :     ENDDO
     152             :   ENDDO
     153           5 : END subroutine
     154             : 
     155           5 : subroutine write_dos(eigdos,hdf_id)
     156             : #ifdef CPP_HDF
     157             :     use HDF5
     158             :     use m_banddos_io
     159             : #endif
     160             :     class(t_eigdos),INTENT(INOUT):: eigdos
     161             : #ifdef CPP_HDF
     162             :     integer(HID_T),intent(in) ::hdf_id
     163             : #else
     164             :     integer,       intent(in) ::hdf_id
     165             : #endif
     166             :     integer:: jspin,i,ind,id, n
     167             :     character(len=100)::filename
     168           5 :     real,allocatable:: dos_grid(:)
     169             :     LOGICAL l_printTextDOS
     170             : 
     171           5 :     l_printTextDOS = .TRUE.
     172             : 
     173             : #ifdef CPP_HDF
     174          87 :     DO n=1,eigdos%get_num_weights()
     175          82 :       print *, "writedos:",n,eigdos%get_num_weights()
     176          87 :       call writedosData(hdf_ID,eigdos%name_of_dos,eigdos%get_dos_grid(),eigdos%get_weight_name(n),eigdos%dos(:,:,n))
     177             :     enddo
     178           5 :     IF(eigdos%get_num_weights().GT.40) THEN
     179           1 :        WRITE(*,*) 'Number of weights in ', TRIM(ADJUSTL(eigdos%name_of_dos)),' DOS too large for simple text output.'
     180           1 :        WRITE(*,*) 'Output only in banddos.hdf file.'
     181             :        l_printTextDOS = .FALSE.
     182             :     END IF
     183             : #endif
     184             : 
     185           1 :     IF (.NOT.l_printTextDOS) RETURN
     186           4 :     if (.not.allocated(eigdos%dos)) return
     187          16 :     if (size(eigdos%dos)==0) return
     188           9 :     DO jspin=1,eigdos%get_spins()
     189           5 :       write(filename,"(a,a,i0)") trim(eigdos%name_of_dos),".",jspin
     190           5 :       open(999,file=filename)
     191          48 :       write(999,"(999a21)") "#energy",(eigdos%get_weight_name(id),id=1,eigdos%get_num_weights())
     192          48 :       write(*,"(999a21)") filename,(eigdos%get_weight_name(id),id=1,eigdos%get_num_weights())
     193        6615 :       dos_grid=eigdos%get_dos_grid()
     194        6610 :       DO i=1,size(dos_grid)
     195       63413 :         write(999,"(999(e20.8,1x))") dos_grid(i)*hartree_to_ev_const,(eigdos%dos(i,jspin,id)/hartree_to_ev_const,id=1,eigdos%get_num_weights())
     196             :       ENDDO
     197           5 :       close(999)
     198           9 :       write(*,*) "done:",filename
     199             :     ENDDO
     200           4 :   END subroutine
     201             : 
     202           2 :   subroutine write_band(eigdos,kpts,title,cell,hdf_id,efermi,banddos)
     203             :     use m_types_kpts
     204             :     use m_types_cell
     205             :     use m_gnuplot_BS
     206             :     use m_types_banddos
     207             : #ifdef CPP_HDF
     208             :      use HDF5
     209             :      use m_banddos_io
     210             : #endif
     211             :     class(t_eigdos),INTENT(INOUT):: eigdos
     212             :     type(t_kpts),intent(in)      :: kpts
     213             :     type(t_banddos),INTENT(IN)   :: banddos
     214             :     type(t_cell),intent(in)      :: cell
     215             :     CHARACTER(LEN=*), INTENT(IN) :: title
     216             :     real,intent(in)              :: efermi
     217             : #ifdef CPP_HDF
     218             :     integer(HID_T),intent(in) ::hdf_id
     219             :     INTEGER::n
     220             : 
     221             : #else
     222             :     integer,intent(in):: hdf_id !not used
     223             : #endif
     224             : 
     225             :     integer:: jspin,i,k
     226           2 :     real,allocatable  :: ev(:,:,:),kx(:)
     227             :     real              :: vkr(3),vkr_prev(3)
     228             :     character(len=100)::file
     229           6 :     allocate(kx(kpts%nkpt))
     230             : #ifdef CPP_HDF
     231          28 :     DO n=1,eigdos%get_num_weights()
     232          28 :       call writebandData(hdf_id,kpts,eigdos%name_of_dos,eigdos%get_weight_name(n),eigdos%get_weight_eig(n),eigdos%get_eig())
     233             :     enddo
     234             : #endif
     235           2 :     if (eigdos%name_of_dos.ne."Local") then
     236             : #ifndef CPP_HDF
     237             :       print *,"WARNING, only very basic BS written without HDF"
     238             : #endif
     239           0 :       return
     240             :     endif
     241           5 :     DO jspin=1,eigdos%get_spins()
     242           3 :       write(file,"(a,i0)") "bands.",jspin
     243           3 :       open(18,file=file)
     244        4279 :       ev=eigdos%get_eig()
     245           3 :       kx(1) = 0.0
     246          39 :       vkr_prev=matmul(kpts%bk(:,1),cell%bmat)
     247          44 :       DO k = 2, kpts%nkpt
     248         533 :         vkr=matmul(kpts%bk(:,k),cell%bmat)
     249         164 :         kx(k)=kx(k-1)+ sqrt(dot_product(vkr-vkr_prev,vkr-vkr_prev))
     250          44 :         vkr_prev=vkr
     251             :       ENDDO
     252         254 :       DO i = 1, minval(eigdos%get_neig())
     253        2461 :         DO k = 1, kpts%nkpt
     254        2458 :           write(18,'(999f15.9)') kx(k),(ev(i,k,jspin)-efermi)*hartree_to_ev_const
     255             :           !-eFermiCorrection
     256             :         ENDDO
     257             :       ENDDO
     258           5 :       CLOSE (18)
     259             :     enddo
     260           2 :     call gnuplot_bs(kpts,title,cell,eigdos%get_spins())
     261           2 :     IF (banddos%unfoldband) call write_gnu_sc(banddos,kpts,title,cell,eigdos%get_spins())
     262           2 :   end subroutine
     263             : 
     264           5 :   subroutine write_EVData(eigdos,hdf_id)
     265             : #ifdef CPP_HDF
     266             :      use HDF5
     267             :      use m_banddos_io
     268             : #endif
     269             :      class(t_eigdos),INTENT(INOUT):: eigdos
     270             : #ifdef CPP_HDF
     271             :      integer(HID_T),intent(in) ::hdf_id
     272             :      INTEGER::n
     273             : #else
     274             :      integer,intent(in):: hdf_id !not used
     275             : #endif
     276             : 
     277             : #ifdef CPP_HDF
     278          87 :      DO n = 1, eigdos%get_num_weights()
     279          87 :         CALL writeEVData(hdf_id,eigdos%name_of_dos,eigdos%get_weight_name(n),eigdos%get_eig(),eigdos%get_weight_eig(n))
     280             :      END DO
     281             : #endif
     282           5 :   end subroutine
     283             : 
     284           5 :   subroutine t_eigdos_make_dos(eigdos,kpts,input,banddos,efermi)
     285             :     use m_types_banddos
     286             :     use m_types_input
     287             :     use m_dosbin
     288             :     use m_ptdos
     289             :     use m_tetra_dos
     290             :     use m_dostetra
     291             :     use m_types_kpts
     292             : 
     293             :     class(t_eigdos),intent(inout):: eigdos
     294             :     type(t_banddos),intent(in)   :: banddos
     295             :     type(t_kpts),intent(in)      :: kpts
     296             :     type(t_input),intent(in)     :: input
     297             :     real,intent(in)              :: efermi
     298             : 
     299             :     integer ::n
     300             :     real :: emin,emax
     301             : 
     302           5 :     emin=min(banddos%e1_dos,banddos%e2_dos)-efermi
     303           5 :     emax=max(banddos%e1_dos,banddos%e2_dos)-efermi
     304           5 :     if (allocated(eigdos%dos)) return
     305             : 
     306          25 :     allocate(eigdos%dos(banddos%ndos_points,eigdos%get_spins(),eigdos%get_num_weights()))
     307             :     !Generate DOS grid
     308           5 :     if (allocated(eigdos%dos_grid)) deallocate(eigdos%dos_grid)
     309          15 :     allocate(eigdos%dos_grid(banddos%ndos_points))
     310        6610 :     DO n=1,banddos%ndos_points
     311        6610 :       eigdos%dos_grid(n)=emin+(emax-emin)/(banddos%ndos_points-1.0)*(n-1.0)
     312             :     ENDDO
     313             : 
     314          87 :     DO n=1,eigdos%get_num_weights()
     315          82 :       print *,eigdos%name_of_dos,n,eigdos%get_num_weights()
     316           5 :       SELECT CASE(input%bz_integration)
     317             : 
     318             :       CASE(BZINT_METHOD_HIST, BZINT_METHOD_GAUSS)
     319          14 :         call dos_bin(input%jspins,kpts%wtkpt,eigdos%dos_grid,eigdos%get_eig(),eigdos%get_weight_eig(n),eigdos%dos(:,:,n),efermi)
     320             :       CASE(BZINT_METHOD_TRIA)
     321          68 :         IF(input%film) THEN
     322           0 :           CALL ptdos(input%jspins,kpts,eigdos%dos_grid,eigdos%get_eig(),eigdos%get_weight_eig(n),eigdos%dos(:,:,n),efermi)
     323             :         ELSE
     324             :           CALL tetra_dos(input%jspins,kpts,eigdos%dos_grid,eigdos%get_neig(),eigdos%get_eig(),&
     325          68 :                          eigdos%get_weight_eig(n),eigdos%dos(:,:,n),efermi)
     326             :         ENDIF
     327             :       CASE(BZINT_METHOD_TETRA)
     328          82 :         CALL dostetra(kpts,input,eigdos%dos_grid,eigdos%get_eig(),eigdos%get_weight_eig(n),eigdos%dos(:,:,n),efermi)
     329             :       END SELECT
     330             : 
     331             :     end do
     332             :   END subroutine
     333             : 
     334             : 
     335             : 
     336           0 :   subroutine write_raw(this,id)
     337             :     class(t_eigdos),INTENT(IN):: this
     338             :     INTEGER,INTENT(IN)        :: id
     339             : 
     340             : 
     341           0 :   end subroutine
     342             : 
     343             : 
     344             : 
     345           0 : END MODULE m_types_eigdos

Generated by: LCOV version 1.14