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-2018 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

PARENTS

CHILDREN

SOURCE

24 #if defined HAVE_CONFIG_H
25 #include "config.h"
26 #endif
27 
28 #include "abi_common.h"
29 
30 module m_xcdata
31 
32  use defs_basis
33  use m_errors
34  use libxc_functionals
35 
36  use defs_abitypes, only : dataset_type
37 
38  implicit none
39 
40  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

PARENTS

      calc_vhxc_me,invars2

CHILDREN

      get_xclevel

SOURCE

336 subroutine get_auxc_ixc(auxc_ixc,ixc)
337 
338 
339 !This section has been created automatically by the script Abilint (TD).
340 !Do not modify the following lines by hand.
341 #undef ABI_FUNC
342 #define ABI_FUNC 'get_auxc_ixc'
343 !End of the abilint section
344 
345  implicit none
346 
347 !Arguments ------------------------------------
348 !scalars
349  integer, intent(in) :: ixc
350  integer, intent(out) :: auxc_ixc
351 
352 !Local variables-------------------------------
353  integer :: usefock,xclevel
354 !integer :: gga_id(2)
355 
356 ! *************************************************************************
357 
358  auxc_ixc=11
359 !Native hybrid functionals from ABINIT
360  if (ixc==40.or.ixc==41.or.ixc==42) then
361    auxc_ixc = 11
362 !Hybrid functionals from libxc
363  else if (ixc<0) then
364    call get_xclevel(ixc,xclevel,usefock)
365    if(usefock==1)then
366      auxc_ixc=11
367 !    if (libxc_functionals_gga_from_hybrid(hybrid_id=ixc,gga_id=gga_id)) then
368 !      auxc_ixc=-gga_id(1)*1000-gga_id(2)
369 !    endif
370    end if
371  end if
372 
373 end subroutine get_auxc_ixc
374 
375 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

PARENTS

      invars2,m_xcdata,setup_sigma

CHILDREN

      get_xclevel

SOURCE

246 subroutine get_xclevel(ixc,xclevel,usefock)
247 
248 
249 !This section has been created automatically by the script Abilint (TD).
250 !Do not modify the following lines by hand.
251 #undef ABI_FUNC
252 #define ABI_FUNC 'get_xclevel'
253 !End of the abilint section
254 
255  implicit none
256 
257 !Arguments ------------------------------------
258 !scalars
259  integer, intent(in) :: ixc
260  integer, intent(out) :: xclevel
261  integer, intent(out), optional :: usefock
262 
263 !Local variables-------------------------------
264  integer :: ii,isiz,jj
265  character(len=500) :: message
266 
267 ! *************************************************************************
268 
269  xclevel=0
270  if( ( 1<=ixc .and. ixc<=10).or.(30<=ixc .and. ixc<=39).or.(ixc==50) )xclevel=1 ! LDA
271  if( (11<=ixc .and. ixc<=19).or.(23<=ixc .and. ixc<=29).or. ixc==1402000)xclevel=2 ! GGA
272  if( 20<=ixc .and. ixc<=22 )xclevel=3 ! ixc for TDDFT kernel tests
273  if(present(usefock))then
274    usefock=0
275    if( ixc>=40 .and. ixc<=42 )usefock=1 ! Hartree-Fock or internal hybrid functionals
276  endif
277  if( ixc>=41 .and. ixc<=42)xclevel=2 ! ixc for internal hybrids using GGA
278  if (ixc<0) then                                  ! libXC: metaGGA and hybrid functionals
279    xclevel=1
280    do isiz=1,2
281 !    ixc has ABINIT sign convention
282 !    ii has Libxc sign convention
283      if (isiz==1) ii=-ixc/1000
284      if (isiz==2) ii=-ixc-ii*1000
285      jj=libxc_functionals_family_from_id(ii)
286      if (jj==XC_FAMILY_GGA    .or.jj==XC_FAMILY_MGGA) xclevel=2
287      if (jj==XC_FAMILY_HYB_GGA.or.jj==XC_FAMILY_HYB_MGGA) then
288        xclevel=2
289        if(present(usefock))then
290          usefock=1
291        endif
292        if (.not.libxc_functionals_gga_from_hybrid(hybrid_id=ii)) then
293          write(message, '(a,i8,3a,i8,2a,2i8,2a)' )&
294 &         'ixc=',ixc,' (libXC hybrid functional) is presently not allowed.',ch10,&
295 &         'XC_FAMILY_HYB_GGA=',XC_FAMILY_HYB_GGA,ch10,&
296 &         'ii,jj=',ii,jj,ch10,&
297 &         'Action: try another hybrid functional.'
298          MSG_ERROR(message)
299        end if
300      end if
301    end do
302  end if
303 
304 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)]
  [usekden = 1 if the XC functional depends on the kinetic energy density]
  [vdw_xc = Choice of van-der-Waals density functional]
  [xc_tb09_c = Parameter for Tran-Blaha functional]

OUTPUT

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

SIDE EFFECTS

PARENTS

      calc_vhxc_me,energy,m_kxc,nonlinear,nres2vres,odamix,prcref,prcref_PMA
      respfn,rhotov,scfcv,setvtr,xchybrid_ncpp_cc

CHILDREN

      get_xclevel

SOURCE

150 subroutine xcdata_init(xcdata,auxc_ixc,dtset,hyb_mixing,intxc,ixc,nelect,nspden,tphysel,usekden,vdw_xc,xc_tb09_c,xc_denpos)
151 
152 
153 !This section has been created automatically by the script Abilint (TD).
154 !Do not modify the following lines by hand.
155 #undef ABI_FUNC
156 #define ABI_FUNC 'xcdata_init'
157 !End of the abilint section
158 
159  implicit none
160 
161 !Arguments ------------------------------------
162 !scalars
163  integer, intent(in),optional :: auxc_ixc,intxc,ixc,nspden,usekden,vdw_xc
164  real(dp),intent(in),optional :: hyb_mixing,nelect,tphysel,xc_denpos,xc_tb09_c
165  type(dataset_type), intent(in),optional :: dtset
166  type(xcdata_type), intent(out) :: xcdata
167 !Local variables-------------------------------
168  integer :: usefock,xclevel
169  character(len=500) :: message
170 
171 ! *************************************************************************
172 
173  if(present(dtset))then
174    xcdata%auxc_ixc=dtset%auxc_ixc
175    xcdata%intxc=dtset%intxc
176    xcdata%ixc=dtset%ixc
177    xcdata%nspden=dtset%nspden
178    xcdata%usekden=dtset%usekden
179    xcdata%vdw_xc=dtset%vdw_xc
180 
181    xcdata%hyb_mixing=abs(dtset%hyb_mixing) ! Warning : the absolute value is needed, because of the singular way
182                                            ! to define the default values for this input variable.
183    xcdata%nelect=dtset%nelect
184    xcdata%tphysel=dtset%tphysel
185    xcdata%xc_denpos=dtset%xc_denpos
186    xcdata%xc_tb09_c=dtset%xc_tb09_c
187 
188  else
189    if(.not.(present(auxc_ixc).and.present(intxc).and.present(ixc).and.&
190 &           present(usekden).and.present(vdw_xc).and.present(hyb_mixing).and.&
191 &           present(nelect).and.present(nspden).and.&
192 &           present(tphysel).and.present(xc_denpos).and.&
193 &           present(xc_tb09_c)))then
194      write(message,'(a)') &
195 &     ' If dtset is not provided, all the other optional arguments must be provided, which is not the case.'
196      MSG_BUG(message)
197    endif
198  endif
199 
200  if(present(auxc_ixc))  xcdata%auxc_ixc=auxc_ixc
201  if(present(intxc))     xcdata%intxc=intxc
202  if(present(ixc))       xcdata%ixc=ixc
203  if(present(nspden))    xcdata%nspden=nspden
204  if(present(usekden))   xcdata%usekden=usekden
205  if(present(vdw_xc))    xcdata%vdw_xc=vdw_xc
206 
207  if(present(hyb_mixing))xcdata%hyb_mixing=hyb_mixing
208  if(present(nelect))    xcdata%nelect=nelect
209  if(present(tphysel))   xcdata%tphysel=tphysel
210  if(present(xc_denpos)) xcdata%xc_denpos=xc_denpos
211  if(present(xc_tb09_c))  xcdata%xc_tb09_c=xc_tb09_c
212 
213 !Compute xclevel
214  call get_xclevel(xcdata%ixc,xclevel,usefock=usefock)
215  xcdata%xclevel=xclevel
216  xcdata%usefock=usefock
217 
218 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

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