File: dpotxc.f90

package info (click to toggle)
elkcode 5.4.24-2
  • links: PTS, VCS
  • area: main
  • in suites: buster, sid
  • size: 12,840 kB
  • sloc: f90: 48,415; fortran: 22,457; perl: 965; makefile: 384; sh: 369; python: 105; ansic: 67
file content (136 lines) | stat: -rw-r--r-- 4,683 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136

! Copyright (C) 2012 J. K. Dewhurst, S. Sharma and E. K. U. Gross.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

subroutine dpotxc
use modmain
use modphonon
implicit none
! local variables
integer idm,is,ias
integer nr,nri,nrc,nrci
integer ir,np
! allocatable arrays
real(8), allocatable :: fxcmt(:,:,:,:),fxcir(:,:,:)
complex(8), allocatable :: dvmt(:),dbmt(:,:),zfmt(:)
! compute the exchange-correlation kernel
if (spinpol) then
  allocate(fxcmt(npmtmax,natmtot,4,4),fxcir(ngtot,4,4))
  call genspfxcr(.false.,fxcmt,fxcir)
else
  allocate(fxcmt(npmtmax,natmtot,1,1),fxcir(ngtot,1,1))
  call genfxcr(.false.,fxcmt,fxcir)
end if
allocate(dvmt(npmtmax))
if (spinpol) allocate(dbmt(npmtmax,3),zfmt(npcmtmax))
!---------------------------------------!
!     muffin-tin potential and field    !
!---------------------------------------!
! note: muffin-tin functions are in spherical coordinates
do ias=1,natmtot
  is=idxis(ias)
  nr=nrmt(is)
  nri=nrmti(is)
  nrc=nrcmt(is)
  nrci=nrcmti(is)
  np=npmt(is)
! charge-charge contribution to potential derivative
  dvmt(1:np)=fxcmt(1:np,ias,1,1)*drhomt(1:np,ias)
! spin-polarised case
  if (spinpol) then
    if (ncmag) then
! non-collinear
! add charge-spin contribution to potential derivative
      dvmt(1:np)=dvmt(1:np) &
       +fxcmt(1:np,ias,1,2)*dmagmt(1:np,ias,1) &
       +fxcmt(1:np,ias,1,3)*dmagmt(1:np,ias,2) &
       +fxcmt(1:np,ias,1,4)*dmagmt(1:np,ias,3)
! spin-charge contribution to B-field derivative
      dbmt(1:np,1)=fxcmt(1:np,ias,1,2)*drhomt(1:np,ias)
      dbmt(1:np,2)=fxcmt(1:np,ias,1,3)*drhomt(1:np,ias)
      dbmt(1:np,3)=fxcmt(1:np,ias,1,4)*drhomt(1:np,ias)
! add spin-spin contribution to B-field derivative
! (note: fxc is stored as an upper triangular matrix)
      dbmt(1:np,1)=dbmt(1:np,1) &
       +fxcmt(1:np,ias,2,2)*dmagmt(1:np,ias,1) &
       +fxcmt(1:np,ias,2,3)*dmagmt(1:np,ias,2) &
       +fxcmt(1:np,ias,2,4)*dmagmt(1:np,ias,3)
      dbmt(1:np,2)=dbmt(1:np,2) &
       +fxcmt(1:np,ias,2,3)*dmagmt(1:np,ias,1) &
       +fxcmt(1:np,ias,3,3)*dmagmt(1:np,ias,2) &
       +fxcmt(1:np,ias,3,4)*dmagmt(1:np,ias,3)
      dbmt(1:np,3)=dbmt(1:np,3) &
       +fxcmt(1:np,ias,2,4)*dmagmt(1:np,ias,1) &
       +fxcmt(1:np,ias,3,4)*dmagmt(1:np,ias,2) &
       +fxcmt(1:np,ias,4,4)*dmagmt(1:np,ias,3)
    else
! collinear
! add charge-spin contribution to potential derivative
      dvmt(1:np)=dvmt(1:np)+fxcmt(1:np,ias,1,4)*dmagmt(1:np,ias,1)
! spin-charge contribution to B-field derivative
      dbmt(1:np,1)=fxcmt(1:np,ias,1,4)*drhomt(1:np,ias)
! add spin-spin contribution to B-field derivative
      dbmt(1:np,1)=dbmt(1:np,1)+fxcmt(1:np,ias,4,4)*dmagmt(1:np,ias,1)
    end if
  end if
! convert potential derivative to spherical harmonics
  call zfsht(nr,nri,dvmt,dvsmt(:,ias))
! convert magnetic field derivative to spherical harmonics on coarse mesh
  do idm=1,ndmag
    call zfmtftoc(nr,nri,dbmt(:,idm),zfmt)
    call zfsht(nrc,nrci,zfmt,dbsmt(:,ias,idm))
  end do
end do
!------------------------------------------!
!     interstitial potential and field     !
!------------------------------------------!
! charge-charge contribution to potential derivative
do ir=1,ngtot
  dvsir(ir)=fxcir(ir,1,1)*drhoir(ir)
end do
! spin-polarised case
if (spinpol) then
  if (ncmag) then
! non-collinear
    do ir=1,ngtot
! add charge-spin contribution to potential derivative
      dvsir(ir)=dvsir(ir) &
       +fxcir(ir,1,2)*dmagir(ir,1) &
       +fxcir(ir,1,3)*dmagir(ir,2) &
       +fxcir(ir,1,4)*dmagir(ir,3)
! spin-charge contribution to B-field derivative
      dbsir(ir,1)=fxcir(ir,1,2)*drhoir(ir)
      dbsir(ir,2)=fxcir(ir,1,3)*drhoir(ir)
      dbsir(ir,3)=fxcir(ir,1,4)*drhoir(ir)
! add spin-spin contribution to B-field derivative
      dbsir(ir,1)=dbsir(ir,1) &
       +fxcir(ir,2,2)*dmagir(ir,1) &
       +fxcir(ir,2,3)*dmagir(ir,2) &
       +fxcir(ir,2,4)*dmagir(ir,3)
      dbsir(ir,2)=dbsir(ir,2) &
       +fxcir(ir,2,3)*dmagir(ir,1) &
       +fxcir(ir,3,3)*dmagir(ir,2) &
       +fxcir(ir,3,4)*dmagir(ir,3)
      dbsir(ir,3)=dbsir(ir,3) &
       +fxcir(ir,2,4)*dmagir(ir,1) &
       +fxcir(ir,3,4)*dmagir(ir,2) &
       +fxcir(ir,4,4)*dmagir(ir,3)
    end do
  else
! collinear
    do ir=1,ngtot
! add charge-spin contribution to potential derivative
      dvsir(ir)=dvsir(ir)+fxcir(ir,1,4)*dmagir(ir,1)
! spin-charge contribution to B-field derivative
      dbsir(ir,1)=fxcir(ir,1,4)*drhoir(ir)
! add spin-spin contribution to B-field derivative
      dbsir(ir,1)=dbsir(ir,1)+fxcir(ir,4,4)*dmagir(ir,1)
    end do
  end if
end if
deallocate(fxcmt,fxcir,dvmt)
if (spinpol) deallocate(dbmt,zfmt)
return
end subroutine