TABLE OF CONTENTS


ABINIT/m_xcdata [ Modules ]

[ Top ] [ Modules ]

NAME

  m_xcdata

FUNCTION

  This module provides the definition of
  the xcdata_type used to drive the computation of the XC energy, potential, kernel, etc.

COPYRIGHT

  Copyright (C) 2017-2022 ABINIT group (XG)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

NOTES

SOURCE

19 #if defined HAVE_CONFIG_H
20 #include "config.h"
21 #endif
22 
23 #include "abi_common.h"
24 
25 module m_xcdata
26 
27  use defs_basis
28  use m_errors
29  use libxc_functionals
30  use m_dtset, only : dataset_type
31  use m_drivexc, only : size_dvxc
32 
33  implicit none
34 
35  private

m_xcdata/get_auxc_ixc [ Functions ]

[ Top ] [ m_xcdata ] [ Functions ]

NAME

  get_auxc_ixc

FUNCTION

  Returns the ixc of an auxiliary XC functional to be used instead of the input ixc
  For most of the functionals, there is no need of an auxiliary functional, in which case auxc_ixc=0
  For hybrid functionals, on the contrary, some speedup can be achieved by using such an auxiliary functional
  Note that this XC functional intend to replace the whole ixc functional. Generally speakin, it should be
  mistaken for the GGA part of the hybrid functional (that is for exchange only, actually).

  At present, always return ixc=1, but this might change in the future ...

INPUTS

  ixc= index of exchange-correlation functional

OUTPUT

  auxc_ixc= 0 if no need of an auxiliary functional, otherwise, returns the ixc of an auxiliary functional.

SIDE EFFECTS

SOURCE

298 subroutine get_auxc_ixc(auxc_ixc,ixc)
299 
300 !Arguments ------------------------------------
301 !scalars
302  integer, intent(in) :: ixc
303  integer, intent(out) :: auxc_ixc
304 
305 !Local variables-------------------------------
306  integer :: usefock,xclevel
307 !integer :: gga_id(2)
308 
309 ! *************************************************************************
310 
311  auxc_ixc=11
312 !Native hybrid functionals from ABINIT
313  if (ixc==40.or.ixc==41.or.ixc==42) then
314    auxc_ixc = 11
315 !Hybrid functionals from libxc
316  else if (ixc<0) then
317    call get_xclevel(ixc,xclevel,usefock)
318    if(usefock==1)then
319      auxc_ixc=11
320 !    if (libxc_functionals_gga_from_hybrid(hybrid_id=ixc,gga_id=gga_id)) then
321 !      auxc_ixc=-gga_id(1)*1000-gga_id(2)
322 !    endif
323    end if
324  end if
325 
326 end subroutine get_auxc_ixc
327 
328 end module m_xcdata

m_xcdata/get_xclevel [ Functions ]

[ Top ] [ m_xcdata ] [ Functions ]

NAME

  get_xclevel

FUNCTION

  Compute xclevel.

INPUTS

  ixc= index of exchange-correlation functional

OUTPUT

  [usefock = 1 if the XC functional needs the Fock operator]
  xclevel= 0 if no XC functional except possibly Fock; 1 if LDA; 2 if GGA ; 3 for TDDFT kernel tests

SIDE EFFECTS

SOURCE

221 subroutine get_xclevel(ixc,xclevel,usefock)
222 
223 !Arguments ------------------------------------
224 !scalars
225  integer, intent(in) :: ixc
226  integer, intent(out) :: xclevel
227  integer, intent(out), optional :: usefock
228 
229 !Local variables-------------------------------
230  integer :: ii,isiz,jj
231  character(len=500) :: message
232 
233 ! *************************************************************************
234 
235  xclevel=0 ; if(present(usefock)) usefock=0
236  if( ( 1<=ixc .and. ixc<=10).or.(30<=ixc .and. ixc<=39).or.(ixc==50) )xclevel=1 ! LDA
237  if( (11<=ixc .and. ixc<=19).or.(23<=ixc .and. ixc<=29).or. ixc==1402000)xclevel=2 ! GGA
238  if( 20<=ixc .and. ixc<=22 )xclevel=3 ! ixc for TDDFT kernel tests
239  if(present(usefock))then
240    if( ixc>=40 .and. ixc<=42 )usefock=1 ! Hartree-Fock or internal hybrid functionals
241  endif
242  if( ixc>=31 .and. ixc<=35)xclevel=2 ! ixc for internal fake mGGA
243  if( ixc>=41 .and. ixc<=42)xclevel=2 ! ixc for internal hybrids using GGA
244  if (ixc<0) then                     ! libXC: metaGGA and hybrid functionals
245    xclevel=1
246    do isiz=1,2
247 !    ixc has ABINIT sign convention
248 !    ii has Libxc sign convention
249      if (isiz==1) ii=-ixc/1000
250      if (isiz==2) ii=-ixc-ii*1000
251      if (ii<=0) cycle
252      jj=libxc_functionals_family_from_id(ii)
253      if (jj==XC_FAMILY_GGA    .or.jj==XC_FAMILY_MGGA) xclevel=2
254      if (jj==XC_FAMILY_HYB_GGA.or.jj==XC_FAMILY_HYB_MGGA) xclevel=2
255      if (present(usefock)) then
256        if (libxc_functionals_is_hybrid_from_id(ii)) usefock=1
257        if (usefock==1) then
258          if (.not.libxc_functionals_gga_from_hybrid(hybrid_id=ii)) then
259            write(message, '(a,i8,3a,2i8,2a)' )&
260 &           'ixc=',ixc,' (libXC hybrid functional) is presently not allowed.',ch10,&
261 &           'ii,jj=',ii,jj,ch10,&
262 &           'Action: try another hybrid functional.'
263            ABI_ERROR(message)
264          end if
265        end if
266      end if
267    end do
268  end if
269 
270 end subroutine get_xclevel

m_xcdata/xcdata_init [ Functions ]

[ Top ] [ m_xcdata ] [ Functions ]

NAME

  xcdata_init

FUNCTION

  Init the structure. Mostly copy input variables, except compute and usefock and xclevel.

INPUTS

  [dtset = the dataset from which the other input variables are taken, if they are not present]
  [auxc_ixc = possibly the index of the auxiliary xc functional, otherwise 0.]
  [hyb_mixing = parameter for mixing Fock exchange in native PBEx functionals]
  [intxc = 1 if the XC functional has to be interpolated on a more refined mesh than the FFT one]
  [ixc= index of exchange-correlation functional]
  [nelect = Number of electrons in the cell (for Fermi-Amaldi only)]
  [tphysel = Physical temperature (for temperature-dependent functional)]
  [vdw_xc = Choice of van-der-Waals density functional]

OUTPUT

  xcdata <type(xcdata_type)>= the data to calculate exchange-correlation are initialized

SIDE EFFECTS

SOURCE

140 subroutine xcdata_init(xcdata,auxc_ixc,dtset,hyb_mixing,intxc,ixc,nelect,nspden,tphysel,&
141 &                      vdw_xc,xc_denpos)
142 
143 !Arguments ------------------------------------
144 !scalars
145  integer, intent(in),optional :: auxc_ixc,intxc,ixc,nspden,vdw_xc
146  real(dp),intent(in),optional :: hyb_mixing,nelect,tphysel,xc_denpos
147  type(dataset_type), intent(in),optional :: dtset
148  type(xcdata_type), intent(out) :: xcdata
149 !Local variables-------------------------------
150  integer :: nspden_updn
151  character(len=500) :: msg
152 
153 ! *************************************************************************
154 
155  if(present(dtset))then
156    xcdata%auxc_ixc=dtset%auxc_ixc
157    xcdata%intxc=dtset%intxc
158    xcdata%ixc=dtset%ixc
159    xcdata%nspden=dtset%nspden
160    xcdata%vdw_xc=dtset%vdw_xc
161 
162    xcdata%hyb_mixing=abs(dtset%hyb_mixing) ! Warning : the absolute value is needed, because of the singular way
163                                            ! to define the default values for this input variable.
164    xcdata%nelect=dtset%nelect
165    xcdata%tphysel=merge(dtset%tphysel,dtset%tsmear,dtset%tphysel>tol8.and.dtset%occopt/=3.and.dtset%occopt/=9)
166 
167    xcdata%xc_denpos=dtset%xc_denpos
168 
169  else
170    if(.not.(present(auxc_ixc).and.present(intxc).and.present(ixc).and.&
171 &           present(vdw_xc).and.present(hyb_mixing).and.&
172 &           present(nelect).and.present(nspden).and.&
173 &           present(tphysel).and.present(xc_denpos)))then
174      msg='If dtset is not provided, all the other optional arguments must be provided, which is not the case!'
175      ABI_BUG(msg)
176    endif
177  endif
178 
179  if(present(auxc_ixc))  xcdata%auxc_ixc=auxc_ixc
180  if(present(intxc))     xcdata%intxc=intxc
181  if(present(ixc))       xcdata%ixc=ixc
182  if(present(nspden))    xcdata%nspden=nspden
183  if(present(vdw_xc))    xcdata%vdw_xc=vdw_xc
184 
185  if(present(hyb_mixing))xcdata%hyb_mixing=hyb_mixing
186  if(present(nelect))    xcdata%nelect=nelect
187  if(present(tphysel))   xcdata%tphysel=tphysel
188  if(present(xc_denpos)) xcdata%xc_denpos=xc_denpos
189 
190 !Compute xclevel
191  call get_xclevel(xcdata%ixc,xcdata%xclevel,usefock=xcdata%usefock)
192 
193 !Compute usegradient,uselaplacian,usekden
194  nspden_updn=min(xcdata%nspden,2)
195  call size_dvxc(xcdata%ixc,1,nspden_updn,usegradient=xcdata%usegradient,&
196 &               uselaplacian=xcdata%uselaplacian,usekden=xcdata%usekden)
197 
198 end subroutine xcdata_init

m_xcdata/xcdata_type [ Types ]

[ Top ] [ m_xcdata ] [ Types ]

NAME

  xcdata_type

FUNCTION

   This object stores the input variables (and derived parameters) needed to compute the exchange-correlation functional,
   not simply to define it.

NOTES

SOURCE

 51  type, public :: xcdata_type
 52 
 53 ! Integer scalars
 54 
 55   integer :: auxc_ixc
 56     ! Choice of auxiliary exchange-correlation functional. See input variable documentation
 57     ! If 0, there is no auxiliary xc functional, the one corresponding to ixc has to be used.
 58 
 59   integer :: intxc
 60     ! 1 if the XC functional has to be interpolated on a more refined mesh than the FFT one.
 61     ! 0 stick to the original FFT mesh
 62 
 63   integer :: ixc
 64     ! Choice of exchange-correlation functional. See input variable documentation
 65 
 66   integer :: nspden
 67     ! Number of spin components of the density
 68 
 69   integer :: usefock
 70     ! 1 if the XC functional includes a (possibly screened) Fock contribution
 71 
 72   integer :: usegradient
 73     ! 1 if the XC functional depends on the density gradient
 74 
 75   integer :: uselaplacian
 76     ! 1 if the XC functional depends on the density laplacian
 77 
 78   integer :: usekden
 79     ! 1 if the XC functional depends on the kinetic energy density
 80 
 81   integer :: vdw_xc
 82     ! Choice of van-der-Waals density functional. See input variable documentation
 83 
 84   integer :: xclevel
 85     ! Determined from ixc
 86     ! 0 if no XC functional
 87     ! 1 if LDA-type XC functional
 88     ! 2 if GGA-type XC functional
 89     ! 3 if for TDDFT kernel
 90 
 91 ! Real scalars
 92 
 93   real(dp) :: hyb_mixing
 94     ! Parameter for mixing Fock exchange in native PBEx functionals
 95 
 96   real(dp) :: nelect
 97     ! Number of electrons in the cell (for Fermi-Amaldi only)
 98 
 99   real(dp) :: tphysel
100     ! Physical temperature (for temperature-dependent functional)
101 
102   real(dp) :: xc_denpos
103     ! density positivity value
104 
105  end type xcdata_type
106 
107 !----------------------------------------------------------------------
108 
109  public :: xcdata_init                ! Initialize the object.
110  public :: get_xclevel                ! Get the xclevel from ixc (as well as usefock)
111  public :: get_auxc_ixc               ! Get the auxiliary xc functional (if it exists)
112 
113 contains