1
2! Copyright (C) 2012 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 dpotxc
7use modmain
8use modphonon
9implicit none
10! local variables
11integer idm,is,ias
12integer nr,nri,nrc,nrci
13integer ir,np
14! allocatable arrays
15real(8), allocatable :: fxcmt(:,:,:,:),fxcir(:,:,:)
16complex(8), allocatable :: dvmt(:),dbmt(:,:),zfmt(:)
17! compute the exchange-correlation kernel
18if (spinpol) then
19  allocate(fxcmt(npmtmax,natmtot,4,4),fxcir(ngtot,4,4))
20  call genspfxcr(.false.,fxcmt,fxcir)
21else
22  allocate(fxcmt(npmtmax,natmtot,1,1),fxcir(ngtot,1,1))
23  call genfxcr(.false.,fxcmt,fxcir)
24end if
25allocate(dvmt(npmtmax))
26if (spinpol) allocate(dbmt(npmtmax,3),zfmt(npcmtmax))
27!---------------------------------------!
28!     muffin-tin potential and field    !
29!---------------------------------------!
30! note: muffin-tin functions are in spherical coordinates
31do ias=1,natmtot
32  is=idxis(ias)
33  nr=nrmt(is)
34  nri=nrmti(is)
35  nrc=nrcmt(is)
36  nrci=nrcmti(is)
37  np=npmt(is)
38! charge-charge contribution to potential derivative
39  dvmt(1:np)=fxcmt(1:np,ias,1,1)*drhomt(1:np,ias)
40! spin-polarised case
41  if (spinpol) then
42    if (ncmag) then
43! non-collinear
44! add charge-spin contribution to potential derivative
45      dvmt(1:np)=dvmt(1:np) &
46       +fxcmt(1:np,ias,1,2)*dmagmt(1:np,ias,1) &
47       +fxcmt(1:np,ias,1,3)*dmagmt(1:np,ias,2) &
48       +fxcmt(1:np,ias,1,4)*dmagmt(1:np,ias,3)
49! spin-charge contribution to B-field derivative
50      dbmt(1:np,1)=fxcmt(1:np,ias,1,2)*drhomt(1:np,ias)
51      dbmt(1:np,2)=fxcmt(1:np,ias,1,3)*drhomt(1:np,ias)
52      dbmt(1:np,3)=fxcmt(1:np,ias,1,4)*drhomt(1:np,ias)
53! add spin-spin contribution to B-field derivative
54! (note: fxc is stored as an upper triangular matrix)
55      dbmt(1:np,1)=dbmt(1:np,1) &
56       +fxcmt(1:np,ias,2,2)*dmagmt(1:np,ias,1) &
57       +fxcmt(1:np,ias,2,3)*dmagmt(1:np,ias,2) &
58       +fxcmt(1:np,ias,2,4)*dmagmt(1:np,ias,3)
59      dbmt(1:np,2)=dbmt(1:np,2) &
60       +fxcmt(1:np,ias,2,3)*dmagmt(1:np,ias,1) &
61       +fxcmt(1:np,ias,3,3)*dmagmt(1:np,ias,2) &
62       +fxcmt(1:np,ias,3,4)*dmagmt(1:np,ias,3)
63      dbmt(1:np,3)=dbmt(1:np,3) &
64       +fxcmt(1:np,ias,2,4)*dmagmt(1:np,ias,1) &
65       +fxcmt(1:np,ias,3,4)*dmagmt(1:np,ias,2) &
66       +fxcmt(1:np,ias,4,4)*dmagmt(1:np,ias,3)
67    else
68! collinear
69! add charge-spin contribution to potential derivative
70      dvmt(1:np)=dvmt(1:np)+fxcmt(1:np,ias,1,4)*dmagmt(1:np,ias,1)
71! spin-charge contribution to B-field derivative
72      dbmt(1:np,1)=fxcmt(1:np,ias,1,4)*drhomt(1:np,ias)
73! add spin-spin contribution to B-field derivative
74      dbmt(1:np,1)=dbmt(1:np,1)+fxcmt(1:np,ias,4,4)*dmagmt(1:np,ias,1)
75    end if
76  end if
77! convert potential derivative to spherical harmonics
78  call zfsht(nr,nri,dvmt,dvsmt(:,ias))
79! convert magnetic field derivative to spherical harmonics on coarse mesh
80  do idm=1,ndmag
81    call zfmtftoc(nrc,nrci,dbmt(:,idm),zfmt)
82    call zfsht(nrc,nrci,zfmt,dbsmt(:,ias,idm))
83  end do
84end do
85!------------------------------------------!
86!     interstitial potential and field     !
87!------------------------------------------!
88! charge-charge contribution to potential derivative
89do ir=1,ngtot
90  dvsir(ir)=fxcir(ir,1,1)*drhoir(ir)
91end do
92! spin-polarised case
93if (spinpol) then
94  if (ncmag) then
95! non-collinear
96    do ir=1,ngtot
97! add charge-spin contribution to potential derivative
98      dvsir(ir)=dvsir(ir) &
99       +fxcir(ir,1,2)*dmagir(ir,1) &
100       +fxcir(ir,1,3)*dmagir(ir,2) &
101       +fxcir(ir,1,4)*dmagir(ir,3)
102! spin-charge contribution to B-field derivative
103      dbsir(ir,1)=fxcir(ir,1,2)*drhoir(ir)
104      dbsir(ir,2)=fxcir(ir,1,3)*drhoir(ir)
105      dbsir(ir,3)=fxcir(ir,1,4)*drhoir(ir)
106! add spin-spin contribution to B-field derivative
107      dbsir(ir,1)=dbsir(ir,1) &
108       +fxcir(ir,2,2)*dmagir(ir,1) &
109       +fxcir(ir,2,3)*dmagir(ir,2) &
110       +fxcir(ir,2,4)*dmagir(ir,3)
111      dbsir(ir,2)=dbsir(ir,2) &
112       +fxcir(ir,2,3)*dmagir(ir,1) &
113       +fxcir(ir,3,3)*dmagir(ir,2) &
114       +fxcir(ir,3,4)*dmagir(ir,3)
115      dbsir(ir,3)=dbsir(ir,3) &
116       +fxcir(ir,2,4)*dmagir(ir,1) &
117       +fxcir(ir,3,4)*dmagir(ir,2) &
118       +fxcir(ir,4,4)*dmagir(ir,3)
119    end do
120  else
121! collinear
122    do ir=1,ngtot
123! add charge-spin contribution to potential derivative
124      dvsir(ir)=dvsir(ir)+fxcir(ir,1,4)*dmagir(ir,1)
125! spin-charge contribution to B-field derivative
126      dbsir(ir,1)=fxcir(ir,1,4)*drhoir(ir)
127! add spin-spin contribution to B-field derivative
128      dbsir(ir,1)=dbsir(ir,1)+fxcir(ir,4,4)*dmagir(ir,1)
129    end do
130  end if
131end if
132deallocate(fxcmt,fxcir,dvmt)
133if (spinpol) deallocate(dbmt,zfmt)
134end subroutine
135
136