Commit dd1053a5 authored by Daniel Wortmann's avatar Daniel Wortmann

Reordered the hybrid code.

Now fleur.F90 calls a subroutine that does most of the hybrid stuff.
In eigen.F90 only the non-local potential is added to the Hamiltonian
parent 0f56b9dc
......@@ -13,7 +13,7 @@ MODULE m_apws
!*********************************************************************
CONTAINS
SUBROUTINE apws(dimension,input,noco,kpts,&
nk,cell,l_zref,n_size,jspin, bkpt,lapw,matind,nred)
nk,cell,l_zref,n_size,jspin,bkpt,lapw,matind,nred)
USE m_types
USE m_sort
......
......@@ -30,17 +30,14 @@ CONTAINS
USE m_usetup
USE m_pot_io
USE m_eigen_diag
USE m_eigen_hf_init
USE m_eigen_hf_setup
USE m_hsfock
USE m_add_vnonlocal
USE m_subvxc
USE m_hsefunctional
USE m_hybridmix , ONLY: amix_pbe0,amix_hf
!USE m_hsefunctional
USE m_util
USE m_icorrkeys
USE m_io_hybrid
!USE m_icorrkeys
USE m_eig66_io, ONLY : open_eig, write_eig, close_eig,read_eig
USE m_xmlOutput
USE m_io_hybrid
#ifdef CPP_MPI
USE m_mpi_bc_pot
#endif
......@@ -107,10 +104,8 @@ CONTAINS
! .. variables for HF or hybrid functional calculation ..
!
TYPE(t_hybdat) :: hybdat
INTEGER :: comm(kpts%nkpt),irank2(kpts%nkpt),isize2(kpts%nkpt)
REAL, ALLOCATABLE :: eig_irr(:,:),vr0(:,:,:)
#ifdef CPP_MPI
INTEGER :: sndreqd,sndreq(mpi%isize*kpts%nkpt)
#endif
......@@ -155,10 +150,6 @@ CONTAINS
8120 FORMAT (' IT=',i4,' ISEC1=',i4,' reduced diagonalization')
l_wu = .TRUE.
END IF
!
! load potential from file pottot (=unit 8)
!
ALLOCATE ( vr0(atoms%jmtd,atoms%ntype,DIMENSION%jspd) ) ; vr0 = 0
!IF (mpi%irank.EQ.0) THEN
! CALL readPotential(stars,vacuum,atoms,sphhar,input,sym,POT_ARCHIVE_TYPE_TOT_const,&
......@@ -180,8 +171,7 @@ CONTAINS
CALL lodpot(mpi,atoms,sphhar,obsolete,vacuum,&
input, v%mt,v%vacz, enpara_in, enpara)
!
CALL eigen_hf_init(hybrid,kpts,atoms,input,dimension,hybdat,irank2,isize2,l_real)
!---> set up and solve the eigenvalue problem
!---> loop over energy windows
......@@ -268,8 +258,7 @@ CONTAINS
!---> loop over k-points: each can be a separate task
DO jsp = 1,nspins
CALL eigen_HF_setup(hybrid,input,sym,kpts,dimension,atoms,mpi,noco,cell,oneD,results,jsp,eig_id,&
hybdat,irank2,it,l_real,vr0,eig_irr)
!
!---> set up k-point independent t(l'm',lm) matrices
......@@ -315,17 +304,8 @@ CONTAINS
#endif
k_loop:DO nk = mpi%n_start,kpts%nkpt,mpi%n_stride
#if defined(CPP_MPI)&&defined(CPP_NEVER)
IF ( hybrid%l_calhf ) THEN
! jump to next k-point if this process is not present in communicator
IF ( comm(nk) == MPI_COMM_NULL ) CYCLE
! allocate buffer for communication of the results
IF ( irank2(nk) /= 0 ) CALL work_dist_reserve_buffer( hybdat%nbands(nk) )
END IF
#endif
nrec = kpts%nkpt*(jsp-1) + nk
nrec = mpi%n_size*(nrec-1) + mpi%n_rank + 1
!
!---> set up lapw list
!
......@@ -362,15 +342,13 @@ CONTAINS
call olap%from_packed(l_real,dimension%nbasfcn,hamovlp%b_r,hamovlp%b_c)
call write_olap(olap,nrec)
CALL hsfock(nk,atoms,hybrid,lapw,DIMENSION,kpts,jsp,input,hybdat,eig_irr,&
sym,cell,noco,results,it,maxval(hybdat%nobd),xcpot,&
mpi,irank2(nk),isize2(nk),comm(nk), hamovlp)
if (hybrid%l_addhf) CALL add_Vnonlocal(nk,hybrid,dimension, kpts,jsp,results,xcpot,hamovlp)
IF ( irank2(nk) /= 0 ) CYCLE
IF( hybrid%l_subvxc ) THEN
CALL subvxc(lapw,kpts%bk(:,nk),DIMENSION,input,jsp,vr0,atoms,ud,hybrid,enpara%el0,enpara%ello0,&
CALL subvxc(lapw,kpts%bk(:,nk),DIMENSION,input,jsp,v%mt(:,0,:,:),atoms,ud,hybrid,enpara%el0,enpara%ello0,&
sym, atoms%nlotot,kveclo, cell,sphhar, stars, xcpot,mpi,&
oneD, hamovlp,vx)
END IF
......@@ -433,20 +411,6 @@ CONTAINS
CALL MPI_BARRIER(mpi%MPI_COMM,ierr)
#endif
#if defined(CPP_MPI)&&defined(CPP_NEVER)
IF ( hybrid%l_calhf ) THEN
IF ( isize2(nk) == 1 ) THEN
WRITE(*,'(a,i6,a,i6,a)') 'HF: kpt ', nk, ' was done by rank ', mpi%irank, '.'
ELSE
WRITE(*,'(a,i6,a,i6,a,i6,a)')&
'HF: kpt ', nk, ' was done by rank ', mpi%irank, ' and ', isize2(nk)-1, ' more.'
END IF
! ELSE
! WRITE(*,'(a,i6,a,i6,a)') ' kpt ', nk, ' was done by rank ', irank, '.'
END IF
!# else
! WRITE (*,*) 'pe: ',irank,' wrote ',nrec
# endif
CALL timestop("EV output")
!#ifdef CPP_MPI
IF (l_real) THEN
......@@ -460,11 +424,6 @@ ENDIF
DEALLOCATE (td%tuu,td%tud,td%tdu,td%tdd)
DEALLOCATE (td%ind,td%tuulo,td%tdulo)
DEALLOCATE (td%tuloulo)
#ifdef CPP_NEVER
IF ( hybrid%l_calhf ) THEN
DEALLOCATE ( eig_irr,kveclo_eig )
END IF
#endif
END DO ! spin loop ends
DEALLOCATE( vs_mmp )
DEALLOCATE (matind)
......@@ -473,30 +432,11 @@ ENDIF
ELSE
DEALLOCATE(hamOvlp%a_c,hamOvlp%b_c)
ENDIF
#ifdef CPP_NEVER
IF( hybrid%l_calhf ) THEN
DEALLOCATE( fac,sfac,gauntarr )
DEALLOCATE( nindxc,core1,core2,nbasm,eig_c )
END IF
#endif
#if defined(CPP_MPI)&&defined(CPP_NEVER)
IF ( hybrid%l_calhf ) DEALLOCATE (nkpt_EIBZ)
#endif
! hf: write out radial potential vr0
IF (hybrid%l_hybrid.OR.hybrid%l_calhf) THEN
OPEN(unit=120,file='vr0',form='unformatted')
DO isp=1,DIMENSION%jspd
DO nn=1,atoms%ntype
DO i=1,atoms%jmtd
WRITE(120) vr0(i,nn,isp)
END DO
END DO
END DO
CLOSE(120)
ENDIF
#ifdef CPP_MPI
CALL MPI_BARRIER(mpi%MPI_COMM,ierr)
#endif
......
......@@ -209,7 +209,7 @@ CONTAINS
ALLOCATE (hybrid_temp%nindx(0:atoms%lmaxd,atoms%ntype))
ALLOCATE (hybrid_temp%select1(4,atoms%ntype),hybrid_temp%lcutm1(atoms%ntype))
ALLOCATE (hybrid_temp%select2(4,atoms%ntype),hybrid_temp%lcutm2(atoms%ntype),hybrid_temp%lcutwf(atoms%ntype))
ALLOCATE (hybrid_temp%lcutwf(atoms%ntype))
CALL rw_inp('r',atoms_temp,obsolete_temp,vacuum_temp,input_temp,stars_temp,sliceplot_temp,&
banddos_temp,cell_temp,sym_temp,xcpot_temp,noco_temp,Jij_temp,oneD_temp,hybrid_temp,&
......
MODULE m_icorrkeys
INTEGER, PARAMETER :: icorr_exx = -3
INTEGER, PARAMETER :: icorr_pbe0 = 12
INTEGER, PARAMETER :: icorr_hf = -2
!HSE, hybrid exchange functional: J. Chem. Phys. 118, 8207 (2003)
INTEGER, PARAMETER :: icorr_hse = 13
! only local part of HSE
INTEGER, PARAMETER :: icorr_hseloc = 14
! hybrid functional similar to HSE but with variable screening and mixing parameter
INTEGER, PARAMETER :: icorr_vhse = 15
INTEGER, PARAMETER :: icorr_exx = -3
INTEGER, PARAMETER :: icorr_pbe0 = 12
INTEGER, PARAMETER :: icorr_hf = -2
!HSE, hybrid exchange functional: J. Chem. Phys. 118, 8207 (2003)
INTEGER, PARAMETER :: icorr_hse = 13
! only local part of HSE
INTEGER, PARAMETER :: icorr_hseloc = 14
! hybrid functional similar to HSE but with variable screening and mixing parameter
INTEGER, PARAMETER :: icorr_vhse = 15
CONTAINS
FUNCTION get_exchange_weight(icorr) result(a_ex)
USE m_hybridmix
USE m_judft
IMPLICIT NONE
INTEGER,INTENT(IN)::icorr
REAL a_ex
SELECT CASE (icorr)
CASE (icorr_pbe0)
a_ex = amix_pbe0
CASE (icorr_hf )
a_ex = amix_hf
CASE ( icorr_hse)
a_ex = aMix_HSE
CASE (icorr_vhse )
a_ex = aMix_VHSE()
CASE DEFAULT
call judft_error('xc functional can not be identified')
END SELECT
END
END MODULE m_icorrkeys
......@@ -397,56 +397,41 @@ MODULE m_types
TYPE t_hybrid
LOGICAL :: l_hybrid
LOGICAL :: l_subvxc
LOGICAL :: l_calhf
LOGICAL :: l_addhf
INTEGER :: ewaldlambda
INTEGER :: lexp
INTEGER :: bands1 !Only read in
INTEGER :: bands2 !Only read in
INTEGER :: nbasp
INTEGER :: nbasp
INTEGER :: maxlcutm1
INTEGER :: maxindxm1
INTEGER :: maxbasm1
INTEGER :: maxlcutm2
INTEGER :: maxindxm2
INTEGER :: maxbasm2
INTEGER :: maxindxp1
INTEGER :: maxindxp2
INTEGER :: maxgptm
INTEGER :: maxgptm1
INTEGER :: maxgptm2
INTEGER :: maxindx
INTEGER :: maxlmindx
INTEGER :: gptmd
INTEGER,ALLOCATABLE :: nindx(:,:)
INTEGER,ALLOCATABLE :: select1(:,:)
INTEGER,ALLOCATABLE :: lcutm1(:)
INTEGER,ALLOCATABLE :: select2(:,:)
INTEGER,ALLOCATABLE :: lcutm2(:)
INTEGER,ALLOCATABLE :: nindxm1(:,:)
INTEGER,ALLOCATABLE :: nindxm2(:,:)
INTEGER,ALLOCATABLE :: gptm(:,:)
INTEGER,ALLOCATABLE :: ngptm1(:)
INTEGER,ALLOCATABLE :: pgptm1(:,:)
INTEGER,ALLOCATABLE :: ngptm2(:)
INTEGER,ALLOCATABLE :: pgptm2(:,:)
INTEGER,ALLOCATABLE :: ngptm (:)
INTEGER,ALLOCATABLE :: pgptm (:,:)
INTEGER,ALLOCATABLE :: nindxp1(:,:)
INTEGER,ALLOCATABLE :: nindxp2(:,:)
INTEGER,ALLOCATABLE :: lcutwf(:)
INTEGER,ALLOCATABLE :: map(:,:)
INTEGER,ALLOCATABLE :: tvec(:,:,:)
REAL :: radshmin
INTEGER , ALLOCATABLE :: nbasm(:)
REAL :: gcutm1
REAL :: gcutm2
REAL :: tolerance1 !only read in
REAL :: tolerance2 !only read in
REAL ,ALLOCATABLE :: ddist(:)
REAL ,ALLOCATABLE :: basm1(:,:,:,:)
REAL ,ALLOCATABLE :: basm2(:,:,:,:)
COMPLEX,ALLOCATABLE :: d_wgn2(:,:,:,:)
LOGICAL :: l_subvxc
LOGICAL :: l_calhf
LOGICAL,ALLOCATABLE :: l_exxc(:,:)
INTEGER,ALLOCATABLE :: ne_eig(:),nbands(:),nobd(:) !alloc in eigen_HF_init
END TYPE t_hybrid
TYPE prodtype
......@@ -460,8 +445,6 @@ MODULE m_types
REAL, ALLOCATABLE :: core1(:,:,:,:),core2(:,:,:,:) !alloc in eigen_HF_init
REAL, ALLOCATABLE :: eig_c(:,:,:) !alloc in eigen_HF_init
INTEGER , ALLOCATABLE:: kveclo_eig(:,:) !alloc in eigen_HF_setup
INTEGER,ALLOCATABLE :: ne_eig(:),nbands(:),nobd(:) !alloc in eigen_HF_init
INTEGER , ALLOCATABLE:: nbasm(:) !alloc in eigen_HF_init
INTEGER :: maxfac
REAL, ALLOCATABLE :: sfac(:),fac(:) !alloc in eigen_HF_init
REAL, ALLOCATABLE :: gauntarr(:,:,:,:,:,:) !alloc in eigen_HF_init
......
set(fleur_F90 ${fleur_F90}
hybrid/abcrot.F90 hybrid/exchange_core.F90 hybrid/mixedbasis.F90 hybrid/symmetrizeh.F90
hybrid/checkolap.F90 hybrid/exchange_val_hf.F90 hybrid/read_core.F90 hybrid/wavefproducts.F90
hybrid/coulombmatrix.F90 hybrid/gen_wavf.F90 hybrid/spmvec.F90
hybrid/eigen_HF_init.F90 hybrid/hsfock.F90 hybrid/subvxc.F90
hybrid/eigen_HF_setup.F90 hybrid/kp_perturbation.F90 hybrid/symm_hf.F90
hybrid/coulombmatrix.F90 hybrid/gen_wavf.F90 hybrid/spmvec.F90 hybrid/hybrid.F90
hybrid/HF_init.F90 hybrid/hsfock.F90 hybrid/subvxc.F90 hybrid/add_Vnonlocal.F90
hybrid/HF_setup.F90 hybrid/kp_perturbation.F90 hybrid/symm_hf.F90
)
MODULE m_hf_init
!
! preparations for HF and hybrid functional calculation
!
CONTAINS
SUBROUTINE hf_init(hybrid,kpts,atoms,input,DIMENSION,hybdat,irank2,isize2,l_real)
USE m_types
USE m_read_core
USE m_util
USE m_io_hybrid
IMPLICIT NONE
TYPE(t_hybrid),INTENT(INOUT) :: hybrid
TYPE(t_kpts),INTENT(IN) :: kpts
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_input),INTENT(IN) :: input
TYPE(t_dimension),INTENT(IN) :: DIMENSION
INTEGER,INTENT(OUT) :: irank2(:),isize2(:)
TYPE(t_hybdat),INTENT(OUT) :: hybdat
LOGICAL,INTENT(IN) :: l_real
INTEGER:: itype,ieq,l,m,i,nk,l1,l2,m1,m2,ok
!initialize hybdat%gridf for radial integration
CALL intgrf_init(atoms%ntype,atoms%jmtd,atoms%jri,atoms%dx,atoms%rmsh,hybdat%gridf)
!Alloc variables
ALLOCATE(hybdat%lmaxc(atoms%ntype))
ALLOCATE(hybrid%ne_eig(kpts%nkpt),hybrid%nbands(kpts%nkpt))
ALLOCATE(hybrid%nobd(kpts%nkptf))
ALLOCATE(hybdat%bas1(atoms%jmtd,hybrid%maxindx,0:atoms%lmaxd,atoms%ntype))
ALLOCATE(hybdat%bas2(atoms%jmtd,hybrid%maxindx,0:atoms%lmaxd,atoms%ntype))
ALLOCATE(hybdat%bas1_MT(hybrid%maxindx,0:atoms%lmaxd,atoms%ntype))
ALLOCATE(hybdat%drbas1_MT(hybrid%maxindx,0:atoms%lmaxd,atoms%ntype))
!sym%tau = oneD%ods%tau
! preparations for core states
CALL core_init(dimension,input,atoms, hybdat%lmaxcd,hybdat%maxindxc)
ALLOCATE( hybdat%nindxc(0:hybdat%lmaxcd,atoms%ntype), stat = ok )
IF( ok .ne. 0 ) STOP 'eigen_hf: failure allocation hybdat%nindxc'
ALLOCATE( hybdat%core1(atoms%jmtd,hybdat%maxindxc,0:hybdat%lmaxcd,atoms%ntype), stat=ok )
IF( ok .ne. 0 ) STOP 'eigen_hf: failure allocation core1'
ALLOCATE( hybdat%core2(atoms%jmtd,hybdat%maxindxc,0:hybdat%lmaxcd,atoms%ntype), stat=ok )
IF( ok .ne. 0 ) STOP 'eigen_hf: failure allocation core2'
ALLOCATE( hybdat%eig_c(hybdat%maxindxc,0:hybdat%lmaxcd,atoms%ntype), stat=ok )
IF( ok .ne. 0 ) STOP 'eigen_hf: failure allocation hybdat%eig_c'
hybdat%nindxc = 0 ; hybdat%core1 = 0 ; hybdat%core2 = 0 ; hybdat%eig_c = 0
! pre-calculate gaunt coefficients
hybdat%maxfac = max(2*atoms%lmaxd+hybrid%maxlcutm1+1,2*hybdat%lmaxcd+2*atoms%lmaxd+1)
ALLOCATE ( hybdat%fac( 0:hybdat%maxfac),hybdat%sfac( 0:hybdat%maxfac),stat=ok)
IF( ok .ne. 0 ) STOP 'eigen_hf: failure allocation fac,hybdat%sfac'
hybdat%fac(0) = 1
hybdat%sfac(0) = 1
DO i=1,hybdat%maxfac
hybdat%fac(i) = hybdat%fac(i-1)*i ! hybdat%fac(i) = i!
hybdat%sfac(i) = hybdat%sfac(i-1)*sqrt(i*1d0) ! hybdat%sfac(i) = sqrt(i!)
END DO
ALLOCATE ( hybdat%gauntarr( 2, 0:atoms%lmaxd, 0:atoms%lmaxd, 0:hybrid%maxlcutm1, -atoms%lmaxd:atoms%lmaxd ,-hybrid%maxlcutm1:hybrid%maxlcutm1 ),stat=ok)
IF( ok .ne. 0 ) STOP 'eigen: failure allocation hybdat%gauntarr'
hybdat%gauntarr = 0
DO l2 = 0,atoms%lmaxd
DO l1 = 0,atoms%lmaxd
DO l = abs(l1-l2),min(l1+l2,hybrid%maxlcutm1)
DO m = -l,l
DO m1 = -l1,l1
m2 = m1 + m ! Gaunt condition -m1+m2-m = 0
IF(abs(m2).le.l2) hybdat%gauntarr(1,l1,l2,l,m1,m) = gaunt(l1,l2,l,m1,m2,m,hybdat%maxfac,hybdat%fac,hybdat%sfac)
m2 = m1 - m ! switch role of l2-index
IF(abs(m2).le.l2) hybdat%gauntarr(2,l1,l2,l,m1,m) = gaunt(l2,l1,l,m2,m1,m,hybdat%maxfac,hybdat%fac,hybdat%sfac)
END DO
END DO
END DO
END DO
END DO
irank2 = 0
isize2 = 1
!skip_kpt = .false.
END SUBROUTINE hf_init
END MODULE m_hf_init
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_hybridmix
USE m_icorrkeys
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
! - local arrays -
REAL :: div_vv(hybrid%nbands(nk))
TYPE(t_mat) :: olap,tmp,v_x,z
COMPLEX :: exch(dimension%neigd,dimension%neigd)
!
! initialize weighting factor for HF exchange part
!
a_ex=get_exchange_weight(xcpot%icorr)
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
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) ###'
WRITE(6,'(A)') new_line('n') // ' k-point '// 'band tail pole input%total(valence+core)'
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)
! in the case of a spin-unpolarized calculation the factor 2 is added in eigen_hf.F TODO
if (.not.v_x%l_real) v_x%data_c=conjg(v_x%data_c)
exch = 0
print *,"sizes:",shape(z%data_r),shape(v_x%data_r)
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)')&
& kpts%bkf(:,nk),iband, (REAL(exch(iband,iband))-div_vv(iband))*(-27.211608),&
& div_vv(iband)*(-27.211608),REAL(exch(iband,iband))*(-27.211608)
END IF
END DO
END SUBROUTINE add_vnonlocal
END MODULE m_add_vnonlocal
......@@ -62,7 +62,7 @@
REAL :: q(3)
REAL :: integrand(atoms%jmtd)
REAL :: bkpt(3)
REAL :: rarr(maxval(hybdat%nbands))
REAL :: rarr(maxval(hybrid%nbands))
REAL :: rtaual(3)
REAL , ALLOCATABLE :: olapcb(:)
REAL , ALLOCATABLE :: olapcv_avg(:,:,:,:),olapcv_max(:,:,:,:)
......@@ -119,7 +119,7 @@
END DO
IF ( mpi%irank == 0 ) WRITE(6,'(/A)') ' Overlap <core|basis>'
ALLOCATE( olapcb(hybrid%maxindx),olapcv(maxval(hybdat%nbands),nkpti),&
ALLOCATE( olapcb(hybrid%maxindx),olapcv(maxval(hybrid%nbands),nkpti),&
& olapcv_avg( -hybdat%lmaxcd:hybdat%lmaxcd,hybdat%maxindxc,0:hybdat%lmaxcd,atoms%ntype),&
& olapcv_max( -hybdat%lmaxcd:hybdat%lmaxcd,hybdat%maxindxc,0:hybdat%lmaxcd,atoms%ntype),&
& olapcv_loc(2,-hybdat%lmaxcd:hybdat%lmaxcd,hybdat%maxindxc,0:hybdat%lmaxcd,atoms%ntype) )
......@@ -154,13 +154,13 @@
DO j=1,hybrid%nindx(l,itype)
lm = lm + 1
olapcv(:,:) = olapcv(:,:) + &
& olapcb(j)*cmt(:maxval(hybdat%nbands),lm,iatom,:nkpti)
& olapcb(j)*cmt(:maxval(hybrid%nbands),lm,iatom,:nkpti)
END DO
rdum = sum ( abs(olapcv(:,:))**2 )
rdum1 = maxval ( abs(olapcv(:,:)) )
iarr = maxloc ( abs(olapcv(:,:)) )
olapcv_avg( m,i,l,itype) = &
& sqrt( rdum / nkpti / sum(hybdat%nbands(:nkpti)) * nkpti )
& sqrt( rdum / nkpti / sum(hybrid%nbands(:nkpti)) * nkpti )
olapcv_max( m,i,l,itype) = rdum1
olapcv_loc(:,m,i,l,itype) = iarr
END DO
......@@ -233,9 +233,9 @@
IF ( mpi%irank == 0 ) WRITE(6,'(/A)') &
& 'Mismatch of wave functions at the MT-sphere boundaries'
ALLOCATE (carr1(maxval(hybdat%nbands),(atoms%lmaxd+1)**2))
ALLOCATE (carr2(maxval(hybdat%nbands),(atoms%lmaxd+1)**2))
ALLOCATE (carr3(maxval(hybdat%nbands),(atoms%lmaxd+1)**2))
ALLOCATE (carr1(maxval(hybrid%nbands),(atoms%lmaxd+1)**2))
ALLOCATE (carr2(maxval(hybrid%nbands),(atoms%lmaxd+1)**2))
ALLOCATE (carr3(maxval(hybrid%nbands),(atoms%lmaxd+1)**2))
DO ikpt = 1,nkpti
call read_z(z(ikpt),ikpt)
END DO
......@@ -280,7 +280,7 @@
cdum = 4*pi_const*img**l/sqrt(cell%omtil) * sphbes(l) * cexp
DO m = -l,l
lm = lm + 1
DO iband = 1,hybdat%nbands(ikpt)
DO iband = 1,hybrid%nbands(ikpt)
if (z(1)%l_real) THEN
carr2(iband,lm) = carr2(iband,lm) + cdum * z(ikpt)%data_r(igpt,iband) * y(lm)
Else
......@@ -300,7 +300,7 @@
DO n = 1,hybrid%nindx(l,itype)
lm1 = lm1 + 1
rdum = hybdat%bas1(atoms%jri(itype),n,l,itype) / atoms%rmt(itype)
DO iband = 1,hybdat%nbands(ikpt)
DO iband = 1,hybrid%nbands(ikpt)
carr3(iband,lm) = carr3(iband,lm) + cmt(iband,lm1,iatom,ikpt) * rdum
END DO
END DO
......
......@@ -31,7 +31,7 @@
! For the PW contribution we have to construct plane waves within the MT spheres with the help
! of spherical Bessel functions. The value lexp (LEXP in gwinp) is the corresponding cutoff.
!
MODULE m_coulomb
MODULE m_coulombmatrix
CONTAINS
......@@ -1886,7 +1886,7 @@ CONTAINS
!
CALL getshells(ptsh,nptsh,radsh,nshell,rrad,cell%bmat,first)
! minimum nonzero reciprocal-shell radius (needed in routines concerning the non-local hartree-fock exchange)
hybrid%radshmin = radsh(2)
!hybrid%radshmin = radsh(2)
!
! Fourier-space sum
!
......@@ -2231,4 +2231,4 @@ CONTAINS
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
END MODULE m_coulomb
END MODULE m_coulombmatrix
module m_eigen_hf_init
!
! preparations for HF and hybrid functional calculation
!
contains
subroutine eigen_hf_init(hybrid,kpts,atoms,input,dimension,hybdat,irank2,isize2,l_real)
USE m_types
USE m_read_core
USE m_util
USE m_io_hybrid
implicit none
TYPE(t_hybrid),INTENT(INOUT) :: hybrid
TYPE(t_kpts),INTENT(IN) :: kpts
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_input),INTENT(IN) :: input
TYPE(t_dimension),INTENT(IN) :: dimension
INTEGER,INTENT(OUT) :: irank2(:),isize2(:)
TYPE(t_hybdat),INTENT(OUT) :: hybdat
LOGICAL,INTENT(IN) :: l_real
LOGICAL,SAVE :: init_vex=.true.
INTEGER,SAVE :: nohf_it
integer:: itype,ieq,l,m,i,nk,l1,l2,m1,m2,ok
IF( .NOT. hybrid%l_hybrid ) THEN
hybrid%l_calhf = .false.
ELSE IF( (all(hybrid%ddist .lt. 1E-5) .or. init_vex .or. nohf_it >= 50 ) .and. input%imix .gt. 10) THEN
hybrid%l_calhf = .true.
init_vex = .false.
nohf_it = 0
ELSE IF( input%imix .lt. 10 ) THEN
hybrid%l_calhf = .true.
init_vex = .true.
nohf_it = 0
ELSE
hybrid%l_calhf = .false.
nohf_it = nohf_it + 1
END IF
IF( hybrid%l_calhf ) THEN
!initialize hybdat%gridf for radial integration