Line data Source code
1 : !!****m* ABINIT/m_xchybrid
2 : !! NAME
3 : !! m_xchybrid
4 : !!
5 : !! FUNCTION
6 : !!
7 : !! COPYRIGHT
8 : !! Copyright (C) 2015-2026 ABINIT group (FA,MT,FJ)
9 : !! This file is distributed under the terms of the
10 : !! GNU General Public License, see ~abinit/COPYING
11 : !! or http://www.gnu.org/copyleft/gpl.txt .
12 : !!
13 : !! SOURCE
14 :
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
38 : !!***
39 :
40 : public :: xchybrid_ncpp_cc
41 : !!***
42 :
43 : contains
44 : !!***
45 :
46 : !!****f* ABINIT/xchybrid_ncpp_cc
47 : !! NAME
48 : !! xchybrid_ncpp_cc
49 : !!
50 : !! FUNCTION
51 : !! XC Hybrid Norm-Conserving PseudoPotential Core Correction:
52 : !! Relevant only for Norm-Conserving PseudoPotentials (NCPP) and for hybrid functionals.
53 : !! Compute the correction to the XC energy/potential due to the lack of core wave-functions:
54 : !! As Fock exchange cannot be computed for core-core and core-valence interactions, these
55 : !! contribution have to be also substracted from the GGA exchange-correlation.
56 : !!
57 : !! INPUTS
58 : !! dtset <type(dataset_type)>= all input variables in this dataset
59 : !! mpi_enreg= information about MPI parallelization
60 : !! nfft= number of fft grid points.
61 : !! ngfft(1:3)= integer fft box dimensions, see getng for ngfft(4:8).
62 : !! n1xccc=dimension of xccc1d ; 0 if no XC core correction is used
63 : !! n3xccc=dimension of the xccc3d array (0 or nfft).
64 : !! optstr= calculate corrected vxc if optstr=1
65 : !! rhor(nfft,nspden)= electron density in real space in electrons/bohr**3
66 : !! rprimd(3,3)= dimensional primitive translations for real space in Bohr.
67 : !! xccc3d(n3xccc)=3D core electron density for XC core correction (bohr^-3)
68 : !! xcccrc(ntypat)=XC core correction cutoff radius (bohr) for each atom type
69 : !! xccc1d(n1xccc*(1-usepaw),6,ntypat)=1D core charge function and five derivatives,
70 : !! for each type of atom, from psp (used in Norm-conserving only)
71 : !! xred(3,natom)=reduced dimensionless atomic coordinates
72 : !!
73 : !! OUTPUT
74 : !!
75 : !! SIDE EFFECTS
76 : !! bigexc= exchange correlation energy
77 : !! bigsxc= exchange correlation entropy energy (for finite-temperature xc functionals)
78 : !! grxc= correction to the forces
79 : !! strsxc(6)= exchange correlation contribution to stress tensor
80 : !! vxc= exchange correlation potential
81 : !! vxcavg= unit cell average of Vxc
82 : !!
83 : !! NOTES
84 : !! The final expression of the XC potential (that should be added to alpha*VFock[Psi_val]) is:
85 : !! Vxc=Vx[rho_core+rho_val] - alpha*Vx[rho_val] + Vc[rho_core+rho_val]
86 : !! To accomodate libXC convention, Vxc is computed as follows:
87 : !! Vxc=Vxc_libXC[rho_val] + Vxc_gga[rho_core+rho_val] - Vxc_gga[rho_val]
88 : !! Note that this is equivalent to
89 : !! Vxc=Vx_libXC[rho_val] + Vxc_gga[rho_core+rho_val] - Vx_gga[rho_val]
90 : !! but needs one less call to libxc
91 : !!
92 : !! SOURCE
93 :
94 986 : subroutine xchybrid_ncpp_cc(dtset,bigexc,bigsxc,mpi_enreg,nfft,ngfft,n3xccc,rhor,rprimd,strsxc,&
95 886 : & vxcavg,xccc3d,vxc,grxc,xcccrc,xccc1d,xred,n1xccc,optstr)
96 :
97 : !Arguments -------------------------------------------------------------
98 : !scalars
99 : integer,intent(in) :: nfft,n3xccc
100 : integer,optional,intent(in) :: n1xccc,optstr
101 : real(dp),intent(out) :: bigexc,vxcavg,bigsxc
102 : type(dataset_type),intent(in) :: dtset
103 : type(MPI_type),intent(in) :: mpi_enreg
104 : !arrays
105 : integer,intent(in) :: ngfft(18)
106 : real(dp),intent(in) :: rhor(nfft,dtset%nspden),rprimd(3,3),xccc3d(n3xccc)
107 : real(dp),intent(out) :: strsxc(6)
108 : real(dp),optional,intent(in) :: xcccrc(dtset%ntypat),xred(3,dtset%natom),xccc1d(:,:,:)
109 : real(dp),optional,intent(out) :: grxc(3,dtset%natom),vxc(nfft,dtset%nspden)
110 :
111 : !Local variables -------------------------------------------------------
112 : !scalars
113 : integer :: ixc_gga,libxc_gga_initialized,ndim,nkxc,n3xccc_null,option,optstr_loc,usexcnhat
114 : real(dp) :: bigexc_corr,ucvol,vxcavg_corr,bigsxc_corr
115 : character(len=500) :: msg
116 : type(xcdata_type) :: xcdata_gga,xcdata_hybrid
117 : logical :: calcgrxc,nmxc
118 : !arrays
119 : integer :: gga_id(2)
120 : real(dp) :: nhat(1,0),nhatgr(1,1,0),strsxc_corr(6),gmet(3,3),gprimd(3,3),rmet(3,3)
121 493 : real(dp),allocatable :: kxc_dum(:,:),vxc_corr(:,:),xccc3d_null(:),dyfrx2_dum(:,:,:)
122 1479 : type(libxc_functional_type) :: xc_funcs_gga(2)
123 :
124 : ! *************************************************************************
125 :
126 : DBG_ENTER("COLL")
127 :
128 : !Not relevant for PAW
129 493 : if (dtset%usepaw==1) return
130 493 : if(n3xccc==0) return
131 493 : calcgrxc=(present(grxc).and.present(n1xccc).and.present(xcccrc).and.present(xred).and.present(xccc1d))
132 493 : optstr_loc=0
133 493 : if(present(optstr)) optstr_loc=1
134 : !Not applicable for electron-positron
135 493 : if (dtset%positron>0) then
136 0 : msg='NCPP+Hybrid functionals not applicable for electron-positron calculations!'
137 0 : ABI_ERROR(msg)
138 : end if
139 :
140 : !Select the GGA on which the hybrid functional is based on
141 : !or return if not applicable
142 493 : if (dtset%ixc==41.or.dtset%ixc==42) then
143 199 : ixc_gga = 11
144 294 : else if (dtset%ixc<0) then
145 294 : if (libxc_functionals_gga_from_hybrid(gga_id=gga_id)) then
146 294 : ixc_gga=-gga_id(2)*1000-gga_id(1)
147 : else
148 : return
149 : end if
150 : else
151 : return
152 : end if
153 :
154 : !Define xcdata_hybrid as well as xcdata_gga
155 493 : call xcdata_init(xcdata_hybrid,dtset=dtset)
156 493 : call xcdata_init(xcdata_gga,dtset=dtset,auxc_ixc=0,ixc=ixc_gga)
157 493 : libxc_gga_initialized=0 ; nmxc=.false.
158 :
159 493 : nkxc=0;ndim=0;usexcnhat=0;n3xccc_null=0
160 986 : ABI_MALLOC(kxc_dum,(nfft,nkxc))
161 1972 : ABI_MALLOC(vxc_corr,(nfft,dtset%nspden))
162 :
163 493 : if (present(vxc).and.optstr_loc==0) then
164 : !Initialize args for rhotoxc
165 363 : option=0 ! XC only
166 363 : ABI_MALLOC(xccc3d_null,(n3xccc_null))
167 : ! Compute Vxc^Hybrid(rho_val)
168 : call rhotoxc(bigexc,bigsxc,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
169 : & n3xccc_null,option,rhor,rprimd,usexcnhat,vxc,vxcavg,xccc3d_null,xcdata_hybrid,&
170 363 : & strsxc=strsxc)
171 :
172 : ! Initialize GGA functional
173 363 : if (ixc_gga<0) then
174 219 : call libxc_functionals_init(ixc_gga,dtset%nspden,xc_functionals=xc_funcs_gga)
175 219 : libxc_gga_initialized=1
176 : end if
177 :
178 : !Add Vxc^GGA(rho_core+rho_val)
179 363 : if (ixc_gga<0) then
180 : call rhotoxc(bigexc_corr,bigsxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
181 : & n3xccc,option,rhor,rprimd,usexcnhat,vxc_corr,vxcavg_corr,xccc3d,xcdata_gga,&
182 219 : & xc_funcs=xc_funcs_gga,strsxc=strsxc_corr)
183 : else
184 : call rhotoxc(bigexc_corr,bigsxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
185 : & n3xccc,option,rhor,rprimd,usexcnhat,vxc_corr,vxcavg_corr,xccc3d,xcdata_gga,&
186 144 : & strsxc=strsxc_corr)
187 : end if
188 363 : bigexc=bigexc+bigexc_corr
189 363 : bigsxc=bigsxc+bigsxc_corr
190 3848990 : vxc(:,:)=vxc(:,:)+vxc_corr(:,:)
191 363 : vxcavg=vxcavg+vxcavg_corr
192 2541 : strsxc(:)=strsxc(:)+strsxc_corr(:)
193 :
194 : !Substract Vxc^GGA(rho_val)
195 363 : if (ixc_gga<0) then
196 : call rhotoxc(bigexc_corr,bigsxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
197 : & n3xccc_null,option,rhor,rprimd,usexcnhat,vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga,&
198 219 : & xc_funcs=xc_funcs_gga,strsxc=strsxc_corr)
199 : else
200 : call rhotoxc(bigexc_corr,bigsxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
201 : & n3xccc_null,option,rhor,rprimd,usexcnhat,vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga,&
202 144 : & strsxc=strsxc_corr)
203 : end if
204 363 : bigexc=bigexc-bigexc_corr
205 363 : bigsxc=bigsxc-bigsxc_corr
206 3848990 : vxc(:,:)=vxc(:,:)-vxc_corr(:,:)
207 363 : vxcavg=vxcavg-vxcavg_corr
208 2541 : strsxc(:)=strsxc(:)-strsxc_corr(:)
209 :
210 : !Release memory
211 726 : ABI_FREE(xccc3d_null)
212 : end if
213 :
214 493 : if (calcgrxc) then
215 :
216 100 : call metric(gmet,gprimd,-1,rmet,rprimd,ucvol)
217 : ! Initialize GGA functional
218 100 : if (ixc_gga<0 .and. libxc_gga_initialized==0) then
219 46 : call libxc_functionals_init(ixc_gga,dtset%nspden,xc_functionals=xc_funcs_gga)
220 46 : libxc_gga_initialized=1
221 : end if
222 300 : ABI_MALLOC(dyfrx2_dum,(3,3,dtset%natom))
223 300 : ABI_MALLOC(xccc3d_null,(n3xccc))
224 100 : option=1
225 : ! calculate xccc3d in this case
226 : call mkcore(strsxc_corr,dyfrx2_dum,grxc,mpi_enreg,dtset%natom,nfft,dtset%nspden,dtset%ntypat,ngfft(1),n1xccc,ngfft(2),&
227 100 : & ngfft(3),option,rprimd,dtset%typat,ucvol,vxc_corr,xcccrc,xccc1d,xccc3d_null,xred)
228 : !Add Vxc^GGA(rho_core+rho_val)
229 100 : if (ixc_gga<0) then
230 : call rhotoxc(bigexc_corr,bigsxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
231 : & n3xccc,option,rhor,rprimd,usexcnhat,vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga,&
232 46 : & strsxc=strsxc_corr,xc_funcs=xc_funcs_gga)
233 : else
234 : call rhotoxc(bigexc_corr,bigsxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
235 : & n3xccc,option,rhor,rprimd,usexcnhat,vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga,&
236 54 : & strsxc=strsxc_corr)
237 : end if
238 100 : option=2
239 : call mkcore(strsxc_corr,dyfrx2_dum,grxc,mpi_enreg,dtset%natom,nfft,dtset%nspden,dtset%ntypat,ngfft(1),n1xccc,ngfft(2),&
240 100 : & ngfft(3),option,rprimd,dtset%typat,ucvol,vxc_corr,xcccrc,xccc1d,xccc3d_null,xred)
241 100 : ABI_FREE(dyfrx2_dum)
242 100 : ABI_FREE(xccc3d_null)
243 : end if
244 :
245 493 : if(optstr_loc==1) then
246 : ! Initialize GGA functional
247 30 : if (ixc_gga<0 .and. libxc_gga_initialized==0) then
248 19 : call libxc_functionals_init(ixc_gga,dtset%nspden,xc_functionals=xc_funcs_gga)
249 19 : libxc_gga_initialized=1
250 : end if
251 : !calculate Vxc^GGA(rho_core+rho_val)
252 30 : option=0
253 30 : if (ixc_gga<0) then
254 : call rhotoxc(bigexc_corr,bigsxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
255 : & n3xccc,option,rhor,rprimd,usexcnhat,vxc,vxcavg_corr,xccc3d,xcdata_gga,&
256 19 : & strsxc=strsxc_corr,xc_funcs=xc_funcs_gga)
257 : else
258 : call rhotoxc(bigexc_corr,bigsxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
259 : & n3xccc,option,rhor,rprimd,usexcnhat,vxc,vxcavg_corr,xccc3d,xcdata_gga,&
260 11 : & strsxc=strsxc_corr)
261 : end if
262 : end if
263 :
264 : ! Suppress the temporary used xc functional
265 493 : if(libxc_gga_initialized==1) then
266 284 : call libxc_functionals_end(xc_functionals=xc_funcs_gga)
267 : end if
268 493 : ABI_FREE(vxc_corr)
269 493 : ABI_FREE(kxc_dum)
270 :
271 : DBG_EXIT("COLL")
272 :
273 986 : end subroutine xchybrid_ncpp_cc
274 : !!***
275 :
276 : end module m_xchybrid
277 : !!***
|