LCOV - code coverage report
Current view: top level - ldau - umtx.f90 (source / functions) Hit Total Coverage
Test: combined.info Lines: 50 52 96.2 %
Date: 2019-09-08 04:53:50 Functions: 1 1 100.0 %

          Line data    Source code
       1             : MODULE m_umtx
       2             :   USE m_juDFT
       3             :   !*********************************************************************
       4             :   !* The calculation of the "U"-contribution to Hartree-Fock matrix.   *
       5             :   !*-------------------------------------------------------------------*
       6             :   !* Extension to multiple U per atom type by G.M. 2017                *
       7             :   !*********************************************************************
       8             : CONTAINS
       9          28 :   SUBROUTINE umtx(atoms,f0,f2,f4,f6,&
      10          28 :                   u)
      11             : 
      12             :     USE m_constants
      13             :     USE m_sgaunt
      14             :     USE m_types
      15             :     IMPLICIT NONE
      16             : 
      17             :     TYPE(t_atoms), INTENT(IN)  :: atoms
      18             :     REAL,          INTENT(IN)  :: f0(atoms%n_u),f2(atoms%n_u),f4(atoms%n_u),f6(atoms%n_u)
      19             :     REAL,          INTENT(OUT) :: u(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,&
      20             :                                     -lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,atoms%n_u)
      21             : 
      22             :     INTEGER, PARAMETER         :: lmaxw=3,lmmaxw1=(2*lmaxw+2)**2
      23             : 
      24             :     INTEGER i,j,k,l,m,mk,nfk,itype,i_u
      25             :     INTEGER m1,m2,m3,m4,lm1,lm2,lm3,lm4,kf
      26             :     REAL    uk,uq,avu,avj,cgk1,cgk2,tol
      27          56 :     REAL    fk(lmaxU_const+1,atoms%n_u)
      28          28 :     REAL,   ALLOCATABLE :: c(:,:,:)
      29             :     !
      30          28 :     tol = 1.0e-14
      31             :     !
      32             :     ! transformation to Hr-units:
      33             :     !
      34         140 :     DO i_u = 1, atoms%n_u
      35         112 :        itype = atoms%lda_u(i_u)%atomType
      36         112 :        l = atoms%lda_u(i_u)%l
      37         112 :        fk(1,i_u) = f0(i_u) / hartree_to_ev_const
      38         112 :        fk(2,i_u) = f2(i_u) / hartree_to_ev_const
      39         112 :        fk(3,i_u) = f4(i_u) / hartree_to_ev_const
      40         140 :        IF (l.EQ.3) THEN
      41           0 :           fk(4,i_u) = f6(i_u) / hartree_to_ev_const
      42         112 :        ELSE IF (l.GT.3) THEN
      43           0 :           CALL juDFT_error("LDA+U for p, d or f-states!", calledby="umtx")
      44             :        END IF
      45             :     END DO
      46             :     !
      47             :     ! evaluate Gaunt parameter
      48             :     !
      49          28 :     ALLOCATE(c(0:2*lmaxw+1,lmmaxw1,lmmaxw1))
      50        3612 :     DO k = 1,lmmaxw1
      51      231196 :        DO j = 1,lmmaxw1
      52     1951488 :           DO i = 0,2*lmaxw+1
      53     1032192 :              c(i,j,k) = 0.0
      54             :           END DO
      55             :        END DO
      56             :     END DO
      57             : 
      58          28 :     CALL sgaunt(lmaxw,lmmaxw1,lmaxU_const,c)
      59             : 
      60         140 :     DO i_u = 1, atoms%n_u                     !!! over U parameters
      61         112 :        l = atoms%lda_u(i_u)%l
      62         112 :        kf = 2*l
      63         560 :        DO m1 = -l,l
      64         448 :           lm1 = l*(l+1)+m1+1
      65        2464 :           DO m2 = -l,l
      66        1904 :              lm2 = l*(l+1)+m2+1
      67       10864 :              DO m3 = -l,l
      68        8512 :                 lm3 = l*(l+1)+m3+1
      69       49952 :                 DO m4 = -l,l
      70       39536 :                    lm4 = l*(l+1)+m4+1
      71       39536 :                    uk = 0.0e0
      72      153608 :                    DO k=0,kf,2
      73      114072 :                       uq = 0.e0
      74      666288 :                       DO mk=-k,k
      75      552216 :                          IF (mk.NE.m1-m3)  CYCLE
      76       74648 :                          cgk1 = c(k/2,lm1,lm3)
      77       74648 :                          IF (ABS(cgk1).LT.tol) CYCLE
      78       74648 :                          IF (mk.NE.m4-m2)  CYCLE
      79       11928 :                          cgk2 = c(k/2,lm4,lm2)
      80       11928 :                          IF (ABS(cgk2).LT.tol) CYCLE
      81      654360 :                          uq = uq+cgk1*cgk2
      82             :                       END DO                   ! mk
      83      114072 :                       IF (ABS(uq).LT.tol) CYCLE
      84       11928 :                       nfk=k/2+1
      85      153608 :                       uk=uk+uq*fk(nfk,i_u)*4*pi_const/(2*k+1)
      86             :                    END DO                     ! k
      87       48048 :                    u(m1,m2,m3,m4,i_u)=uk
      88             :                 END DO                       ! m4 etc.
      89             :              END DO
      90             :           END DO
      91             :        END DO
      92             :        avu=0.e0
      93             :        avj=0.e0
      94             : 
      95             :        DO i = -l,l
      96             :           DO j = -l,l
      97             :              avu = avu+u(i,j,i,j,i_u)
      98             :              avj = avj+(u(i,j,i,j,i_u)-u(i,j,j,i,i_u))
      99             :           END DO
     100             :        END DO
     101         112 :        avu = avu/(2*l+1)/(2*l+1)
     102         112 :        avj = avj/(2*l+1)/(2*l)
     103         140 :        avj = avu-avJ
     104             :        !        WRITE (6,*) 'U-matr:'
     105             :        !        IF (l.eq.2) WRITE (6,111) ((u(i,j,i,j,i_u),i=-l,l),j=-l,l)
     106             :        !        IF (l.eq.3) WRITE (6,211) ((u(i,j,i,j,i_u),i=-l,l),j=-l,l)
     107             :        !        WRITE (6,*) 'J-matr:'
     108             :        !        IF (l.eq.2) WRITE (6,111) ((u(i,j,j,i,i_u),i=-l,l),j=-l,l)
     109             :        !        IF (l.eq.3) WRITE (6,211) ((u(i,j,j,i,i_u),i=-l,l),j=-l,l)
     110             :        !         PRINT*,'U-av:',avu*hartree_to_ev_const
     111             :        !         PRINT*,'J-av:',avj*hartree_to_ev_const
     112             : 111    FORMAT (5f8.4)
     113             : 211    FORMAT (7f8.4)
     114             : 112    FORMAT (10e20.10)
     115             :        !c         WRITE (9,112) ((((u(i,j,k,m,i_u),i=-l,l),j=-l,l),k=-l,l),m=-l,l)
     116             :     END DO
     117          28 :     DEALLOCATE (c)
     118             : 
     119          28 :   END SUBROUTINE umtx
     120             : END MODULE m_umtx

Generated by: LCOV version 1.13