File: init4.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 (200 lines) | stat: -rw-r--r-- 6,521 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
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200

! Copyright (C) 2012 J. K. Dewhurst, S. Sharma 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 init4
use modmain
use modphonon
use modpw
use modvars
implicit none
! local variables
integer ik,ihk,jspn
integer l1,l2,l3,m1,m2,m3
integer lm1,lm2,lm3
real(8) vl(3),vc(3)
! external functions
real(8) gaunt
external gaunt

!---------------------------!
!     H+k-vector arrays     !
!---------------------------!
if ((task.eq.135).or.(task.eq.170).or.(task.eq.171).or.(task.eq.172).or. &
 (task.eq.173)) then
  if (task.eq.135) hkmax=0.5d0*gmaxvr-epslat
  call findngkmax(nkpt,vkc,nspnfv,vqcss,ngvec,vgc,hkmax,nhkmax)
! allocate the H+k-vector arrays
  if (allocated(nhk)) deallocate(nhk)
  allocate(nhk(nspnfv,nkpt))
  if (allocated(ihkig)) deallocate(ihkig)
  allocate(ihkig(nhkmax,nspnfv,nkpt))
  if (allocated(vhkl)) deallocate(vhkl)
  allocate(vhkl(3,nhkmax,nspnfv,nkpt))
  if (allocated(vhkc)) deallocate(vhkc)
  allocate(vhkc(3,nhkmax,nspnfv,nkpt))
  if (allocated(hkc)) deallocate(hkc)
  allocate(hkc(nhkmax,nspnfv,nkpt))
  if (allocated(tphkc)) deallocate(tphkc)
  allocate(tphkc(2,nhkmax,nspnfv,nkpt))
  if (allocated(sfachk)) deallocate(sfachk)
  allocate(sfachk(nhkmax,natmtot,nspnfv,nkpt))
! initialise H+k-vectors arrays
  do ik=1,nkpt
    do jspn=1,nspnfv
      vl(:)=vkl(:,ik)
      vc(:)=vkc(:,ik)
! spin-spiral case
      if (spinsprl) then
        if (jspn.eq.1) then
          vl(:)=vl(:)+0.5d0*vqlss(:)
          vc(:)=vc(:)+0.5d0*vqcss(:)
        else
          vl(:)=vl(:)-0.5d0*vqlss(:)
          vc(:)=vc(:)-0.5d0*vqcss(:)
        end if
      end if
! generate the H+k-vectors
      call gengkvec(ngvec,ivg,vgc,vl,vc,hkmax,nhkmax,nhk(jspn,ik), &
       ihkig(:,jspn,ik),vhkl(:,:,jspn,ik),vhkc(:,:,jspn,ik))
! generate the spherical coordinates of the H+k-vectors
      do ihk=1,nhk(jspn,ik)
        call sphcrd(vhkc(:,ihk,jspn,ik),hkc(ihk,jspn,ik),tphkc(:,ihk,jspn,ik))
      end do
! generate structure factors for H+k-vectors
      call gensfacgp(nhk(jspn,ik),vhkc(:,:,jspn,ik),nhkmax,sfachk(:,:,jspn,ik))
    end do
  end do
! write to VARIABLES.OUT
  call writevars('hkmax',rv=hkmax)
  call writevars('nhk',nv=nspnfv*nkpt,iva=nhk)
  do ik=1,nkpt
    do jspn=1,nspnfv
      call writevars('ihkig',l=jspn,m=ik,nv=nhk(jspn,ik),iva=ihkig(:,jspn,ik))
    end do
  end do
end if

!-----------------------------!
!     G+k+q-vector arrays     !
!-----------------------------!
if ((task.eq.205).or.(task.eq.240)) then
  if (allocated(vkql)) deallocate(vkql)
  allocate(vkql(3,nkptnr))
  if (allocated(vkqc)) deallocate(vkqc)
  allocate(vkqc(3,nkptnr))
  if (allocated(ngkq)) deallocate(ngkq)
  allocate(ngkq(nspnfv,nkptnr))
  if (allocated(igkqig)) deallocate(igkqig)
  allocate(igkqig(ngkmax,nspnfv,nkptnr))
  if (allocated(vgkql)) deallocate(vgkql)
  allocate(vgkql(3,ngkmax,nspnfv,nkptnr))
  if (allocated(vgkqc)) deallocate(vgkqc)
  allocate(vgkqc(3,ngkmax,nspnfv,nkptnr))
  if (allocated(gkqc)) deallocate(gkqc)
  allocate(gkqc(ngkmax,nspnfv,nkptnr))
  if (allocated(tpgkqc)) deallocate(tpgkqc)
  allocate(tpgkqc(2,ngkmax,nspnfv,nkptnr))
  if (allocated(sfacgkq)) deallocate(sfacgkq)
  allocate(sfacgkq(ngkmax,natmtot,nspnfv,nkptnr))
end if

!---------------------------!
!     G+q-vector arrays     !
!---------------------------!
if ((task.eq.205).or.(task.eq.240)) then
  if (allocated(vgqc)) deallocate(vgqc)
  allocate(vgqc(3,ngtot))
  if (allocated(gqc)) deallocate(gqc)
  allocate(gqc(ngtot))
  if (allocated(gclgq)) deallocate(gclgq)
  allocate(gclgq(ngvec))
  if (allocated(jlgqrmt)) deallocate(jlgqrmt)
  allocate(jlgqrmt(0:lnpsd,ngvec,nspecies))
  if (allocated(ylmgq)) deallocate(ylmgq)
  allocate(ylmgq(lmmaxo,ngtot))
  if (allocated(sfacgq)) deallocate(sfacgq)
  allocate(sfacgq(ngvec,natmtot))
  if (allocated(ffacgq)) deallocate(ffacgq)
  allocate(ffacgq(ngtot,nspecies))
  if (allocated(dcfunig)) deallocate(dcfunig)
  allocate(dcfunig(ngtot))
  if (allocated(dcfunir)) deallocate(dcfunir)
  allocate(dcfunir(ngtot))
end if

!-----------------------------------------------------------------!
!     phonon density functional perturbation theory variables     !
!-----------------------------------------------------------------!
if (task.eq.205) then
  if (allocated(drhomt)) deallocate(drhomt)
  allocate(drhomt(npmtmax,natmtot))
  if (allocated(drhoir)) deallocate(drhoir)
  allocate(drhoir(ngtot))
  if (allocated(dmagmt)) deallocate(dmagmt)
  if (allocated(dmagir)) deallocate(dmagir)
  if (spinpol) then
    allocate(dmagmt(npmtmax,natmtot,ndmag))
    allocate(dmagir(ngtot,ndmag))
  end if
  if (allocated(dvclmt)) deallocate(dvclmt)
  allocate(dvclmt(npmtmax,natmtot))
  if (allocated(dvclir)) deallocate(dvclir)
  allocate(dvclir(ngtot))
  if (allocated(zvnmt)) deallocate(zvnmt)
  allocate(zvnmt(npmtmax))
  if (allocated(dvsmt)) deallocate(dvsmt)
  allocate(dvsmt(npmtmax,natmtot))
  if (allocated(dvsir)) deallocate(dvsir)
  allocate(dvsir(ngtot))
  if (allocated(gvsmt)) deallocate(gvsmt)
  allocate(gvsmt(npmtmax))
  if (allocated(dvsig)) deallocate(dvsig)
  allocate(dvsig(ngtot))
  if (allocated(dbsmt)) deallocate(dbsmt)
  if (allocated(dbsir)) deallocate(dbsir)
  if (spinpol) then
    allocate(dbsmt(npcmtmax,natmtot,ndmag))
    allocate(dbsir(ngtot,ndmag))
  end if
  if (allocated(dsocfr)) deallocate(dsocfr)
  if (spinorb) then
    allocate(dsocfr(nrcmtmax,natmtot))
  end if
  if (allocated(dhaa)) deallocate(dhaa)
  allocate(dhaa(lmmaxo,apwordmax,0:lmaxapw,apwordmax,0:lmaxapw,natmtot))
  if (allocated(dhloa)) deallocate(dhloa)
  allocate(dhloa(lmmaxo,apwordmax,0:lmaxapw,nlomax,natmtot))
  if (allocated(dhlolo)) deallocate(dhlolo)
  allocate(dhlolo(lmmaxo,nlomax,nlomax,natmtot))
! allocate and generate real Gaunt coefficient array
  if (allocated(gntyyy)) deallocate(gntyyy)
  allocate(gntyyy(lmmaxapw,lmmaxo,lmmaxapw))
  do l1=0,lmaxapw
    do m1=-l1,l1
      lm1=idxlm(l1,m1)
      do l2=0,lmaxo
        do m2=-l2,l2
          lm2=idxlm(l2,m2)
          do l3=0,lmaxapw
            do m3=-l3,l3
              lm3=idxlm(l3,m3)
              gntyyy(lm1,lm2,lm3)=gaunt(l1,l2,l3,m1,m2,m3)
            end do
          end do
        end do
      end do
    end do
  end do
  if (allocated(devalfv)) deallocate(devalfv)
  allocate(devalfv(nstfv,nspnfv,nkptnr))
  if (allocated(devalsv)) deallocate(devalsv)
  allocate(devalsv(nstsv,nkptnr))
  if (allocated(doccsv)) deallocate(doccsv)
  allocate(doccsv(nstsv,nkptnr))
end if

return
end subroutine