add_Vnonlocal.F90 6.65 KB
Newer Older
1 2 3 4 5 6
!--------------------------------------------------------------------------------
! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
! This file is part of FLEUR and available as free software under the conditions
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------

Daniel Wortmann's avatar
Daniel Wortmann committed
7
MODULE m_add_vnonlocal
8
   USE m_judft
Daniel Wortmann's avatar
Daniel Wortmann committed
9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42
! c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c
!     This module is the driver routine for the calculation of the Hartree    c
!     Fock exchange term by using the mixed basis set.                        c
!                                                                             c
!     hsfock                                                                  c
!         |                                                                   c
!         |- symm.F:                                                          c
!         |  calculates the irreducible representation                        c
!         |                                                                   c
!         |- wavefproducts.F:                 s      s*                       c
!         |  computes the repsentation of phi    phi       in the mixed basis c
!         |                                  n,k    n',k+q                    c
!         |                                                                   c
!         |- exchange.F:                                                      c
!         |  calculates valence-valence part of the exchange matrix (mat_ex), c
!         |                                                                   c
!         |- exchange_core.F                                                  c
!         |  calculate valence-core contribution                              c
!                                                                             c
!     variables:                                                              c
!         kpts%nkptf   :=   number of kpoints                                      c
!         kpts%nkpt   :=   number of irreducible kpoints                          c
!         nbands  :=   number of bands for which the exchange matrix (mat_ex) c
!                      in the space of the wavefunctions is calculated        c
!         te_hfex :=   hf exchange contribution to the total energy           c
!         mnobd   :=   maximum number of occupied bands                       c
!         parent  :=   parent(ikpt) points to the symmetry equivalent point   c
!                      under the little group of kpoint nk                    c
!         symop   :=   symop(ikpt) points to the symmetry operation, which    c
!                      maps parent(ikpt) on ikpt                              c
!                                                                             c
!                                                                             c
!                                               M.Betzinger (09/07)           c
! c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c
43
   CONTAINS
Daniel Wortmann's avatar
Daniel Wortmann committed
44

45
   SUBROUTINE add_vnonlocal(nk,hybrid,dimension,kpts,jsp,results,xcpot,hmat)
Daniel Wortmann's avatar
Daniel Wortmann committed
46

47 48
      USE m_symm_hf,            ONLY: symm_hf
      USE m_util,               ONLY: intgrf,intgrf_init
Daniel Wortmann's avatar
Daniel Wortmann committed
49 50 51 52
      USE m_exchange_valence_hf
      USE m_exchange_core
      USE m_symmetrizeh
      USE m_wrapper
53
      USE m_hsefunctional,      ONLY: exchange_vccvHSE,exchange_ccccHSE
Daniel Wortmann's avatar
Daniel Wortmann committed
54 55 56
      USE m_types
      USE m_io_hybrid

57
      IMPLICIT NONE
Daniel Wortmann's avatar
Daniel Wortmann committed
58

59 60 61 62 63 64 65 66 67
      TYPE(t_results),       INTENT(INOUT) :: results
      CLASS(t_xcpot),        INTENT(IN)    :: xcpot
      TYPE(t_dimension),     INTENT(IN)    :: dimension
      TYPE(t_hybrid),        INTENT(INOUT) :: hybrid
      TYPE(t_kpts),          INTENT(IN)    :: kpts
      TYPE(t_mat),           INTENT(INOUT) :: hmat
   
      INTEGER,               INTENT(IN)    :: jsp 
      INTEGER,               INTENT(IN)    :: nk
Daniel Wortmann's avatar
Daniel Wortmann committed
68

69 70 71
      ! local scalars
      INTEGER                 :: n,nn,iband
      REAL                    :: a_ex
Daniel Wortmann's avatar
Daniel Wortmann committed
72
      TYPE(t_mat)             :: olap,tmp,v_x,z
73 74
      COMPLEX                 :: exch(dimension%neigd,dimension%neigd)

Daniel Wortmann's avatar
Daniel Wortmann committed
75
      ! initialize weighting factor for HF exchange part
76
      a_ex=xcpot%get_exchange_weight()      
Daniel Wortmann's avatar
Daniel Wortmann committed
77
      
78
      v_x%l_real=hmat%l_real
79
      v_x%matsize1=dimension%nbasfcn
Daniel Wortmann's avatar
Daniel Wortmann committed
80

81 82 83 84 85 86 87 88 89 90 91 92 93 94
      CALL read_v_x(v_x,kpts%nkpt*(jsp-1) + nk)
      ! add non-local x-potential to the hamiltonian hmat
      DO n=1,v_x%matsize1
         DO nn=1,n           
            IF (hmat%l_real) THEN
               hmat%data_r(n,nn) = hmat%data_r(n,nn) - a_ex*v_x%data_r(n,nn)
            ELSE
               hmat%data_c(n,nn) = hmat%data_c(n,nn) - a_ex*v_x%data_c(n,nn)
            ENDIF
         END DO
      END DO
      ! calculate HF energy
      IF(hybrid%l_calhf) THEN
         WRITE(6,'(A)') new_line('n')//new_line('n')//' ###     '// '        diagonal HF exchange elements (eV)              ###'
Daniel Wortmann's avatar
Daniel Wortmann committed
95
          
96 97
         WRITE(6,'(A)') new_line('n') // '         k-point      '// 'band          tail           pole       total(valence+core)'
      END IF
Daniel Wortmann's avatar
Daniel Wortmann committed
98

99 100 101 102
      ! read in lower triangle part of overlap matrix from direct acces file olap
      CALL olap%alloc(hmat%l_real,dimension%nbasfcn)
      CALL read_olap(olap,kpts%nkpt*(jsp-1) + nk)
      IF (.NOT.olap%l_real) olap%data_c=conjg(olap%data_c)
Daniel Wortmann's avatar
Daniel Wortmann committed
103
       
104
      CALL z%alloc(olap%l_real,dimension%nbasfcn,dimension%neigd)
Daniel Wortmann's avatar
Daniel Wortmann committed
105
       
106 107
      CALL read_z(z,nk) !what about spin?
      WRITE(*,*) 'What about spin (in add_Vnonlocal)?'
Daniel Wortmann's avatar
Daniel Wortmann committed
108
       
109 110 111 112 113 114 115
      ! calculate exchange contribution of current k-point nk to total energy (te_hfex)
      ! in the case of a spin-unpolarized calculation the factor 2 is added in eigen.F90 
      IF (.NOT.v_x%l_real) v_x%data_c=conjg(v_x%data_c) 
      exch = 0
      z%matsize1=MIN(z%matsize1,v_x%matsize2)

      CALL v_x%multiply(z,tmp)
Daniel Wortmann's avatar
Daniel Wortmann committed
116

117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
      DO iband = 1, hybrid%nbands(nk)
         IF (z%l_real) THEN
            exch(iband,iband) = dot_product(z%data_r(:z%matsize1,iband),tmp%data_r(:,iband))
         ELSE
            exch(iband,iband) = dot_product(z%data_r(:z%matsize1,iband),tmp%data_r(:,iband))
         END IF
         IF(iband.LE.hybrid%nobd(nk)) THEN
            results%te_hfex%valence = results%te_hfex%valence -a_ex*results%w_iks(iband,nk,jsp)*exch(iband,iband)
         END IF
         IF(hybrid%l_calhf) THEN
            WRITE(6, '(      ''  ('',F5.3,'','',F5.3,'','',F5.3,'')'',I4,4X,3F15.5)')&
                    kpts%bkf(:,nk),iband, (REAL(exch(iband,iband))-hybrid%div_vv(iband,nk,jsp))*(-27.211608),&
                    hybrid%div_vv(iband,nk,jsp)*(-27.211608),REAL(exch(iband,iband))*(-27.211608)
         END IF
      END DO
   END SUBROUTINE add_vnonlocal
Daniel Wortmann's avatar
Daniel Wortmann committed
133

134
END MODULE m_add_vnonlocal