File: bandstr.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 (182 lines) | stat: -rw-r--r-- 5,394 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
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182

! Copyright (C) 2002-2005 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.

!BOP
! !ROUTINE: bandstr
! !INTERFACE:
subroutine bandstr
! !USES:
use modmain
use modomp
! !DESCRIPTION:
!   Produces a band structure along the path in reciprocal space which connects
!   the vertices in the array {\tt vvlp1d}. The band structure is obtained from
!   the second-variational eigenvalues and is written to the file {\tt BAND.OUT}
!   with the Fermi energy set to zero. If required, band structures are plotted
!   to files {\tt BAND\_Sss\_Aaaaa.OUT} for atom {\tt aaaa} of species {\tt ss},
!   which include the band characters for each $l$ component of that atom in
!   columns 4 onwards. Column 3 contains the sum over $l$ of the characters.
!   Vertex location lines are written to {\tt BANDLINES.OUT}.
!
! !REVISION HISTORY:
!   Created June 2003 (JKD)
!EOP
!BOC
implicit none
! local variables
integer ik,ist,ispn,is,ia,ias
integer lmax,lmmax,l,m,lm,iv,nthd
real(8) emin,emax,sum
character(256) fname
! allocatable arrays
real(8), allocatable :: evalfv(:,:),e(:,:)
! low precision for band character array saves memory
real(4), allocatable :: bc(:,:,:,:)
complex(8), allocatable :: dmat(:,:,:,:,:),apwalm(:,:,:,:,:)
complex(8), allocatable :: evecfv(:,:,:),evecsv(:,:)
! initialise universal variables
call init0
call init1
! allocate array for storing the eigenvalues
allocate(e(nstsv,nkpt))
! maximum angular momentum for band character
lmax=min(3,lmaxo)
lmmax=(lmax+1)**2
if (task.eq.21) then
  allocate(bc(0:lmax,natmtot,nstsv,nkpt))
end if
! read density and potentials from file
call readstate
! Fourier transform Kohn-Sham potential to G-space
call genvsig
! read Fermi energy from file
call readfermi
! find the new linearisation energies
call linengy
! generate the APW and local-orbital radial functions and integrals
call genapwlofr
! generate the spin-orbit coupling radial functions
call gensocfr
emin=1.d5
emax=-1.d5
! begin parallel loop over k-points
call omp_hold(nkpt,nthd)
!$OMP PARALLEL DEFAULT(SHARED) &
!$OMP PRIVATE(evalfv,evecfv,evecsv) &
!$OMP PRIVATE(ist,dmat,apwalm,ispn) &
!$OMP PRIVATE(ias,l,m,lm,sum) &
!$OMP NUM_THREADS(nthd)
!$OMP DO
do ik=1,nkpt
  allocate(evalfv(nstfv,nspnfv))
  allocate(evecfv(nmatmax,nstfv,nspnfv))
  allocate(evecsv(nstsv,nstsv))
!$OMP CRITICAL(bandstr_1)
  write(*,'("Info(bandstr): ",I6," of ",I6," k-points")') ik,nkpt
!$OMP END CRITICAL(bandstr_1)
! solve the first- and second-variational eigenvalue equations
  call eveqn(ik,evalfv,evecfv,evecsv)
  do ist=1,nstsv
! subtract the Fermi energy
    e(ist,ik)=evalsv(ist,ik)-efermi
!$OMP CRITICAL(bandstr_2)
    emin=min(emin,e(ist,ik))
    emax=max(emax,e(ist,ik))
!$OMP END CRITICAL(bandstr_2)
  end do
! compute the band characters if required
  if (task.eq.21) then
    allocate(dmat(lmmax,nspinor,lmmax,nspinor,nstsv))
    allocate(apwalm(ngkmax,apwordmax,lmmaxapw,natmtot,nspnfv))
! 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
! average band character over spin and m for all atoms
    do ias=1,natmtot
! generate the diagonal of the density matrix
      call gendmat(.true.,.true.,0,lmax,ias,ngk(:,ik),apwalm,evecfv,evecsv, &
       lmmax,dmat)
      do ist=1,nstsv
        do l=0,lmax
          sum=0.d0
          do m=-l,l
            lm=idxlm(l,m)
            do ispn=1,nspinor
              sum=sum+dble(dmat(lm,ispn,lm,ispn,ist))
            end do
          end do
          bc(l,ias,ist,ik)=real(sum)
        end do
      end do
    end do
    deallocate(dmat,apwalm)
  end if
  deallocate(evalfv,evecfv,evecsv)
! end loop over k-points
end do
!$OMP END DO
!$OMP END PARALLEL
call omp_free(nthd)
emax=emax+(emax-emin)*0.5d0
emin=emin-(emax-emin)*0.5d0
! output the band structure
if (task.eq.20) then
  open(50,file='BAND.OUT',form='FORMATTED')
  do ist=1,nstsv
    do ik=1,nkpt
      write(50,'(2G18.10)') dpp1d(ik),e(ist,ik)
    end do
    write(50,'("     ")')
  end do
  close(50)
  write(*,*)
  write(*,'("Info(bandstr):")')
  write(*,'(" band structure plot written to BAND.OUT")')
else
  do is=1,nspecies
    do ia=1,natoms(is)
      ias=idxas(ia,is)
      write(fname,'("BAND_S",I2.2,"_A",I4.4,".OUT")') is,ia
      open(50,file=trim(fname),form='FORMATTED')
      do ist=1,nstsv
        do ik=1,nkpt
! sum band character over l
          sum=0.d0
          do l=0,lmax
            sum=sum+bc(l,ias,ist,ik)
          end do
          write(50,'(2G18.10,8F12.6)') dpp1d(ik),e(ist,ik),sum, &
           (bc(l,ias,ist,ik),l=0,lmax)
        end do
        write(50,'("     ")')
      end do
      close(50)
    end do
  end do
  write(*,*)
  write(*,'("Info(bandstr):")')
  write(*,'(" band structure plot written to BAND_Sss_Aaaaa.OUT")')
  write(*,'(" for all species and atoms")')
end if
write(*,*)
write(*,'(" Fermi energy is at zero in plot")')
! output the vertex location lines
open(50,file='BANDLINES.OUT',form='FORMATTED')
do iv=1,nvp1d
  write(50,'(2G18.10)') dvp1d(iv),emin
  write(50,'(2G18.10)') dvp1d(iv),emax
  write(50,'("     ")')
end do
close(50)
write(*,*)
write(*,'(" vertex location lines written to BANDLINES.OUT")')
deallocate(e)
if (task.eq.21) deallocate(bc)
return
end subroutine
!EOC