TABLE OF CONTENTS
ABINIT/size_dvxc [ Functions ]
NAME
size_dvxc
FUNCTION
Give the size of the array dvxc(npts,ndvxc) and the second dimension of the d2vxc(npts,nd2vxc) needed for the allocations depending on the routine which is called from the drivexc routine
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (TD) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt. This routine has been written from rhotoxc (DCA, XG, GMR, MF, GZ)
INPUTS
[add_tfw]= optional flag controling the addition of Weiszacker gradient correction to Thomas-Fermi XC energy ixc= choice of exchange-correlation scheme order=gives the maximal derivative of Exc computed. 1=usual value (return exc and vxc) 2=also computes the kernel (return exc,vxc,kxc) -2=like 2, except (to be described) 3=also computes the derivative of the kernel (return exc,vxc,kxc,k3xc) [xc_funcs(2)]= <type(libxc_functional_type)>
OUTPUT
ndvxc size of the array dvxc(npts,ndvxc) for allocation ngr2 size of the array grho2_updn(npts,ngr2) for allocation nd2vxc size of the array d2vxc(npts,nd2vxc) for allocation nvxcdgr size of the array dvxcdgr(npts,nvxcdgr) for allocation
PARENTS
m_pawxc,rhotoxc
CHILDREN
SOURCE
40 #if defined HAVE_CONFIG_H 41 #include "config.h" 42 #endif 43 44 #include "abi_common.h" 45 46 47 subroutine size_dvxc(ixc,ndvxc,ngr2,nd2vxc,nspden,nvxcdgr,order,& 48 & add_tfw,xc_funcs) ! Optional 49 50 use defs_basis 51 use m_profiling_abi 52 use libxc_functionals 53 54 !This section has been created automatically by the script Abilint (TD). 55 !Do not modify the following lines by hand. 56 #undef ABI_FUNC 57 #define ABI_FUNC 'size_dvxc' 58 !End of the abilint section 59 60 implicit none 61 62 !Arguments---------------------- 63 integer, intent(in) :: ixc,nspden,order 64 integer, intent(out) :: ndvxc,nd2vxc,ngr2,nvxcdgr 65 logical, intent(in), optional :: add_tfw 66 type(libxc_functional_type),intent(in),optional :: xc_funcs(2) 67 68 !Local variables---------------- 69 logical :: add_tfw_,isgga,ismgga,is_hybrid 70 71 ! ************************************************************************* 72 73 add_tfw_=.false.;if (present(add_tfw)) add_tfw_=add_tfw 74 isgga=.false. ; ismgga=.false. ; is_hybrid=.false. 75 if(ixc<0)then 76 if(present(xc_funcs))then 77 isgga=libxc_functionals_isgga(xc_functionals=xc_funcs) 78 ismgga=libxc_functionals_ismgga(xc_functionals=xc_funcs) 79 is_hybrid=libxc_functionals_is_hybrid(xc_functionals=xc_funcs) 80 else 81 isgga=libxc_functionals_isgga() 82 ismgga=libxc_functionals_ismgga() 83 is_hybrid=libxc_functionals_is_hybrid() 84 end if 85 end if 86 87 ngr2=0 88 nvxcdgr=0 89 ndvxc=0 90 nd2vxc=0 91 92 !Dimension for the gradient of the density (only allocated for GGA or mGGA) 93 if ((ixc>=11.and.ixc<=17).or.(ixc>=23.and.ixc<=24).or.ixc==26.or.ixc==27.or. & 94 & (ixc>=31.and.ixc<=34).or.(ixc==41.or.ixc==42).or.ixc==1402000.or.(add_tfw_)) ngr2=2*min(nspden,2)-1 95 if (ixc<0.and.isgga.or.ismgga.or.is_hybrid) ngr2=2*min(nspden,2)-1 96 97 !A-Only Exc and Vxc 98 !======================================================================================= 99 if (order**2 <= 1) then 100 if (((ixc>=11 .and. ixc<=15) .and. ixc/=13) .or. (ixc>=23 .and. ixc<=24) .or. & 101 & (ixc==41 .or. ixc==42) .or. ixc==1402000) nvxcdgr=3 102 if (ixc==16.or.ixc==17.or.ixc==26.or.ixc==27) nvxcdgr=2 103 if (ixc<0) nvxcdgr=3 104 if (ixc>=31 .and. ixc<=34) nvxcdgr=3 !Native fake metaGGA functionals (for testing purpose only) 105 if (add_tfw_) nvxcdgr=3 106 107 ! B- Exc+Vxc and other derivatives 108 ! ======================================================================================= 109 else 110 111 ! Definition of ndvxc and nvxcdgr, 2nd dimension of the arrays of 2nd-order derivatives 112 ! ------------------------------------------------------------------------------------- 113 if (ixc==1 .or. ixc==21 .or. ixc==22 .or. (ixc>=7 .and. ixc<=10) .or. ixc==13) then 114 ! Routine xcspol: new Teter fit (4/93) to Ceperley-Alder data, with spin-pol option routine xcspol 115 ! Routine xcpbe, with different options (optpbe) and orders (order) 116 ndvxc=min(nspden,2)+1 117 else if (ixc>=2 .and. ixc<=6) then 118 ! Perdew-Zunger fit to Ceperly-Alder data (no spin-pol) !routine xcpzca 119 ! Teter fit (4/91) to Ceperley-Alder values (no spin-pol) !routine xctetr 120 ! Wigner xc (no spin-pol) !routine xcwign 121 ! Hedin-Lundqvist xc (no spin-pol) !routine xchelu 122 ! X-alpha (no spin-pol) !routine xcxalp 123 ndvxc=1 124 else if (ixc==12 .or. ixc==24) then 125 ! Routine xcpbe, with optpbe=-2 and different orders (order) 126 ndvxc=8 127 nvxcdgr=3 128 else if ((ixc>=11 .and. ixc<=15 .and. ixc/=13) .or. ixc==23 .or. ixc==41 .or. ixc==42) then 129 ! Routine xcpbe, with different options (optpbe) and orders (order) 130 ndvxc=15 131 nvxcdgr=3 132 else if(ixc==16 .or. ixc==17 .or. ixc==26 .or. ixc==27 ) then 133 nvxcdgr=2 134 else if (ixc==50) then 135 ndvxc=1 ! IIT xc (no spin-pol) 136 else if (ixc==1402000) then 137 ndvxc=15 138 nvxcdgr=3 139 else if (ixc<0) then 140 if(libxc_functionals_isgga() .or. libxc_functionals_ismgga() .or. libxc_functionals_is_hybrid()) then 141 ndvxc=15 142 else 143 ndvxc=3 144 end if 145 nvxcdgr=3 146 end if 147 if (add_tfw_) nvxcdgr=3 148 149 ! Definition of nd2vxc, 2nd dimension of the array of 3rd-order derivatives 150 ! ------------------------------------------------------------------------------------- 151 if (order==3) then 152 if (ixc==3) nd2vxc=1 ! Non spin polarized LDA case 153 if ((ixc>=7 .and. ixc<=10) .or. (ixc==13)) nd2vxc=3*min(nspden,2)-2 154 ! Following line to be corrected when the calculation of d2vxcar is implemented for these functionals 155 if ((ixc>=11 .and. ixc<=15 .and. ixc/=13) .or. (ixc==23.and.ixc<=24) .or. (ixc==41.or.ixc==42)) nd2vxc=1 156 if(ixc==1402000)nd2vxc=3*min(nspden,2)-2 157 if ((ixc<0.and.(.not.(libxc_functionals_isgga().or.libxc_functionals_ismgga().or.libxc_functionals_is_hybrid() )))) & 158 & nd2vxc=3*min(nspden,2)-2 159 end if 160 161 end if 162 163 end subroutine size_dvxc