1
2! Copyright (C) 2011 J. K. Dewhurst, S. Sharma and E. K. U. Gross.
3! This file is distributed under the terms of the GNU General Public License.
4! See the file COPYING for license details.
5
6subroutine xc_c_tb09
7use modmain
8use libxcifc
9implicit none
10! local variables
11integer is,ias,i
12integer nr,nri,ir
13real(8), parameter :: alpha=-0.012d0, beta=1.023d0
14real(8) t1
15! allocatable arrays
16real(8), allocatable :: grfmt(:,:,:),grfir(:,:)
17real(8), allocatable :: rfmt(:,:),rfir(:)
18real(8), allocatable :: rfmt1(:),rfmt2(:,:)
19! external functions
20real(8), external :: rfint
21! check that the Tran-Blaha functional is being used (A. Shyichuk)
22if (xctype(2).ne.XC_MGGA_X_TB09) return
23! if Tran-Blaha constant has been read in return
24if (tc_tb09) return
25! compute the gradient of the density
26allocate(grfmt(npmtmax,natmtot,3),grfir(ngtot,3))
27call gradrf(rhomt,rhoir,grfmt,grfir)
28allocate(rfmt(npmtmax,natmtot),rfmt1(npmtmax),rfmt2(npmtmax,3))
29do ias=1,natmtot
30  is=idxis(ias)
31  nr=nrmt(is)
32  nri=nrmti(is)
33! convert muffin-tin density to spherical coordinates
34  call rbsht(nr,nri,rhomt(:,ias),rfmt1)
35! convert muffin-tin gradient to spherical coordinates
36  do i=1,3
37    call rbsht(nr,nri,grfmt(:,ias,i),rfmt2(:,i))
38  end do
39! integrand in muffin-tin
40  do i=1,npmt(is)
41    t1=sqrt(rfmt2(i,1)**2+rfmt2(i,2)**2+rfmt2(i,3)**2)
42    rfmt1(i)=t1/rfmt1(i)
43  end do
44! convert to spherical harmonics
45  call rfsht(nr,nri,rfmt1,rfmt(:,ias))
46end do
47deallocate(grfmt,rfmt1,rfmt2)
48! integrand in interstitial
49allocate(rfir(ngtot))
50do ir=1,ngtot
51  t1=sqrt(grfir(ir,1)**2+grfir(ir,2)**2+grfir(ir,3)**2)
52  rfir(ir)=t1/rhoir(ir)
53end do
54! integrate over the unit cell
55t1=rfint(rfmt,rfir)
56! set the constant
57c_tb09=alpha+beta*sqrt(abs(t1)/omega)
58deallocate(grfir,rfmt,rfir)
59end subroutine
60
61