LCOV - code coverage report
Current view: top level - dos - tetra_dos.f (source / functions) Hit Total Coverage
Test: combined.info Lines: 75 75 100.0 %
Date: 2019-09-08 04:53:50 Functions: 1 1 100.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_tetrados
       8             : c----------------------------------------------------------------------
       9             : c
      10             : c This subroutine evaluates the density of states (g) by the linear
      11             : c tetrahedron method on a energy grid (e) of 'ned' points.
      12             : c
      13             : c ev()          ... eigenvalues 
      14             : c qal()         ... partial charges
      15             : c ntetra)       ... number of tetrahedrons
      16             : c itetra(1-4,nt)... index of k-points forming tetrahedron nt
      17             : c voltet(nt)    ... volume of tetrahedron nt
      18             : c omega_bz      ... volume of irreducible part of BZ
      19             : c
      20             : c                                                      gb 2000
      21             : c----------------------------------------------------------------------
      22             :       CONTAINS
      23           3 :       SUBROUTINE tetra_dos(
      24             :      >                     lmax,ntype,neigd,ned,ntetra,nkpt,
      25           6 :      >                     itetra,efermi,voltet,energy,nevk,
      26           6 :      >                     ev,qal,
      27           3 :      <                     g)
      28             : 
      29             :       IMPLICIT NONE
      30             : c
      31             : C     ..Scalar Arguments ..
      32             :       INTEGER, INTENT (IN)  :: ntype,neigd,ned,lmax
      33             :       INTEGER, INTENT (IN)  :: ntetra,nkpt
      34             :       REAL,    INTENT (IN)  :: efermi
      35             : C     ..
      36             : C     .. Array Arguments ..
      37             :       INTEGER, INTENT (IN)    :: itetra(4,6*nkpt),nevk(nkpt)
      38             :       REAL,    INTENT (IN)    :: voltet(6*nkpt),energy(ned)
      39             :       REAL,    INTENT (IN)    :: qal(lmax*ntype+3,neigd,nkpt)
      40             :       REAL,    INTENT (INOUT) :: ev(neigd,nkpt)
      41             :       REAL,    INTENT (OUT)   :: g(ned,lmax*ntype+3)
      42             : C     ..
      43             : C     .. Local Variables ..
      44             :       INTEGER i,j,neig,nk,ntp,ne,ns,nc,ntet
      45             :       REAL    ener,efer,w
      46             : C     ..
      47             : C     .. Local Arrays ..
      48           6 :       REAL weight(4),eval(4),ecmax(neigd),term(ned)
      49           6 :       REAL wpar(4,ntype,neigd,nkpt),wparint(neigd,nkpt)
      50           6 :       REAL wght(neigd,nkpt)
      51             : 
      52          63 :       DO nk = 1,nkpt
      53           3 :          ev(nevk(nk)+1:neigd,nk) = 1.0e10 
      54             :       ENDDO
      55          63 :       wpar(:,:,:,:) = 0.0
      56          63 :       wparint=0.0
      57             : 
      58         111 :       DO neig = 1,neigd
      59          54 :          ecmax(neig) = -1.0e25
      60        1137 :          DO nk = 1,nkpt
      61        1080 :             wght(neig,nk) = 0.0e0
      62        1134 :             IF ( ev(neig,nk).GT.ecmax(neig) ) ecmax(neig) = ev(neig,nk)
      63             :          ENDDO
      64             :       ENDDO
      65             : c
      66             : c  check for energy degeneracies in tetrahedrons
      67             : c
      68         135 :       DO ntet = 1,ntetra
      69        4887 :         DO neig = 1,neigd
      70        9636 :           DO i = 1,3
      71       21384 :             DO j = i+1,4
      72       14256 :               IF (abs(ev(neig,itetra(i,ntet)) -
      73        7128 :      +                ev(neig,itetra(j,ntet))).LT.1.0e-7) THEN
      74             :                    ev(neig,itetra(i,ntet))=
      75           6 :      +             ev(neig,itetra(i,ntet))+i*(1.e-7)*ntet
      76             :                    ev(neig,itetra(j,ntet))=
      77           6 :      +             ev(neig,itetra(j,ntet))-i*(1.e-7)*ntet
      78             :               ENDIF
      79             :             ENDDO
      80             :           ENDDO
      81             :         ENDDO
      82             :       ENDDO
      83           3 :       WRITE (*,*) 'in tetra_dos'  ! do not remove  this statement !
      84             : c
      85             : c calculate weight factors
      86             : c
      87         135 :       DO ntet = 1,ntetra
      88         132 :         w = 6.0*voltet(ntet)
      89        2511 :         DO neig = 1,neigd
      90       21384 :           DO i = 1,4
      91       11880 :             eval(i) = ev(neig,itetra(i,ntet))
      92             :           ENDDO
      93             : 
      94        2508 :           IF (max(eval(1),eval(2),eval(3),eval(4)).LT.9999.9) THEN
      95       21384 :             DO i=1,4
      96        9504 :               weight(i) = 1.0
      97       47520 :               DO j=1,4
      98       47520 :                 IF (i.NE.j) weight(i) = weight(i)*(eval(j)-eval(i))
      99             :               ENDDO
     100             :               wght(neig,itetra(i,ntet)) =
     101       11880 :      +        wght(neig,itetra(i,ntet)) + w/weight(i)
     102             :             ENDDO
     103             :           ENDIF
     104             :         ENDDO
     105             :       ENDDO
     106             : c
     107             : c calculate partial weights
     108             : c
     109         123 :       DO nk=1,nkpt
     110        1083 :         DO neig = 1,nevk(nk)
     111        5460 :           DO ntp = 1,ntype
     112        2160 :             nc = lmax*(ntp-1)
     113             : 
     114       98280 :             DO ntet = 1,ntetra
     115       95040 :               IF (ALL(itetra(:,ntet).ne.nk)) CYCLE
     116             :            
     117       19008 :               eval(1:4) = ev(neig,itetra(1:4,ntet))
     118             : 
     119       21168 :               IF (max(eval(1),eval(2),eval(3),eval(4)).LT.9999.9) THEN
     120      171072 :                 DO i=1,4
     121       76032 :                   weight(i)=1.0
     122      380160 :                   DO j=1,4
     123      380160 :                     IF (i.NE.j) weight(i)=weight(i)*(eval(j)-eval(i))
     124             :                   ENDDO
     125       76032 :                   weight(i)=6.0*voltet(ntet)/weight(i)
     126      380160 :                   DO ns=1,4
     127             :                     wpar(ns,ntp,neig,itetra(i,ntet)) =
     128             :      +              wpar(ns,ntp,neig,itetra(i,ntet)) + 
     129      380160 :      +               0.25*weight(i)*qal(nc+ns,neig,nk)
     130             :                   ENDDO
     131       76032 :                   IF (ntp.EQ.1) wparint(neig,itetra(i,ntet)) = 
     132             :      +                          wparint(neig,itetra(i,ntet)) + 
     133       57024 :      +                0.25*weight(i)*qal(lmax*ntype+1,neig,nk)
     134             :                 ENDDO
     135             :               ENDIF
     136             : 
     137             :             ENDDO
     138             :           ENDDO
     139             :         ENDDO
     140             :       ENDDO
     141             : c
     142             : c---------------------------------------------------
     143             : c evaluate d.o.s.
     144             : c---------------------------------------------------
     145             : c
     146          36 :       g(:,:) = 0.0
     147             : 
     148         123 :       DO nk = 1,nkpt
     149        2223 :         DO neig = 1,neigd
     150             :           
     151        1080 :           ener = ev(neig,nk)
     152        3240 :           DO ntp = 1,ntype
     153       20520 :             DO ns = 1,lmax
     154        8640 :               nc = ns + lmax*(ntp-1)
     155        8640 :               w  = 0.5*wpar(ns,ntp,neig,nk)
     156    11251440 :               DO ne = 1,ned
     157    11240640 :                 term(ne) = energy(ne) - ener
     158    11240640 :                 IF ( energy(ne).GT.ecmax(neig) )
     159     1002400 :      +                term(ne) = ecmax(neig) - ener
     160    11240640 :                 IF ( term(ne).LT.0.0e0 ) term(ne) = 0.0e0
     161    11249280 :                 g(ne,nc) = g(ne,nc) + w * term(ne)**2
     162             :               ENDDO
     163             :             ENDDO
     164             :           ENDDO
     165             : 
     166        1080 :           nc = lmax*ntype+1
     167        1080 :           w = 0.5*wparint(neig,nk)
     168     1406220 :           DO ne = 1,ned
     169     1405080 :             term(ne) = energy(ne) - ener
     170     1405080 :             IF ( energy(ne).GT.ecmax(neig) ) term(ne) = ecmax(neig)-ener
     171     1405080 :             IF ( term(ne).lt.0.0e0 ) term(ne)=0.0e0
     172     1406160 :             g(ne,nc) = g(ne,nc) + w * term(ne)**2
     173             :           ENDDO
     174             : 
     175             :         ENDDO
     176             :       ENDDO
     177             : 
     178           3 :       END SUBROUTINE tetra_dos
     179             :       END MODULE m_tetrados

Generated by: LCOV version 1.13