1 2! Copyright (C) 2015 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 genvmmtsv(wfmt,vmat) 7use modmain 8use moddftu 9use modomp 10implicit none 11! arguments 12complex(8), intent(in) :: wfmt(npcmtmax,natmtot,nspinor,nstsv) 13complex(8), intent(inout) :: vmat(nstsv,nstsv) 14! local variables 15integer ist,jst,ispn,jspn 16integer is,ias,ld,nthd 17integer nrc,nrci,nrco,i 18integer nm,l,lm,npc,npci 19! allocatable arrays 20complex(8), allocatable :: wfmt1(:,:,:),wfmt2(:) 21! external functions 22complex(8), external :: zfmtinp 23if (.not.tvmatmt) return 24ld=lmmaxdm*nspinor 25allocate(wfmt1(npcmtmax,nspinor,nstsv),wfmt2(npcmtmax)) 26! loop over atoms 27do ias=1,natmtot 28 is=idxis(ias) 29 nrc=nrcmt(is) 30 nrci=nrcmti(is) 31 nrco=nrc-nrci 32 npc=npcmt(is) 33 npci=npcmti(is) 34! convert wavefunctions to spherical harmonics 35 do ist=1,nstsv 36 do ispn=1,nspinor 37 call zfsht(nrc,nrci,wfmt(:,ias,ispn,ist),wfmt1(:,ispn,ist)) 38 end do 39 end do 40! loop over second-variational states 41 do jst=1,nstsv 42! loop over spins 43 do ispn=1,nspinor 44 wfmt2(1:npc)=0.d0 45 do jspn=1,nspinor 46 do l=0,lmaxdm 47 if (tvmmt(l,ias)) then 48 nm=2*l+1 49 lm=l**2+1 50 if (l.le.lmaxi) then 51 call zgemm('N','N',nm,nrci,nm,zone,vmatmt(lm,ispn,lm,jspn,ias), & 52 ld,wfmt1(lm,jspn,jst),lmmaxi,zone,wfmt2(lm),lmmaxi) 53 end if 54 i=npci+lm 55 call zgemm('N','N',nm,nrco,nm,zone,vmatmt(lm,ispn,lm,jspn,ias),ld, & 56 wfmt1(i,jspn,jst),lmmaxo,zone,wfmt2(i),lmmaxo) 57 end if 58 end do 59 end do 60! compute the matrix elements 61 call holdthd(nstsv,nthd) 62!$OMP PARALLEL DO DEFAULT(SHARED) & 63!$OMP NUM_THREADS(nthd) 64 do ist=1,nstsv 65 vmat(ist,jst)=vmat(ist,jst)+zfmtinp(nrc,nrci,wrcmt(:,is), & 66 wfmt1(:,ispn,ist),wfmt2) 67 end do 68!$OMP END PARALLEL DO 69 call freethd(nthd) 70 end do 71 end do 72end do 73deallocate(wfmt1,wfmt2) 74end subroutine 75 76