File: gentauk.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 (99 lines) | stat: -rw-r--r-- 2,933 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

! Copyright (C) 2011 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 gentauk(ik)
use modmain
implicit none
! arguments
integer, intent(in) :: ik
! local variables
integer ispn,jspn,nst,ist,jst
integer is,ias,nrc,nrci
integer npc,igk,ifg,i
real(8) wo
complex(8) z1
! automatic arrays
integer idx(nstsv)
! allocatable arrays
complex(8), allocatable :: apwalm(:,:,:,:,:),evecfv(:,:),evecsv(:,:)
complex(8), allocatable :: wfmt(:,:,:,:),wfgp(:,:,:)
complex(8), allocatable :: gzfmt(:,:),zfmt(:),zfft(:)
allocate(apwalm(ngkmax,apwordmax,lmmaxapw,natmtot,nspnfv))
allocate(evecfv(nmatmax,nstfv),evecsv(nstsv,nstsv))
! find the matching coefficients
do ispn=1,nspnfv
  call match(ngk(ispn,ik),gkc(:,ispn,ik),tpgkc(:,:,ispn,ik), &
   sfacgk(:,:,ispn,ik),apwalm(:,:,:,:,ispn))
end do
! get the eigenvectors from file
call getevecfv(filext,ik,vkl(:,ik),vgkl(:,:,:,ik),evecfv)
call getevecsv(filext,ik,vkl(:,ik),evecsv)
! count and index the occupied states
nst=0
do ist=1,nstsv
  if (abs(occsv(ist,ik)).lt.epsocc) cycle
  nst=nst+1
  idx(nst)=ist
end do
! calculate the second-variational wavefunctions for occupied states
allocate(wfmt(npcmtmax,natmtot,nspinor,nst),wfgp(ngkmax,nspinor,nst))
call genwfsv(.true.,.true.,nst,idx,ngridg,igfft,ngk(:,ik),igkig(:,:,ik), &
 apwalm,evecfv,evecsv,wfmt,ngkmax,wfgp)
deallocate(apwalm,evecfv,evecsv)
!-------------------------!
!     muffin-tin part     !
!-------------------------!
allocate(gzfmt(npcmtmax,3),zfmt(npcmtmax))
do ist=1,nst
  jst=idx(ist)
  wo=0.5d0*wkpt(ik)*occsv(jst,ik)
  do ispn=1,nspinor
    do ias=1,natmtot
      is=idxis(ias)
      nrc=nrcmt(is)
      nrci=nrcmti(is)
      npc=npcmt(is)
! compute the gradient of the wavefunction
      call gradzfmt(nrc,nrci,rcmt(:,is),wfmt(:,ias,ispn,ist),npcmtmax,gzfmt)
      do i=1,3
! convert gradient to spherical coordinates
        call zbsht(nrc,nrci,gzfmt(:,i),zfmt)
! add to total in muffin-tin
!$OMP CRITICAL(gentauk_1)
        taumt(1:npc,ias,ispn)=taumt(1:npc,ias,ispn) &
         +wo*(dble(zfmt(1:npc))**2+aimag(zfmt(1:npc))**2)
!$OMP END CRITICAL(gentauk_1)
      end do
    end do
  end do
end do
deallocate(wfmt,gzfmt,zfmt)
!---------------------------!
!     interstitial part     !
!---------------------------!
allocate(zfft(ngtot))
do ist=1,nst
  jst=idx(ist)
  wo=0.5d0*wkpt(ik)*occsv(jst,ik)/omega
  do ispn=1,nspinor
    jspn=jspnfv(ispn)
    do i=1,3
      zfft(:)=0.d0
      do igk=1,ngk(jspn,ik)
        ifg=igfft(igkig(igk,jspn,ik))
        z1=wfgp(igk,ispn,ist)
        zfft(ifg)=vgkc(i,igk,jspn,ik)*cmplx(-aimag(z1),dble(z1),8)
      end do
      call zfftifc(3,ngridg,1,zfft)
!$OMP CRITICAL(gentauk_2)
      tauir(:,ispn)=tauir(:,ispn)+wo*(dble(zfft(:))**2+aimag(zfft(:))**2)
!$OMP END CRITICAL(gentauk_2)
    end do
  end do
end do
deallocate(wfgp,zfft)
return
end subroutine