add_Vnonlocal.F90 6.34 KB
Newer Older
Daniel Wortmann's avatar
Daniel Wortmann committed
1 2 3 4 5 6 7 8 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 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77
MODULE m_add_vnonlocal
  USE m_judft
! 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
      CONTAINS

      SUBROUTINE add_vnonlocal(nk,hybrid,dimension, kpts,jsp,results,xcpot,hamovlp)

      USE m_symm_hf       ,ONLY: symm_hf
      USE m_util          ,ONLY: intgrf,intgrf_init
      USE m_exchange_valence_hf
      USE m_exchange_core
      USE m_symmetrizeh
      USE m_wrapper
      USE m_hsefunctional ,ONLY: exchange_vccvHSE,exchange_ccccHSE
      USE m_types
      USE m_io_hybrid
      IMPLICIT NONE
      TYPE(t_results),INTENT(INOUT)   :: results
      TYPE(t_xcpot),INTENT(IN)        :: xcpot
      TYPE(t_dimension),INTENT(IN)    :: dimension
      TYPE(t_hybrid),INTENT(INOUT)    :: hybrid
      TYPE(t_kpts),INTENT(IN)         :: kpts
   
      INTEGER,INTENT(IN)      :: jsp 
      INTEGER,INTENT(IN)      :: nk
     
      TYPE(t_hamovlp),INTENT(INOUT)::hamovlp

  

!     - local scalars -
      INTEGER                 ::  ic,n,nn
      INTEGER                 ::  iband

      REAL                    ::  a_ex

      TYPE(t_mat)             :: olap,tmp,v_x,z
      COMPLEX                 ::  exch(dimension%neigd,dimension%neigd)
    
      
    
      !
      ! initialize weighting factor for HF exchange part
      !
78
      a_ex=xcpot%get_exchange_weight()
Daniel Wortmann's avatar
Daniel Wortmann committed
79 80
      
      
81 82
      v_x%l_real=hamovlp%l_real
      v_x%matsize1=dimension%nbasfcn
Daniel Wortmann's avatar
Daniel Wortmann committed
83 84 85 86 87 88 89 90 91 92 93

       CALL read_v_x(v_x,kpts%nkpt*(jsp-1) + nk)
      ! add non-local x-potential to the hamiltonian a (in packed storage)
       ic=1
       DO n=1,v_x%matsize1
          DO nn=1,n           
             IF (hamovlp%l_real) THEN
                hamovlp%a_r(ic) = hamovlp%a_r(ic) - a_ex*v_x%data_r(n,nn)
             ELSE
                hamovlp%a_c(ic) = hamovlp%a_c(ic) - a_ex*v_x%data_c(n,nn)
             ENDIF
94
             ic=ic+1
Daniel Wortmann's avatar
Daniel Wortmann committed
95 96 97 98 99 100
          ENDDO
       END DO
       ! calculate HF energy
       IF( hybrid%l_calhf ) THEN
          WRITE(6,'(A)') new_line('n')//new_line('n')//' ###     '// '        diagonal HF exchange elements (eV)              ###'
          
101
          WRITE(6,'(A)') new_line('n') // '         k-point      '// 'band       tail      pole     total(valence+core)'
Daniel Wortmann's avatar
Daniel Wortmann committed
102 103 104 105 106 107 108 109 110 111 112 113 114 115
          
       END IF

       ! read in lower triangle part of overlap matrix from direct acces file olap
       call olap%alloc(hamovlp%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)
       
       call z%alloc(olap%l_real,dimension%nbasfcn,dimension%neigd)
       
       call read_z(z,nk) !what about spin?
       
       
       ! calculate exchange contribution of current k-point nk to total energy (te_hfex)
116
       ! in the case of a spin-unpolarized calculation the factor 2 is added in eigen.F90 
Daniel Wortmann's avatar
Daniel Wortmann committed
117 118
       if (.not.v_x%l_real) v_x%data_c=conjg(v_x%data_c) 
       exch = 0
119
       z%matsize1=MIN(z%matsize1,v_x%matsize2)
Daniel Wortmann's avatar
Daniel Wortmann committed
120 121 122 123 124 125 126 127 128 129 130 131
       call v_x%multiply(z,tmp)
       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))
          endif
          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,3F10.5)')&
132 133
                  &  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)
Daniel Wortmann's avatar
Daniel Wortmann committed
134 135 136 137 138 139 140 141
          END IF
       END DO
     
    
    END SUBROUTINE add_vnonlocal


  END MODULE m_add_vnonlocal