File: gendmat.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 (137 lines) | stat: -rw-r--r-- 4,082 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
135
136
137

! Copyright (C) 2007 F. Bultmark, F. Cricchio, L. Nordstrom and J. K. Dewhurst.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

subroutine gendmat(tspndg,tlmdg,lmin,lmax,ias,ngp,apwalm,evecfv,evecsv,ld,dmat)
use modmain
implicit none
! arguments
logical, intent(in) :: tspndg,tlmdg
integer, intent(in) :: lmin,lmax
integer, intent(in) :: ias
integer, intent(in) :: ngp(nspnfv)
complex(8), intent(in) :: apwalm(ngkmax,apwordmax,lmmaxapw,natmtot,nspnfv)
complex(8), intent(in) :: evecfv(nmatmax,nstfv,nspnfv)
complex(8), intent(in) :: evecsv(nstsv,nstsv)
integer, intent(in) :: ld
complex(8), intent(out) :: dmat(ld,nspinor,ld,nspinor,nstsv)
! local variables
integer ist,ispn,jspn,is,ia
integer nrc,nrci,nrco,irco,irc
integer l,m1,m2,lm1,lm2
integer npc,npci,i1,i2,i,j
real(8) a,b,t1
complex(8) zq(2),z1
! automatic arrays
logical done(nstfv,nspnfv)
real(8) fr1(nrcmtmax),fr2(nrcmtmax)
! allocatable arrays
complex(8), allocatable :: wfmt1(:,:,:),wfmt2(:,:)
! external functions
real(8) fintgt
external fintgt
if (lmin.lt.0) then
  write(*,*)
  write(*,'("Error(gendmat): lmin < 0 : ",I8)') lmin
  write(*,*)
  stop
end if
if (lmax.gt.lmaxo) then
  write(*,*)
  write(*,'("Error(gendmat): lmax > lmaxo : ",2I8)') lmax,lmaxo
  write(*,*)
  stop
end if
! allocate local arrays
allocate(wfmt1(npcmtmax,nstfv,nspnfv),wfmt2(npcmtmax,nspinor))
! species and atom numbers
is=idxis(ias)
ia=idxia(ias)
nrc=nrcmt(is)
nrci=nrcmti(is)
nrco=nrc-nrci
irco=nrci+1
npc=npcmt(is)
npci=npcmti(is)
! de-phasing factor for spin-spirals
if (ssdph) then
  t1=-0.5d0*dot_product(vqcss(:),atposc(:,ia,is))
  zq(1)=cmplx(cos(t1),sin(t1),8)
  zq(2)=conjg(zq(1))
end if
! zero the density matrix
dmat(:,:,:,:,:)=0.d0
done(:,:)=.false.
! begin loop over second-variational states
do j=1,nstsv
  if (tevecsv) then
! generate spinor wavefunction from second-variational eigenvectors
    wfmt2(1:npc,:)=0.d0
    i=0
    do ispn=1,nspinor
      jspn=jspnfv(ispn)
      do ist=1,nstfv
        i=i+1
        z1=evecsv(i,j)
        if (ssdph) z1=z1*zq(ispn)
        if (abs(dble(z1))+abs(aimag(z1)).gt.epsocc) then
          if (.not.done(ist,jspn)) then
            call wavefmt(lradstp,ias,ngp(jspn),apwalm(:,:,:,ias,jspn), &
             evecfv(:,ist,jspn),wfmt1(:,ist,jspn))
            done(ist,jspn)=.true.
          end if
! add to spinor wavefunction
          wfmt2(1:npc,ispn)=wfmt2(1:npc,ispn)+z1*wfmt1(1:npc,ist,jspn)
        end if
      end do
    end do
  else
! spin-unpolarised wavefunction
    call wavefmt(lradstp,ias,ngp,apwalm(:,:,:,ias,1),evecfv(:,j,1),wfmt2)
  end if
  do ispn=1,nspinor
    do jspn=1,nspinor
      if (tspndg.and.(ispn.ne.jspn)) cycle
      do l=lmin,lmax
        do m1=-l,l
          lm1=idxlm(l,m1)
          do m2=-l,l
            lm2=idxlm(l,m2)
            if (tlmdg.and.(lm1.ne.lm2)) cycle
            if (l.le.lmaxi) then
              i1=lm1; i2=lm2
              do irc=1,nrci
                z1=wfmt2(i1,ispn)*conjg(wfmt2(i2,jspn))*r2cmt(irc,is)
                fr1(irc)=dble(z1); fr2(irc)=aimag(z1)
                i1=i1+lmmaxi; i2=i2+lmmaxi
              end do
              do irc=irco,nrc
                z1=wfmt2(i1,ispn)*conjg(wfmt2(i2,jspn))*r2cmt(irc,is)
                fr1(irc)=dble(z1); fr2(irc)=aimag(z1)
                i1=i1+lmmaxo; i2=i2+lmmaxo
              end do
              a=fintgt(-2,nrc,rcmt(1,is),fr1)
              b=fintgt(-2,nrc,rcmt(1,is),fr2)
            else
              i1=npci+lm1; i2=npci+lm2
              do irc=irco,nrc
                z1=wfmt2(i1,ispn)*conjg(wfmt2(i2,jspn))*r2cmt(irc,is)
                fr1(irc)=dble(z1); fr2(irc)=aimag(z1)
                i1=i1+lmmaxo; i2=i2+lmmaxo
              end do
              a=fintgt(-2,nrco,rcmt(irco,is),fr1(irco))
              b=fintgt(-2,nrco,rcmt(irco,is),fr2(irco))
            end if
            dmat(lm1,ispn,lm2,jspn,j)=cmplx(a,b,8)
          end do
        end do
      end do
    end do
  end do
! end loop over second-variational states
end do
deallocate(wfmt1,wfmt2)
return
end subroutine