File: genidxbse.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 (154 lines) | stat: -rw-r--r-- 4,017 bytes parent folder | download | duplicates (2)
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

! Copyright (C) 2010 S. Sharma, J. K. Dewhurst and E. K. U. Gross.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

subroutine genidxbse
use modmain
implicit none
integer ik,jk,a,ntop
integer ist,jst,i,j,k
! allocatable arrays
integer, allocatable :: idx(:)
! check if the BSE extra valence or conduction states are in range
do i=1,nvxbse
  ist=istxbse(i)
  if ((ist.lt.1).or.(ist.gt.nstsv)) then
    write(*,*)
    write(*,'("Error(genidxbse): extra valence state out of range : ",I8)') ist
    write(*,*)
    stop
  end if
end do
do j=1,ncxbse
  jst=jstxbse(j)
  if ((jst.lt.1).or.(jst.gt.nstsv)) then
    write(*,*)
    write(*,'("Error(genidxbse): extra conduction state out of range : ",I8)') &
     jst
    write(*,*)
    stop
  end if
end do
! number of valence states for transitions
nvbse=nvbse0+nvxbse
! number of conduction states for transitions
ncbse=ncbse0+ncxbse
if ((nvbse.le.0).or.(ncbse.le.0)) then
  write(*,*)
  write(*,'("Error(genidxbse): invalid number of valence or conduction &
   &transition states : ",2I8)') nvbse,ncbse
  write(*,*)
  stop
end if
! total number of transitions
nvcbse=nvbse*ncbse
! block size in BSE matrix
nbbse=nvcbse*nkptnr
! BSE matrix size
if (bsefull) then
  nmbse=2*nbbse
else
  nmbse=nbbse
end if
allocate(idx(nstsv))
! allocate global BSE index arrays
if (allocated(istbse)) deallocate(istbse)
allocate(istbse(nvbse,nkptnr))
if (allocated(jstbse)) deallocate(jstbse)
allocate(jstbse(ncbse,nkptnr))
if (allocated(ijkbse)) deallocate(ijkbse)
allocate(ijkbse(nvbse,ncbse,nkptnr))
a=0
! loop over non-reduced k-points
do ik=1,nkptnr
! equivalent reduced k-point
  jk=ivkik(ivk(1,ik),ivk(2,ik),ivk(3,ik))
! index for sorting the eigenvalues into ascending order
  call sortidx(nstsv,evalsv(:,jk),idx)
! find the topmost occupied band
  ntop=nstsv
  do ist=nstsv,1,-1
    if (evalsv(idx(ist),jk).lt.efermi) then
      ntop=ist
      exit
    end if
  end do
  if ((ntop-nvbse0+1).lt.1) then
    write(*,*)
    write(*,'("Error(genidxbse): not enough valence states, reduce nvbse")')
    write(*,*)
    stop
  end if
  if ((ntop+ncbse0).gt.nstsv) then
    write(*,*)
    write(*,'("Error(genidxbse): not enough conduction states")')
    write(*,'(" reduce ncbse or increase nempty")')
    write(*,*)
    stop
  end if
! index from BSE valence states to second-variational state numbers
  do i=1,nvbse0
    istbse(i,ik)=idx(ntop-nvbse0+i)
  end do
! index from BSE conduction states to second-variational state numbers
  do j=1,ncbse0
    jstbse(j,ik)=idx(ntop+j)
  end do
! add extra states to the list
  do i=1,nvxbse
    ist=istxbse(i)
    if (evalsv(ist,jk).gt.efermi) then
      write(*,*)
      write(*,'("Error(genidxbse): extra valence state above Fermi energy : ",&
       &I6)') ist
      write(*,'(" for k-point ",I8)') jk
      write(*,*)
      stop
    end if
    do k=1,nvbse0+i-1
      if (ist.eq.istbse(k,ik)) then
        write(*,*)
        write(*,'("Error(genidxbse): redundant extra valence state : ",I6)') ist
        write(*,'(" for k-point ",I8)') jk
        write(*,*)
        stop
      end if
    end do
    istbse(nvbse0+i,ik)=ist
  end do
  do j=1,ncxbse
    jst=jstxbse(j)
    if (evalsv(jst,jk).lt.efermi) then
      write(*,*)
      write(*,'("Error(genidxbse): extra conduction state below Fermi &
       &energy : ",I6)') jst
      write(*,'(" for k-point ",I8)') jk
      write(*,*)
      stop
    end if
    do k=1,ncbse0+j-1
      if (jst.eq.jstbse(k,ik)) then
        write(*,*)
        write(*,'("Error(genidxbse): redundant extra conduction state : ",&
         &I6)') jst
        write(*,'(" for k-point ",I8)') jk
        write(*,*)
        stop
      end if
    end do
    jstbse(ncbse0+j,ik)=jst
  end do
! index from BSE valence-conduction pair and k-point to location in BSE matrix
  do i=1,nvbse
    do j=1,ncbse
      a=a+1
      ijkbse(i,j,ik)=a
    end do
  end do
! end loop over non-reduced k-points
end do
deallocate(idx)
return
end subroutine