TABLE OF CONTENTS


ABINIT/hybrid_corr [ Functions ]

[ Top ] [ 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.

PARENTS

CHILDREN

      dtset_copy,dtset_free,libxc_functionals_end,libxc_functionals_init
      rhohxc

SOURCE

338 subroutine hybrid_corr(dtset,ixc,nkxc,mpi_enreg,nfft,ngfft,nspden,rhor,rprimd,hybrid_mixing,vxc,enxc)
339 
340 
341 !This section has been created automatically by the script Abilint (TD).
342 !Do not modify the following lines by hand.
343 #undef ABI_FUNC
344 #define ABI_FUNC 'hybrid_corr'
345 !End of the abilint section
346 
347  implicit none
348 
349 !Arguments -------------------------------------------------------------
350 !scalars
351  integer,intent(in) :: ixc,nfft,nspden,nkxc
352  real(dp),intent(in) :: hybrid_mixing
353  real(dp),intent(inout) :: enxc
354  type(MPI_type),intent(in) :: mpi_enreg
355  type(dataset_type),intent(in) :: dtset
356 !arrays
357  integer,intent(in) :: ngfft(18)
358  real(dp),intent(in) :: rhor(nfft,nspden),rprimd(3,3)
359  real(dp),intent(inout) :: vxc(nfft,nspden)
360 
361 !Local variables -------------------------------------------------------
362 !No improved xc quadrature.
363 !No core correction.
364 !Dummy here.
365 !For debugging purposes (see tests below):
366 !integer :: i1,i2,i3,k1,n1,n2,n3
367 !real(dp) :: kx,rho,rhomax,ftest
368 !scalars
369  character(len=500) :: message
370  type(dataset_type) :: dtLocal
371 !arrays
372  real(dp) :: dum(0)
373  real(dp) :: enxc_corr
374  real(dp),allocatable :: kxcr(:,:)
375  real(dp),allocatable :: xccc3d(:),vxc_corr(:,:)
376 
377 !***********************************************************************
378 !For debugging purposes (see tests below):
379 
380 !ftest(i1,n1,k1) = 0._dp+1._dp*cos(k1*two_pi*float(i1)/float(n1))
381 
382 !***********************************************************************
383 
384 !Check input parameters.
385 
386  if (nspden > 2) then
387    message = ' kxc_alda does not work yet for nspden > 2.'
388    MSG_ERROR(message)
389  end if
390 
391 !Copy the input variables from the current dataset to a temporary one
392 !to tune some parameters
393  call dtset_copy(dtLocal, dtset)
394  dtLocal%intxc = 0
395  dtLocal%ixc   = -101
396 
397  ! Reinitialize the libxc module with the overriden values
398  if (dtset%ixc<0) then
399    call libxc_functionals_end()
400  end if
401  if (dtLocal%ixc<0) then
402    call libxc_functionals_init(dtLocal%ixc,dtLocal%nspden)
403  end if
404 
405  call rhohxc(dtLocal,enxc_corr,dum,0,kxcr,mpi_enreg,nfft,dum,dum,0,dum,0,nkxc,0,&
406 & dtset%usepawu==4.or.dtset%usepawu==14,nspden,n3xccc,&
407 & 0,dum,rhor,rprimd,strsxc,1,dum,vxc_corr,dum,xccc3d)
408 
409  vxc(:,:) = vxc(:,:) + hybrid_mixing*vxc_corr(:,:)
410  enxc = enxc + hybrid_mixing*enxc_corr
411 
412  call rhohxc(dtLocal,enxc_corr,dum,0,kxcr,mpi_enreg,dum,dum,dum,0,dum,0,nkxc,0,&
413 & dtset%usepawu==4.or.dtset%usepawu==14,nspden,0,&
414 & 0,dum,rhor,rprimd,strsxc,1,dum,vxc_corr,vxcavg,0)
415 
416  vxc(:,:) = vxc(:,:) - hybrid_mixing*vxc_corr(:,:)
417  enxc = enxc - hybrid_mixing*enxc_corr
418 
419 
420 ! Revert libxc module to the original settings
421  if (dtLocal%ixc<0) then
422    call libxc_functionals_end()
423  end if
424  if (dtset%ixc<0) then
425    call libxc_functionals_init(dtset%ixc,dtset%nspden)
426  end if
427 
428 !Free memory.
429  call dtset_free(dtLocal)
430  ABI_FREE(kxcr)
431  ABI_FREE(xccc3d)
432 
433 end subroutine hybrid_corr

ABINIT/m_xchybrid [ Modules ]

[ Top ] [ Modules ]

NAME

  m_xchybrid

FUNCTION

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 .

PARENTS

CHILDREN

SOURCE

20 #if defined HAVE_CONFIG_H
21 #include "config.h"
22 #endif
23 
24 #include "abi_common.h"
25 
26 module m_xchybrid
27 
28  use defs_basis
29  use m_abicore
30  use m_errors
31  use m_xcdata
32  use libxc_functionals
33 
34  use m_geometry,    only : metric
35  use defs_abitypes, only : MPI_type, dataset_type
36  use m_dtset,       only : dtset_copy, dtset_free
37  use m_rhotoxc,     only : rhotoxc
38  use m_mkcore,      only : mkcore
39 
40  implicit none
41 
42  private

ABINIT/xchybrid_ncpp_cc [ Functions ]

[ Top ] [ 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

PARENTS

      forces,forstr,rhotov,setvtr

CHILDREN

      libxc_functionals_end,libxc_functionals_init,metric,mkcore,rhotoxc
      xcdata_init

SOURCE

105 subroutine xchybrid_ncpp_cc(dtset,enxc,mpi_enreg,nfft,ngfft,n3xccc,rhor,rprimd,strsxc,vxcavg,xccc3d,vxc,grxc,xcccrc,xccc1d,&
106 &                           xred,n1xccc,optstr)
107 
108 
109 !This section has been created automatically by the script Abilint (TD).
110 !Do not modify the following lines by hand.
111 #undef ABI_FUNC
112 #define ABI_FUNC 'xchybrid_ncpp_cc'
113 !End of the abilint section
114 
115  implicit none
116 
117 !Arguments -------------------------------------------------------------
118 !scalars
119  integer,intent(in) :: nfft,n3xccc
120  integer,optional,intent(in) :: n1xccc,optstr
121  real(dp),intent(out) :: enxc,vxcavg
122  type(dataset_type),intent(in) :: dtset
123  type(MPI_type),intent(in) :: mpi_enreg
124 !arrays
125  integer,intent(in) :: ngfft(18)
126  real(dp),intent(in) :: rhor(nfft,dtset%nspden),rprimd(3,3),xccc3d(n3xccc)
127  real(dp),intent(out) :: strsxc(6)
128  real(dp),optional,intent(in) :: xcccrc(dtset%ntypat),xred(3,dtset%natom),xccc1d(:,:,:)
129  real(dp),optional,intent(out) :: grxc(3,dtset%natom),vxc(nfft,dtset%nspden)
130 
131 !Local variables -------------------------------------------------------
132 !scalars
133  integer :: ixc_gga,libxc_gga_initialized,ndim,nkxc,n3xccc_null,option,optstr_loc,usexcnhat
134  real(dp) :: enxc_corr,ucvol,vxcavg_corr
135  character(len=500) :: msg
136  type(xcdata_type) :: xcdata_gga,xcdata_hybrid
137  logical :: calcgrxc,nmxc
138 !arrays
139  integer :: gga_id(2)
140  real(dp) :: nhat(1,0),nhatgr(1,1,0),strsxc_corr(6),gmet(3,3),gprimd(3,3),rmet(3,3)
141  real(dp),allocatable :: kxc_dum(:,:),vxc_corr(:,:),xccc3d_null(:),dyfrx2_dum(:,:,:)
142  type(libxc_functional_type) :: xc_funcs_gga(2)
143 
144 ! *************************************************************************
145 
146  DBG_ENTER("COLL")
147 
148 
149  nmxc=(dtset%usepawu==4).or.(dtset%usepawu==14)
150 
151 !Not relevant for PAW
152  if (dtset%usepaw==1) return
153  if(n3xccc==0) return
154  calcgrxc=(present(grxc).and.present(n1xccc).and.present(xcccrc).and.present(xred).and.present(xccc1d))
155  optstr_loc=0
156  if(present(optstr)) optstr_loc=1
157 !Not applicable for electron-positron
158  if (dtset%positron>0) then
159    msg='NCPP+Hybrid functionals not applicable for electron-positron calculations!'
160    MSG_ERROR(msg)
161  end if
162 
163 !Select the GGA on which the hybrid functional is based on
164 !or return if not applicable
165  if (dtset%ixc==41.or.dtset%ixc==42) then
166    ixc_gga = 11
167  else if (dtset%ixc<0) then
168    if (libxc_functionals_gga_from_hybrid(gga_id=gga_id)) then
169      ixc_gga=-gga_id(2)*1000-gga_id(1)
170    else
171      return
172    end if
173  else
174    return
175  end if
176 
177 !Define xcdata_hybrid as well as xcdata_gga
178  call xcdata_init(xcdata_hybrid,dtset=dtset)
179  call xcdata_init(xcdata_gga,dtset=dtset,auxc_ixc=0,ixc=ixc_gga)
180  libxc_gga_initialized=0
181 
182  nkxc=0;ndim=0;usexcnhat=0;n3xccc_null=0
183  ABI_ALLOCATE(kxc_dum,(nfft,nkxc))
184  ABI_ALLOCATE(vxc_corr,(nfft,dtset%nspden))
185 
186  if (present(vxc).and.optstr_loc==0) then
187 !Initialize args for rhotoxc
188    option=0 ! XC only
189    ABI_ALLOCATE(xccc3d_null,(n3xccc_null))
190 !  Compute Vxc^Hybrid(rho_val)
191    call rhotoxc(enxc,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
192 &   n3xccc_null,option,dtset%paral_kgb,rhor,rprimd,&
193 &   strsxc,usexcnhat,vxc,vxcavg,xccc3d_null,xcdata_hybrid)
194 
195 !  Initialize GGA functional
196    if (ixc_gga<0) then
197      call libxc_functionals_init(ixc_gga,dtset%nspden,xc_functionals=xc_funcs_gga)
198      libxc_gga_initialized=1
199    end if
200 
201 !Add Vxc^GGA(rho_core+rho_val)
202    if (ixc_gga<0) then
203      call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
204 &     n3xccc,option,dtset%paral_kgb,rhor,rprimd,&
205 &     strsxc_corr,usexcnhat,vxc_corr,vxcavg_corr,xccc3d,xcdata_gga,xc_funcs=xc_funcs_gga)
206    else
207      call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
208 &     n3xccc,option,dtset%paral_kgb,rhor,rprimd,&
209 &     strsxc_corr,usexcnhat,vxc_corr,vxcavg_corr,xccc3d,xcdata_gga)
210    end if
211    enxc=enxc+enxc_corr
212    vxc(:,:)=vxc(:,:)+vxc_corr(:,:)
213    vxcavg=vxcavg+vxcavg_corr
214    strsxc(:)=strsxc(:)+strsxc_corr(:)
215 
216 !Substract Vxc^GGA(rho_val)
217    if (ixc_gga<0) then
218      call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
219 &     n3xccc_null,option,dtset%paral_kgb,rhor,rprimd,strsxc_corr,usexcnhat,&
220 &     vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga,xc_funcs=xc_funcs_gga)
221    else
222      call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
223 &     n3xccc_null,option,dtset%paral_kgb,rhor,rprimd,strsxc_corr,usexcnhat,&
224 &     vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga)
225    end if
226    enxc=enxc-enxc_corr
227    vxc(:,:)=vxc(:,:)-vxc_corr(:,:)
228    vxcavg=vxcavg-vxcavg_corr
229    strsxc(:)=strsxc(:)-strsxc_corr(:)
230 
231 !Release memory
232    ABI_DEALLOCATE(xccc3d_null)
233  end if
234 
235  if (calcgrxc) then
236 
237    call metric(gmet,gprimd,-1,rmet,rprimd,ucvol)
238 !  Initialize GGA functional
239    if (ixc_gga<0 .and. libxc_gga_initialized==0) then
240      call libxc_functionals_init(ixc_gga,dtset%nspden,xc_functionals=xc_funcs_gga)
241      libxc_gga_initialized=1
242    end if
243    ABI_ALLOCATE(dyfrx2_dum,(3,3,dtset%natom))
244    ABI_ALLOCATE(xccc3d_null,(n3xccc))
245    option=1
246 ! calculate xccc3d in this case
247    call mkcore(strsxc_corr,dyfrx2_dum,grxc,mpi_enreg,dtset%natom,nfft,dtset%nspden,dtset%ntypat,ngfft(1),n1xccc,ngfft(2),&
248 &   ngfft(3),option,rprimd,dtset%typat,ucvol,vxc_corr,xcccrc,xccc1d,xccc3d_null,xred)
249 !Add Vxc^GGA(rho_core+rho_val)
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,dtset%paral_kgb,&
253 &     rhor,rprimd,strsxc_corr,usexcnhat,vxc_corr,vxcavg_corr,xccc3d_null,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,dtset%paral_kgb,&
257 &     rhor,rprimd,strsxc_corr,usexcnhat,vxc_corr,vxcavg_corr,xccc3d_null,xcdata_gga)
258    end if
259    option=2
260    call mkcore(strsxc_corr,dyfrx2_dum,grxc,mpi_enreg,dtset%natom,nfft,dtset%nspden,dtset%ntypat,ngfft(1),n1xccc,ngfft(2),&
261 &   ngfft(3),option,rprimd,dtset%typat,ucvol,vxc_corr,xcccrc,xccc1d,xccc3d_null,xred)
262    ABI_DEALLOCATE(dyfrx2_dum)
263    ABI_DEALLOCATE(xccc3d_null)
264  end if
265 
266  if(optstr_loc==1) then
267 !  Initialize GGA functional
268    if (ixc_gga<0 .and. libxc_gga_initialized==0) then
269      call libxc_functionals_init(ixc_gga,dtset%nspden,xc_functionals=xc_funcs_gga)
270      libxc_gga_initialized=1
271    end if
272 !calculate Vxc^GGA(rho_core+rho_val)
273    option=0
274    if (ixc_gga<0) then
275      call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
276 &     n3xccc,option,dtset%paral_kgb,rhor,rprimd,&
277 &     strsxc_corr,usexcnhat,vxc,vxcavg_corr,xccc3d,xcdata_gga,xc_funcs=xc_funcs_gga)
278    else
279      call rhotoxc(enxc_corr,kxc_dum,mpi_enreg,nfft,ngfft,nhat,ndim,nhatgr,ndim,nkxc,nkxc,nmxc,&
280 &     n3xccc,option,dtset%paral_kgb,rhor,rprimd,&
281 &     strsxc_corr,usexcnhat,vxc,vxcavg_corr,xccc3d,xcdata_gga)
282    end if
283  end if
284 
285 ! Suppress the temporary used xc functional
286  if(libxc_gga_initialized==1) then
287    call libxc_functionals_end(xc_functionals=xc_funcs_gga)
288  end if
289  ABI_DEALLOCATE(vxc_corr)
290  ABI_DEALLOCATE(kxc_dum)
291 
292  DBG_EXIT("COLL")
293 
294 end subroutine xchybrid_ncpp_cc