File: writeexpmat.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 (95 lines) | stat: -rw-r--r-- 2,582 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

! Copyright (C) 2008 J. K. Dewhurst, S. Sharma and C. Ambrosch-Draxl.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

subroutine writeexpmat
use modmain
implicit none
! local variables
integer nk,ik,jk,i,j
real(8) vgqc(3),gqc
real(8) a,b
! allocatable arrays
real(8), allocatable :: jlgqr(:,:)
complex(8), allocatable :: ylmgq(:),sfacgq(:)
complex(8), allocatable :: expmt(:,:),emat(:,:)
! initialise universal variables
call init0
call init1
call init2
! read in the density and potentials from file
call readstate
! read Fermi energy from file
call readfermi
! find the new linearisation energies
call linengy
! generate the APW radial functions
call genapwfr
! generate the local-orbital radial functions
call genlofr
! generate the phase factor function exp(iq.r) in the muffin-tins
allocate(jlgqr(njcmax,nspecies))
allocate(ylmgq(lmmaxo),sfacgq(natmtot))
allocate(expmt(npcmtmax,natmtot))
ngrf=1
call gengqrf(vecqc,vgqc,gqc,jlgqr,ylmgq,sfacgq)
call genexpmt(1,jlgqr,ylmgq,1,sfacgq,expmt)
expmt(:,:)=omega*expmt(:,:)
deallocate(jlgqr,ylmgq,sfacgq)
! number of k-points to write out
if (kstlist(1,1).le.0) then
  nk=nkpt
else
  nk=nkstlist
end if
open(50,file='EXPIQR.OUT',form='FORMATTED')
write(50,*)
write(50,'("q-vector (lattice coordinates) :")')
write(50,'(3G18.10)') vecql
write(50,'("q-vector (Cartesian coordinates) :")')
write(50,'(3G18.10)') vecqc
write(50,*)
write(50,'(I8," : number of k-points")') nk
write(50,'(I6," : number of states per k-point")') nstsv
allocate(emat(nstsv,nstsv))
do jk=1,nk
  if (kstlist(1,1).le.0) then
    ik=jk
  else
    ik=kstlist(1,jk)
  end if
  if ((ik.le.0).or.(ik.gt.nkpt)) then
    write(*,*)
    write(*,'("Error(writeexpiqr): k-point out of range : ",I8)') ik
    write(*,*)
    stop
  end if
  write(50,*)
  write(50,'(" k-point (lattice coordinates) :")')
  write(50,'(3G18.10)') vkl(:,ik)
  write(50,*)
  write(50,'(" k-point (Cartesian coordinates) :")')
  write(50,'(3G18.10)') vkc(:,ik)
  call genexpmat(vkl(:,ik),expmt,emat)
  do i=1,nstsv
    write(50,*)
    write(50,'(I6," : state i; state j, <...>, |<...>|^2 below")') i
    do j=1,nstsv
      a=dble(emat(i,j))
      b=aimag(emat(i,j))
      write(50,'(I6,3G18.10)') j,a,b,a**2+b**2
    end do
  end do
! end loop over k-points
end do
close(50)
write(*,*)
write(*,'("Info(writeexpmat)")')
write(*,'(" < i,k+q | exp(iq.r) | j,k > matrix elements written to &
 &EXPIQR.OUT")')
write(*,'(" for the q-vector in vecql and all k-points in kstlist")')
deallocate(expmt,emat)
return
end subroutine