TABLE OF CONTENTS
ABINIT/hybrid_corr [ Functions ]
NAME
hybrid_corr
FUNCTION
Compute the correction to the XC potential and energy due to the hybridation
INPUTS
dtset <type(dataset_type)>=all input variables in this dataset ixc = choice of exchange-correlation functional. mpi_enreg=information about MPI parallelization nfft = number of fft grid points. ngfft(1:3) = integer fft box dimensions, see getng for ngfft(4:8). nspden = number of spin-density components. rhor(nfft,nspden) = electron density in real space in electrons/bohr**3 (total in first half and spin-up in second half if nspden = 2). rprimd(3,3) = dimensional primitive translations for real space in Bohr.
TODO
This routine is not used and does no compile
OUTPUT
SIDE EFFECTS
enxc = exchange correlation energy corrected for hybrid functionals vxc = exchange correlation potential corrected for hybrid functionals
WARNINGS
Current restrictions are: a - Spin-polarized case not tested.
SOURCE
308 subroutine hybrid_corr(dtset,ixc,nkxc,mpi_enreg,nfft,ngfft,nspden,rhor,rprimd,hybrid_mixing,vxc,enxc) 309 310 !Arguments ------------------------------------------------------------- 311 !scalars 312 integer,intent(in) :: ixc,nfft,nspden,nkxc 313 real(dp),intent(in) :: hybrid_mixing 314 real(dp),intent(inout) :: enxc 315 type(MPI_type),intent(in) :: mpi_enreg 316 type(dataset_type),intent(in) :: dtset 317 !arrays 318 integer,intent(in) :: ngfft(18) 319 real(dp),intent(in) :: rhor(nfft,nspden),rprimd(3,3) 320 real(dp),intent(inout) :: vxc(nfft,nspden) 321 322 !Local variables ------------------------------------------------------- 323 !No improved xc quadrature. 324 !No core correction. 325 !Dummy here. 326 !For debugging purposes (see tests below): 327 !integer :: i1,i2,i3,k1,n1,n2,n3 328 !real(dp) :: kx,rho,rhomax,ftest 329 !scalars 330 logical :: nmxc 331 character(len=500) :: message 332 type(dataset_type) :: dtLocal 333 !arrays 334 real(dp) :: dum(0) 335 real(dp) :: enxc_corr 336 real(dp),allocatable :: kxcr(:,:) 337 real(dp),allocatable :: xccc3d(:),vxc_corr(:,:) 338 339 !*********************************************************************** 340 !For debugging purposes (see tests below): 341 342 !ftest(i1,n1,k1) = 0._dp+1._dp*cos(k1*two_pi*float(i1)/float(n1)) 343 344 !*********************************************************************** 345 346 !Check input parameters. 347 348 if (nspden > 2) then 349 message = ' kxc_alda does not work yet for nspden > 2.' 350 ABI_ERROR(message) 351 end if 352 353 !Copy the input variables from the current dataset to a temporary one 354 !to tune some parameters 355 call dtset_copy(dtLocal, dtset) 356 dtLocal%intxc = 0 357 dtLocal%ixc = -101 358 nmxc=(dtset%usepaw==1.and.mod(abs(dtset%usepawu),10)==4) 359 360 ! Reinitialize the libxc module with the overriden values 361 if (dtset%ixc<0) then 362 call libxc_functionals_end() 363 end if 364 if (dtLocal%ixc<0) then 365 call libxc_functionals_init(dtLocal%ixc,dtLocal%nspden) 366 end if 367 368 call rhohxc(dtLocal,enxc_corr,dum,0,kxcr,mpi_enreg,nfft,dum,dum,0,dum,0,nkxc,0,& 369 & nmxc,nspden,n3xccc,0,dum,rhor,rprimd,strsxc,1,dum,vxc_corr,dum,xccc3d) 370 371 vxc(:,:) = vxc(:,:) + hybrid_mixing*vxc_corr(:,:) 372 enxc = enxc + hybrid_mixing*enxc_corr 373 374 call rhohxc(dtLocal,enxc_corr,dum,0,kxcr,mpi_enreg,dum,dum,dum,0,dum,0,nkxc,0,& 375 & nmxc,nspden,0,0,dum,rhor,rprimd,strsxc,1,dum,vxc_corr,vxcavg,0) 376 377 vxc(:,:) = vxc(:,:) - hybrid_mixing*vxc_corr(:,:) 378 enxc = enxc - hybrid_mixing*enxc_corr 379 380 381 ! Revert libxc module to the original settings 382 if (dtLocal%ixc<0) then 383 call libxc_functionals_end() 384 end if 385 if (dtset%ixc<0) then 386 call libxc_functionals_init(dtset%ixc,dtset%nspden) 387 end if 388 389 !Free memory. 390 call dtset_free(dtLocal) 391 ABI_FREE(kxcr) 392 ABI_FREE(xccc3d) 393 394 end subroutine hybrid_corr
ABINIT/m_xchybrid [ Modules ]
NAME
m_xchybrid
FUNCTION
COPYRIGHT
Copyright (C) 2015-2024 ABINIT group (FA,MT,FJ) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
15 #if defined HAVE_CONFIG_H 16 #include "config.h" 17 #endif 18 19 #include "abi_common.h" 20 21 module m_xchybrid 22 23 use defs_basis 24 use m_abicore 25 use m_errors 26 use m_xcdata 27 use libxc_functionals 28 use m_dtset 29 30 use m_geometry, only : metric 31 use defs_abitypes, only : MPI_type 32 use m_rhotoxc, only : rhotoxc 33 use m_mkcore, only : mkcore 34 35 implicit none 36 37 private
ABINIT/xchybrid_ncpp_cc [ Functions ]
NAME
xchybrid_ncpp_cc
FUNCTION
XC Hybrid Norm-Conserving PseudoPotential Core Correction: Relevant only for Norm-Conserving PseudoPotentials (NCPP) and for hybrid functionals. Compute the correction to the XC energy/potential due to the lack of core wave-functions: As Fock exchange cannot be computed for core-core and core-valence interactions, these contribution have to be also substracted from the GGA exchange-correlation.
INPUTS
dtset <type(dataset_type)>= all input variables in this dataset mpi_enreg= information about MPI parallelization nfft= number of fft grid points. ngfft(1:3)= integer fft box dimensions, see getng for ngfft(4:8). n1xccc=dimension of xccc1d ; 0 if no XC core correction is used n3xccc=dimension of the xccc3d array (0 or nfft). optstr= calculate corrected vxc if optstr=1 rhor(nfft,nspden)= electron density in real space in electrons/bohr**3 rprimd(3,3)= dimensional primitive translations for real space in Bohr. xccc3d(n3xccc)=3D core electron density for XC core correction (bohr^-3) xcccrc(ntypat)=XC core correction cutoff radius (bohr) for each atom type xccc1d(n1xccc*(1-usepaw),6,ntypat)=1D core charge function and five derivatives, for each type of atom, from psp (used in Norm-conserving only) xred(3,natom)=reduced dimensionless atomic coordinates
OUTPUT
SIDE EFFECTS
enxc= exchange correlation energy grxc= correction to the forces strsxc(6)= exchange correlation contribution to stress tensor vxc= exchange correlation potential vxcavg= unit cell average of Vxc
NOTES
The final expression of the XC potential (that should be added to alpha*VFock[Psi_val]) is: Vxc=Vx[rho_core+rho_val] - alpha*Vx[rho_val] + Vc[rho_core+rho_val] To accomodate libXC convention, Vxc is computed as follows: Vxc=Vxc_libXC[rho_val] + Vxc_gga[rho_core+rho_val] - Vxc_gga[rho_val] Note that this is equivalent to Vxc=Vx_libXC[rho_val] + Vxc_gga[rho_core+rho_val] - Vx_gga[rho_val] but needs one less call to libxc
SOURCE
93 subroutine xchybrid_ncpp_cc(dtset,enxc,mpi_enreg,nfft,ngfft,n3xccc,rhor,rprimd,strsxc,vxcavg,xccc3d,vxc,grxc,xcccrc,xccc1d,& 94 & xred,n1xccc,optstr) 95 96 !Arguments ------------------------------------------------------------- 97 !scalars 98 integer,intent(in) :: nfft,n3xccc 99 integer,optional,intent(in) :: n1xccc,optstr 100 real(dp),intent(out) :: enxc,vxcavg 101 type(dataset_type),intent(in) :: dtset 102 type(MPI_type),intent(in) :: mpi_enreg 103 !arrays 104 integer,intent(in) :: ngfft(18) 105 real(dp),intent(in) :: rhor(nfft,dtset%nspden),rprimd(3,3),xccc3d(n3xccc) 106 real(dp),intent(out) :: strsxc(6) 107 real(dp),optional,intent(in) :: xcccrc(dtset%ntypat),xred(3,dtset%natom),xccc1d(:,:,:) 108 real(dp),optional,intent(out) :: grxc(3,dtset%natom),vxc(nfft,dtset%nspden) 109 110 !Local variables ------------------------------------------------------- 111 !scalars 112 integer :: ixc_gga,libxc_gga_initialized,ndim,nkxc,n3xccc_null,option,optstr_loc,usexcnhat 113 real(dp) :: enxc_corr,ucvol,vxcavg_corr 114 character(len=500) :: msg 115 type(xcdata_type) :: xcdata_gga,xcdata_hybrid 116 logical :: calcgrxc,nmxc 117 !arrays 118 integer :: gga_id(2) 119 real(dp) :: nhat(1,0),nhatgr(1,1,0),strsxc_corr(6),gmet(3,3),gprimd(3,3),rmet(3,3) 120 real(dp),allocatable :: kxc_dum(:,:),vxc_corr(:,:),xccc3d_null(:),dyfrx2_dum(:,:,:) 121 type(libxc_functional_type) :: xc_funcs_gga(2) 122 123 ! ************************************************************************* 124 125 DBG_ENTER("COLL") 126 127 !Not relevant for PAW 128 if (dtset%usepaw==1) return 129 if(n3xccc==0) return 130 calcgrxc=(present(grxc).and.present(n1xccc).and.present(xcccrc).and.present(xred).and.present(xccc1d)) 131 optstr_loc=0 132 if(present(optstr)) optstr_loc=1 133 !Not applicable for electron-positron 134 if (dtset%positron>0) then 135 msg='NCPP+Hybrid functionals not applicable for electron-positron calculations!' 136 ABI_ERROR(msg) 137 end if 138 139 !Select the GGA on which the hybrid functional is based on 140 !or return if not applicable 141 if (dtset%ixc==41.or.dtset%ixc==42) then 142 ixc_gga = 11 143 else if (dtset%ixc<0) then 144 if (libxc_functionals_gga_from_hybrid(gga_id=gga_id)) then 145 ixc_gga=-gga_id(2)*1000-gga_id(1) 146 else 147 return 148 end if 149 else 150 return 151 end if 152 153 !Define xcdata_hybrid as well as xcdata_gga 154 call xcdata_init(xcdata_hybrid,dtset=dtset) 155 call xcdata_init(xcdata_gga,dtset=dtset,auxc_ixc=0,ixc=ixc_gga) 156 libxc_gga_initialized=0 ; nmxc=.false. 157 158 nkxc=0;ndim=0;usexcnhat=0;n3xccc_null=0 159 ABI_MALLOC(kxc_dum,(nfft,nkxc)) 160 ABI_MALLOC(vxc_corr,(nfft,dtset%nspden)) 161 162 if (present(vxc).and.optstr_loc==0) then 163 !Initialize args for rhotoxc 164 option=0 ! XC only 165 ABI_MALLOC(xccc3d_null,(n3xccc_null)) 166 ! Compute Vxc^Hybrid(rho_val) 167 call rhotoxc(enxc,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 168 & n3xccc_null,option,rhor,rprimd,& 169 & strsxc,usexcnhat,vxc,vxcavg,xccc3d_null,xcdata_hybrid) 170 171 ! Initialize GGA functional 172 if (ixc_gga<0) then 173 call libxc_functionals_init(ixc_gga,dtset%nspden,xc_functionals=xc_funcs_gga) 174 libxc_gga_initialized=1 175 end if 176 177 !Add Vxc^GGA(rho_core+rho_val) 178 if (ixc_gga<0) then 179 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 180 & n3xccc,option,rhor,rprimd,& 181 & strsxc_corr,usexcnhat,vxc_corr,vxcavg_corr,xccc3d,xcdata_gga,xc_funcs=xc_funcs_gga) 182 else 183 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 184 & n3xccc,option,rhor,rprimd,& 185 & strsxc_corr,usexcnhat,vxc_corr,vxcavg_corr,xccc3d,xcdata_gga) 186 end if 187 enxc=enxc+enxc_corr 188 vxc(:,:)=vxc(:,:)+vxc_corr(:,:) 189 vxcavg=vxcavg+vxcavg_corr 190 strsxc(:)=strsxc(:)+strsxc_corr(:) 191 192 !Substract Vxc^GGA(rho_val) 193 if (ixc_gga<0) then 194 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 195 & n3xccc_null,option,rhor,rprimd,strsxc_corr,usexcnhat,& 196 & vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga,xc_funcs=xc_funcs_gga) 197 else 198 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 199 & n3xccc_null,option,rhor,rprimd,strsxc_corr,usexcnhat,& 200 & vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga) 201 end if 202 enxc=enxc-enxc_corr 203 vxc(:,:)=vxc(:,:)-vxc_corr(:,:) 204 vxcavg=vxcavg-vxcavg_corr 205 strsxc(:)=strsxc(:)-strsxc_corr(:) 206 207 !Release memory 208 ABI_FREE(xccc3d_null) 209 end if 210 211 if (calcgrxc) then 212 213 call metric(gmet,gprimd,-1,rmet,rprimd,ucvol) 214 ! Initialize GGA functional 215 if (ixc_gga<0 .and. libxc_gga_initialized==0) then 216 call libxc_functionals_init(ixc_gga,dtset%nspden,xc_functionals=xc_funcs_gga) 217 libxc_gga_initialized=1 218 end if 219 ABI_MALLOC(dyfrx2_dum,(3,3,dtset%natom)) 220 ABI_MALLOC(xccc3d_null,(n3xccc)) 221 option=1 222 ! calculate xccc3d in this case 223 call mkcore(strsxc_corr,dyfrx2_dum,grxc,mpi_enreg,dtset%natom,nfft,dtset%nspden,dtset%ntypat,ngfft(1),n1xccc,ngfft(2),& 224 & ngfft(3),option,rprimd,dtset%typat,ucvol,vxc_corr,xcccrc,xccc1d,xccc3d_null,xred) 225 !Add Vxc^GGA(rho_core+rho_val) 226 if (ixc_gga<0) then 227 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 228 & n3xccc,option,& 229 & rhor,rprimd,strsxc_corr,usexcnhat,vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga,xc_funcs=xc_funcs_gga) 230 else 231 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 232 & n3xccc,option,& 233 & rhor,rprimd,strsxc_corr,usexcnhat,vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga) 234 end if 235 option=2 236 call mkcore(strsxc_corr,dyfrx2_dum,grxc,mpi_enreg,dtset%natom,nfft,dtset%nspden,dtset%ntypat,ngfft(1),n1xccc,ngfft(2),& 237 & ngfft(3),option,rprimd,dtset%typat,ucvol,vxc_corr,xcccrc,xccc1d,xccc3d_null,xred) 238 ABI_FREE(dyfrx2_dum) 239 ABI_FREE(xccc3d_null) 240 end if 241 242 if(optstr_loc==1) then 243 ! Initialize GGA functional 244 if (ixc_gga<0 .and. libxc_gga_initialized==0) then 245 call libxc_functionals_init(ixc_gga,dtset%nspden,xc_functionals=xc_funcs_gga) 246 libxc_gga_initialized=1 247 end if 248 !calculate Vxc^GGA(rho_core+rho_val) 249 option=0 250 if (ixc_gga<0) then 251 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 252 & n3xccc,option,rhor,rprimd,& 253 & strsxc_corr,usexcnhat,vxc,vxcavg_corr,xccc3d,xcdata_gga,xc_funcs=xc_funcs_gga) 254 else 255 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 256 & n3xccc,option,rhor,rprimd,& 257 & strsxc_corr,usexcnhat,vxc,vxcavg_corr,xccc3d,xcdata_gga) 258 end if 259 end if 260 261 ! Suppress the temporary used xc functional 262 if(libxc_gga_initialized==1) then 263 call libxc_functionals_end(xc_functionals=xc_funcs_gga) 264 end if 265 ABI_FREE(vxc_corr) 266 ABI_FREE(kxc_dum) 267 268 DBG_EXIT("COLL") 269 270 end subroutine xchybrid_ncpp_cc