File: zftzf.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 (134 lines) | stat: -rw-r--r-- 3,298 bytes parent folder | download
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

! 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 zftzf(ngp,jlgpr,ylmgp,ld,sfacgp,zfmt,zfir,zfgp)
use modmain
implicit none
! arguments
integer, intent(in) :: ngp
real(8), intent(in) :: jlgpr(njcmax,nspecies,ngp)
complex(8), intent(in) :: ylmgp(lmmaxo,ngp)
integer, intent(in) :: ld
complex(8), intent(in) :: sfacgp(ld,natmtot)
complex(8), intent(in) :: zfmt(npcmtmax,natmtot),zfir(ngtc)
complex(8), intent(out) :: zfgp(ngp)
! local variables
integer is,ia,ias,ig
integer nrc,nrci,irco,irc
integer l,m,lm,i,j
real(8) t0,t1,t2
complex(8) z1,z2
! automatic arrays
real(8) fr1(nrcmtmax),fr2(nrcmtmax)
complex(8) ylm(lmmaxo)
! allocatable arrays
complex(8), allocatable :: zfft(:)
! external functions
real(8) splint
external splint
if ((ngp.lt.1).or.(ngp.gt.ngvc)) then
  write(*,*)
  write(*,'("Error(zftzf): ngp out of range : ",I8)') ngp
  write(*,*)
  stop
end if
!-----------------------------------!
!     interstitial contribution     !
!-----------------------------------!
allocate(zfft(ngtc))
! multiply by coarse characteristic function
zfft(:)=zfir(:)*cfrc(:)
! Fourier transform to coarse G-grid
call zfftifc(3,ngdc,-1,zfft)
do ig=1,ngp
  zfgp(ig)=zfft(igfc(ig))
end do
deallocate(zfft)
!---------------------------------!
!     muffin-tin contribution     !
!---------------------------------!
t0=fourpi/omega
do ig=1,ngp
  lm=1
  do l=1,lmaxo
    z1=zilc(l)
    do m=-l,l
      lm=lm+1
      ylm(lm)=z1*ylmgp(lm,ig)
    end do
  end do
  do is=1,nspecies
    nrc=nrcmt(is)
    nrci=nrcmti(is)
    irco=nrci+1
    do ia=1,natoms(is)
      ias=idxas(ia,is)
      i=1
      j=1
! inner part of muffin-tin
      if (lmaxi.eq.1) then
        do irc=1,nrci
          z1=jlgpr(j,is,ig)*zfmt(i,ias)*y00+jlgpr(j+1,is,ig) &
           *(zfmt(i+1,ias)*ylm(2)+zfmt(i+2,ias)*ylm(3)+zfmt(i+3,ias)*ylm(4))
          z1=z1*r2cmt(irc,is)
          fr1(irc)=dble(z1)
          fr2(irc)=aimag(z1)
          i=i+4
          j=j+2
        end do
      else
        do irc=1,nrci
          z1=jlgpr(j,is,ig)*zfmt(i,ias)*y00
          i=i+1
          j=j+1
          lm=1
          do l=1,lmaxi
            lm=lm+1
            z2=zfmt(i,ias)*ylm(lm)
            i=i+1
            do m=1-l,l
              lm=lm+1
              z2=z2+zfmt(i,ias)*ylm(lm)
              i=i+1
            end do
            z1=z1+jlgpr(j,is,ig)*z2
            j=j+1
          end do
          z1=z1*r2cmt(irc,is)
          fr1(irc)=dble(z1)
          fr2(irc)=aimag(z1)
        end do
      end if
! outer part of muffin-tin
      do irc=irco,nrc
        z1=jlgpr(j,is,ig)*zfmt(i,ias)*y00
        i=i+1
        j=j+1
        lm=1
        do l=1,lmaxo
          lm=lm+1
          z2=zfmt(i,ias)*ylm(lm)
          i=i+1
          do m=1-l,l
            lm=lm+1
            z2=z2+zfmt(i,ias)*ylm(lm)
            i=i+1
          end do
          z1=z1+jlgpr(j,is,ig)*z2
          j=j+1
        end do
        z1=z1*r2cmt(irc,is)
        fr1(irc)=dble(z1)
        fr2(irc)=aimag(z1)
      end do
      t1=splint(nrc,rcmt(:,is),fr1)
      t2=splint(nrc,rcmt(:,is),fr2)
      zfgp(ig)=zfgp(ig)+t0*conjg(sfacgp(ig,ias))*cmplx(t1,t2,8)
    end do
  end do
end do
return
end subroutine