1
2! Copyright (C) 2013 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 dwavefmt(lrstp,ias,ngp,ngpq,apwalmq,dapwalm,evecfv,devecfv,dwfmt)
7use modmain
8use modphonon
9implicit none
10! arguments
11integer, intent(in) :: lrstp,ias,ngp,ngpq
12complex(8), intent(in) :: apwalmq(ngkmax,apwordmax,lmmaxapw)
13complex(8), intent(in) :: dapwalm(ngkmax,apwordmax,lmmaxapw)
14complex(8), intent(in) :: evecfv(nmatmax),devecfv(nmatmax)
15real(8), intent(out) :: dwfmt(2,*)
16! local variables
17integer is,ldi,ldo,io,ilo
18integer nrc,nrci,nrco,iro
19integer l,lm,npc,npci,i
20complex(8) z1
21! external functions
22complex(8), external :: zdotu
23is=idxis(ias)
24ldi=2*lmmaxi
25ldo=2*lmmaxo
26iro=nrmti(is)+lrstp
27if (lrstp.eq.1) then
28  nrc=nrmt(is)
29  nrci=nrmti(is)
30  npc=npmt(is)
31  npci=npmti(is)
32else
33  nrc=nrcmt(is)
34  nrci=nrcmti(is)
35  npc=npcmt(is)
36  npci=npcmti(is)
37end if
38nrco=nrc-nrci
39! zero the wavefunction derivative
40dwfmt(:,1:npc)=0.d0
41!-----------------------!
42!     APW functions     !
43!-----------------------!
44do l=0,lmaxo
45  do lm=l**2+1,(l+1)**2
46    i=npci+lm
47    do io=1,apword(l,is)
48      z1=zdotu(ngpq,devecfv,1,apwalmq(:,io,lm),1)
49      if (ias.eq.iasph) then
50        z1=z1+zdotu(ngp,evecfv,1,dapwalm(:,io,lm),1)
51      end if
52      if (abs(dble(z1)).gt.1.d-14) then
53        if (l.le.lmaxi) then
54          call daxpy(nrci,dble(z1),apwfr(1,1,io,l,ias),lrstp,dwfmt(1,lm),ldi)
55        end if
56        call daxpy(nrco,dble(z1),apwfr(iro,1,io,l,ias),lrstp,dwfmt(1,i),ldo)
57      end if
58      if (abs(aimag(z1)).gt.1.d-14) then
59        if (l.le.lmaxi) then
60          call daxpy(nrci,aimag(z1),apwfr(1,1,io,l,ias),lrstp,dwfmt(2,lm),ldi)
61        end if
62        call daxpy(nrco,aimag(z1),apwfr(iro,1,io,l,ias),lrstp,dwfmt(2,i),ldo)
63      end if
64    end do
65  end do
66end do
67!---------------------------------!
68!     local-orbital functions     !
69!---------------------------------!
70do ilo=1,nlorb(is)
71  l=lorbl(ilo,is)
72  do lm=l**2+1,(l+1)**2
73    i=npci+lm
74    z1=devecfv(ngpq+idxlo(lm,ilo,ias))
75    if (abs(dble(z1)).gt.1.d-14) then
76      if (l.le.lmaxi) then
77        call daxpy(nrci,dble(z1),lofr(1,1,ilo,ias),lrstp,dwfmt(1,lm),ldi)
78      end if
79      call daxpy(nrco,dble(z1),lofr(iro,1,ilo,ias),lrstp,dwfmt(1,i),ldo)
80    end if
81    if (abs(aimag(z1)).gt.1.d-14) then
82      if (l.le.lmaxi) then
83        call daxpy(nrci,aimag(z1),lofr(1,1,ilo,ias),lrstp,dwfmt(2,lm),ldi)
84      end if
85      call daxpy(nrco,aimag(z1),lofr(iro,1,ilo,ias),lrstp,dwfmt(2,i),ldo)
86    end if
87  end do
88end do
89end subroutine
90
91