TABLE OF CONTENTS
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.
COPYRIGHT
Copyright (C) 2015-2018 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 .
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
PARENTS
forces,forstr,rhotov,setvtr
CHILDREN
libxc_functionals_end,libxc_functionals_init,metric,mkcore,rhotoxc xcdata_init
SOURCE
61 #if defined HAVE_CONFIG_H 62 #include "config.h" 63 #endif 64 65 #include "abi_common.h" 66 67 subroutine xchybrid_ncpp_cc(dtset,enxc,mpi_enreg,nfft,ngfft,n3xccc,rhor,rprimd,strsxc,vxcavg,xccc3d,vxc,grxc,xcccrc,xccc1d,& 68 & xred,n1xccc,optstr) 69 70 use defs_basis 71 use m_profiling_abi 72 use m_errors 73 use m_xcdata 74 use libxc_functionals 75 76 use m_geometry, only : metric 77 use defs_abitypes, only : MPI_type, dataset_type 78 use m_dtset, only : dtset_copy, dtset_free 79 80 !This section has been created automatically by the script Abilint (TD). 81 !Do not modify the following lines by hand. 82 #undef ABI_FUNC 83 #define ABI_FUNC 'xchybrid_ncpp_cc' 84 use interfaces_56_xc, except_this_one => xchybrid_ncpp_cc 85 !End of the abilint section 86 87 implicit none 88 89 !Arguments ------------------------------------------------------------- 90 !scalars 91 integer,intent(in) :: nfft,n3xccc 92 integer,optional,intent(in) :: n1xccc,optstr 93 real(dp),intent(out) :: enxc,vxcavg 94 type(dataset_type),intent(in) :: dtset 95 type(MPI_type),intent(in) :: mpi_enreg 96 !arrays 97 integer,intent(in) :: ngfft(18) 98 real(dp),intent(in) :: rhor(nfft,dtset%nspden),rprimd(3,3),xccc3d(n3xccc) 99 real(dp),intent(out) :: strsxc(6) 100 real(dp),optional,intent(in) :: xcccrc(dtset%ntypat),xred(3,dtset%natom),xccc1d(:,:,:) 101 real(dp),optional,intent(out) :: grxc(3,dtset%natom),vxc(nfft,dtset%nspden) 102 103 !Local variables ------------------------------------------------------- 104 !scalars 105 integer :: ixc_gga,libxc_gga_initialized,ndim,nkxc,n3xccc_null,option,optstr_loc,usexcnhat 106 real(dp) :: enxc_corr,ucvol,vxcavg_corr 107 character(len=500) :: msg 108 type(xcdata_type) :: xcdata_gga,xcdata_hybrid 109 logical :: calcgrxc,nmxc 110 !arrays 111 integer :: gga_id(2) 112 real(dp) :: nhat(1,0),nhatgr(1,1,0),strsxc_corr(6),gmet(3,3),gprimd(3,3),rmet(3,3) 113 real(dp),allocatable :: kxc_dum(:,:),vxc_corr(:,:),xccc3d_null(:),dyfrx2_dum(:,:,:) 114 type(libxc_functional_type) :: xc_funcs_gga(2) 115 116 ! ************************************************************************* 117 118 DBG_ENTER("COLL") 119 120 121 nmxc=(dtset%usepawu==4).or.(dtset%usepawu==14) 122 123 !Not relevant for PAW 124 if (dtset%usepaw==1) return 125 if(n3xccc==0) return 126 calcgrxc=(present(grxc).and.present(n1xccc).and.present(xcccrc).and.present(xred).and.present(xccc1d)) 127 optstr_loc=0 128 if(present(optstr)) optstr_loc=1 129 !Not applicable for electron-positron 130 if (dtset%positron>0) then 131 msg='NCPP+Hybrid functionals not applicable for electron-positron calculations!' 132 MSG_ERROR(msg) 133 end if 134 135 !Select the GGA on which the hybrid functional is based on 136 !or return if not applicable 137 if (dtset%ixc==41.or.dtset%ixc==42) then 138 ixc_gga = 11 139 else if (dtset%ixc<0) then 140 if (libxc_functionals_gga_from_hybrid(gga_id=gga_id)) then 141 ixc_gga=-gga_id(2)*1000-gga_id(1) 142 else 143 return 144 end if 145 else 146 return 147 end if 148 149 !Define xcdata_hybrid as well as xcdata_gga 150 call xcdata_init(xcdata_hybrid,dtset=dtset) 151 call xcdata_init(xcdata_gga,dtset=dtset,auxc_ixc=0,ixc=ixc_gga) 152 libxc_gga_initialized=0 153 154 nkxc=0;ndim=0;usexcnhat=0;n3xccc_null=0 155 ABI_ALLOCATE(kxc_dum,(nfft,nkxc)) 156 ABI_ALLOCATE(vxc_corr,(nfft,dtset%nspden)) 157 158 if (present(vxc).and.optstr_loc==0) then 159 !Initialize args for rhotoxc 160 option=0 ! XC only 161 ABI_ALLOCATE(xccc3d_null,(n3xccc_null)) 162 ! Compute Vxc^Hybrid(rho_val) 163 call rhotoxc(enxc,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 164 & n3xccc_null,option,dtset%paral_kgb,rhor,rprimd,& 165 & strsxc,usexcnhat,vxc,vxcavg,xccc3d_null,xcdata_hybrid) 166 167 ! Initialize GGA functional 168 if (ixc_gga<0) then 169 call libxc_functionals_init(ixc_gga,dtset%nspden,xc_functionals=xc_funcs_gga) 170 libxc_gga_initialized=1 171 end if 172 173 !Add Vxc^GGA(rho_core+rho_val) 174 if (ixc_gga<0) then 175 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 176 & n3xccc,option,dtset%paral_kgb,rhor,rprimd,& 177 & strsxc_corr,usexcnhat,vxc_corr,vxcavg_corr,xccc3d,xcdata_gga,xc_funcs=xc_funcs_gga) 178 else 179 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 180 & n3xccc,option,dtset%paral_kgb,rhor,rprimd,& 181 & strsxc_corr,usexcnhat,vxc_corr,vxcavg_corr,xccc3d,xcdata_gga) 182 end if 183 enxc=enxc+enxc_corr 184 vxc(:,:)=vxc(:,:)+vxc_corr(:,:) 185 vxcavg=vxcavg+vxcavg_corr 186 strsxc(:)=strsxc(:)+strsxc_corr(:) 187 188 !Substract Vxc^GGA(rho_val) 189 if (ixc_gga<0) then 190 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 191 & n3xccc_null,option,dtset%paral_kgb,rhor,rprimd,strsxc_corr,usexcnhat,& 192 & vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga,xc_funcs=xc_funcs_gga) 193 else 194 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 195 & n3xccc_null,option,dtset%paral_kgb,rhor,rprimd,strsxc_corr,usexcnhat,& 196 & vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga) 197 end if 198 enxc=enxc-enxc_corr 199 vxc(:,:)=vxc(:,:)-vxc_corr(:,:) 200 vxcavg=vxcavg-vxcavg_corr 201 strsxc(:)=strsxc(:)-strsxc_corr(:) 202 203 !Release memory 204 ABI_DEALLOCATE(xccc3d_null) 205 end if 206 207 if (calcgrxc) then 208 209 call metric(gmet,gprimd,-1,rmet,rprimd,ucvol) 210 ! Initialize GGA functional 211 if (ixc_gga<0 .and. libxc_gga_initialized==0) then 212 call libxc_functionals_init(ixc_gga,dtset%nspden,xc_functionals=xc_funcs_gga) 213 libxc_gga_initialized=1 214 end if 215 ABI_ALLOCATE(dyfrx2_dum,(3,3,dtset%natom)) 216 ABI_ALLOCATE(xccc3d_null,(n3xccc)) 217 option=1 218 ! calculate xccc3d in this case 219 call mkcore(strsxc_corr,dyfrx2_dum,grxc,mpi_enreg,dtset%natom,nfft,dtset%nspden,dtset%ntypat,ngfft(1),n1xccc,ngfft(2),& 220 & ngfft(3),option,rprimd,dtset%typat,ucvol,vxc_corr,xcccrc,xccc1d,xccc3d_null,xred) 221 !Add Vxc^GGA(rho_core+rho_val) 222 if (ixc_gga<0) then 223 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 224 & n3xccc,option,dtset%paral_kgb,& 225 & rhor,rprimd,strsxc_corr,usexcnhat,vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga,xc_funcs=xc_funcs_gga) 226 else 227 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 228 & n3xccc,option,dtset%paral_kgb,& 229 & rhor,rprimd,strsxc_corr,usexcnhat,vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga) 230 end if 231 option=2 232 call mkcore(strsxc_corr,dyfrx2_dum,grxc,mpi_enreg,dtset%natom,nfft,dtset%nspden,dtset%ntypat,ngfft(1),n1xccc,ngfft(2),& 233 & ngfft(3),option,rprimd,dtset%typat,ucvol,vxc_corr,xcccrc,xccc1d,xccc3d_null,xred) 234 ABI_DEALLOCATE(dyfrx2_dum) 235 ABI_DEALLOCATE(xccc3d_null) 236 end if 237 238 if(optstr_loc==1) then 239 ! Initialize GGA functional 240 if (ixc_gga<0 .and. libxc_gga_initialized==0) then 241 call libxc_functionals_init(ixc_gga,dtset%nspden,xc_functionals=xc_funcs_gga) 242 libxc_gga_initialized=1 243 end if 244 !calculate Vxc^GGA(rho_core+rho_val) 245 option=0 246 if (ixc_gga<0) then 247 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 248 & n3xccc,option,dtset%paral_kgb,rhor,rprimd,& 249 & strsxc_corr,usexcnhat,vxc,vxcavg_corr,xccc3d,xcdata_gga,xc_funcs=xc_funcs_gga) 250 else 251 call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,& 252 & n3xccc,option,dtset%paral_kgb,rhor,rprimd,& 253 & strsxc_corr,usexcnhat,vxc,vxcavg_corr,xccc3d,xcdata_gga) 254 end if 255 end if 256 257 ! Suppress the temporary used xc functional 258 if(libxc_gga_initialized==1) then 259 call libxc_functionals_end(xc_functionals=xc_funcs_gga) 260 end if 261 ABI_DEALLOCATE(vxc_corr) 262 ABI_DEALLOCATE(kxc_dum) 263 264 DBG_EXIT("COLL") 265 266 end subroutine xchybrid_ncpp_cc