Line data Source code
1 : !--------------------------------------------------------------------------------
2 : ! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
3 : ! This file is part of FLEUR and available as free software under the conditions
4 : ! of the MIT license as expressed in the LICENSE file in more detail.
5 : !--------------------------------------------------------------------------------
6 :
7 : module m_wann_fft5
8 : contains
9 0 : subroutine wann_fft5(
10 : > l_conjugate,
11 : > dim1,dim2,
12 : > filename,title,
13 0 : > rvecnum,rvec,kpoints,
14 : > jspin,nkpts,film,
15 : > l_soc,band_min,band_max,
16 : > neigd,wan90version)
17 :
18 : c*************************************************
19 : c Transform 5-dimensional matrices from
20 : c Bloch representation to Wannier representation.
21 : c
22 : c Frank Freimuth, February 2011
23 : c*************************************************
24 :
25 : use m_constants
26 : use m_wann_read_umatrix
27 :
28 : implicit none
29 : logical, intent(in) :: l_conjugate
30 : integer, intent(in) :: dim1,dim2
31 :
32 : character,intent(in) :: filename*(*)
33 : character,intent(in) :: title*(*)
34 :
35 : integer, intent(in) :: rvecnum
36 : integer, intent(in) :: rvec(:,:)
37 : real, intent(in) :: kpoints(:,:)
38 :
39 : integer, intent(in) :: jspin
40 : integer, intent(in) :: nkpts
41 : logical, intent(in) :: film
42 :
43 : logical, intent(in) :: l_soc
44 : integer,intent(in) :: band_min(2),band_max(2),neigd
45 : integer, intent(in) :: wan90version
46 :
47 : integer :: ikpt
48 : integer :: kpts
49 : logical :: l_file
50 : integer :: num_wann,num_kpts,num_nnmax
51 : integer :: kspin,kkspin
52 : integer :: wann_shift,num_wann2
53 : integer :: i,j,k,m,info,r1,r2,r3,dummy1
54 : integer :: dummy2,dummy3,dummy4,dummy5,dummy6
55 : integer :: hopmin,hopmax,counter,m1,m2
56 : integer :: num_bands2
57 : integer,allocatable :: iwork(:)
58 : real,allocatable :: energy(:,:),ei(:)
59 : real,allocatable :: eigw(:,:),rwork(:)
60 : complex,allocatable :: work(:),vec(:,:)
61 0 : complex,allocatable :: u_matrix(:,:,:)
62 0 : complex,allocatable :: hwann(:,:,:,:,:)
63 0 : complex,allocatable :: hreal(:,:,:,:,:)
64 0 : complex,allocatable :: hsomtx(:,:,:,:,:)
65 0 : complex,allocatable :: hsomtx2(:,:,:,:,:)
66 : complex :: fac,eulav,eulav1
67 : real :: tmp_omi,rdotk,tpi,minenerg,maxenerg
68 : real, allocatable :: minieni(:),maxieni(:)
69 : character :: jobz,uplo
70 : integer :: kpt,band,lee,lwork,lrwork,liwork,n,lda
71 : complex :: value(4)
72 : logical :: um_format
73 : logical :: repro_eig
74 : logical :: l_chk,l_proj
75 : logical :: have_disentangled
76 0 : integer,allocatable :: ndimwin(:)
77 0 : logical,allocatable :: lwindow(:,:)
78 : integer :: chk_unit,nkp,ntmp,ierr
79 : character(len=33) :: header
80 : character(len=20) :: checkpoint
81 : real :: tmp_latt(3,3), tmp_kpt_latt(3,nkpts)
82 : real :: omega_invariant
83 0 : complex,allocatable :: u_matrix_opt(:,:,:)
84 : integer :: num_bands
85 : logical :: l_umdat
86 : real,allocatable :: eigval2(:,:)
87 : real,allocatable :: eigval_opt(:,:)
88 : real :: scale,a,b
89 : character(len=2) :: spinspin12(0:2)
90 : character(len=3) :: spin12(2)
91 : character(len=6) :: filename2
92 : integer :: jp,mp,kk,ii,jj,dir,dir2,rvecind
93 : integer :: spin1,spin2
94 :
95 : data spinspin12/' ','.1' , '.2'/
96 : data spin12/'WF1','WF2'/
97 :
98 0 : call timestart("wann_fft5")
99 0 : tpi=2*pimach()
100 :
101 0 : write(oUnit,*)"nkpts=",nkpts
102 : c*****************************************************
103 : c get num_bands and num_wann from the proj file
104 : c*****************************************************
105 0 : do j=jspin,0,-1
106 0 : inquire(file=trim('proj'//spinspin12(j)),exist=l_file)
107 0 : if(l_file)then
108 0 : filename2='proj'//spinspin12(j)
109 0 : exit
110 : endif
111 : enddo
112 0 : if(l_file)then
113 0 : open (203,file=trim(filename2),status='old')
114 0 : rewind (203)
115 : else
116 0 : stop 'no proj/proj.1/proj.2'
117 : endif
118 0 : read (203,*) num_wann,num_bands
119 0 : close (203)
120 0 : write(oUnit,*)'According to proj there are ',num_bands,' bands'
121 0 : write(oUnit,*)"and ",num_wann," wannier functions."
122 :
123 : c****************************************************************
124 : c read in chk
125 : c****************************************************************
126 0 : num_kpts=nkpts
127 0 : allocate( u_matrix_opt(num_bands,num_wann,nkpts) )
128 0 : allocate( u_matrix(num_wann,num_wann,nkpts) )
129 0 : allocate( lwindow(num_bands,nkpts) )
130 0 : allocate( ndimwin(nkpts) )
131 :
132 : call wann_read_umatrix2(
133 : > nkpts,num_wann,num_bands,
134 : > um_format,jspin,wan90version,
135 : < have_disentangled,
136 : < lwindow(:,:),
137 : < ndimwin(:),
138 : < u_matrix_opt(:,:,:),
139 0 : < u_matrix(:,:,:))
140 0 : num_bands2=num_bands
141 :
142 : c****************************************************
143 : c Read the file "filename".
144 : c****************************************************
145 0 : allocate( hsomtx(dim1,dim2,num_bands,num_bands,nkpts) )
146 0 : open(304,file=filename,form='formatted')
147 0 : read(304,*)
148 0 : read(304,*)
149 0 : if(l_conjugate)then
150 0 : do nkp=1,num_kpts
151 0 : do i=1,num_bands
152 0 : do j=1,num_bands
153 0 : do dir2=1,dim2
154 0 : do dir=1,dim1
155 0 : read(304,*)dummy1,dummy2,dummy3,dummy4,dummy5,a,b
156 0 : hsomtx(dir,dir2,j,i,nkp)=cmplx(a,-b)
157 : enddo !dir
158 : enddo !dir2
159 : enddo !j
160 : enddo !i
161 : enddo !nkp
162 : else
163 0 : do nkp=1,num_kpts
164 0 : do i=1,num_bands
165 0 : do j=1,num_bands
166 0 : do dir2=1,dim2
167 0 : do dir=1,dim1
168 0 : read(304,*)dummy1,dummy2,dummy3,dummy4,dummy5,a,b
169 0 : hsomtx(dir,dir2,j,i,nkp)=cmplx(a,b)
170 : enddo !dir
171 : enddo !dir2
172 : enddo !j
173 : enddo !i
174 : enddo !nkp
175 : endif
176 0 : close(304)
177 :
178 : c****************************************************************
179 : c fft5: transformation to Bloch-like
180 : c****************************************************************
181 0 : allocate( hsomtx2(dim1,dim2,num_wann,num_wann,nkpts) )
182 0 : write(oUnit,*)"fft5: transformation to Bloch-like"
183 :
184 0 : if(have_disentangled) then
185 0 : hsomtx2=0.0
186 0 : do nkp=1,num_kpts
187 0 : do j=1,num_wann
188 0 : do jp=1,num_wann
189 0 : do m=1,ndimwin(nkp)
190 0 : do mp=1,ndimwin(nkp)
191 0 : do dir2=1,dim2
192 0 : do dir=1,dim1
193 : hsomtx2(dir,dir2,jp,j,nkp)=hsomtx2(dir,dir2,jp,j,nkp)+
194 : & conjg(u_matrix_opt(mp,jp,nkp))*
195 : & hsomtx(dir,dir2,mp,m,nkp)*
196 0 : & u_matrix_opt(m,j,nkp)
197 : enddo !dir
198 : enddo !dir2
199 : enddo !mp
200 : enddo !m
201 : enddo !jp
202 : enddo !j
203 : enddo !nkp
204 : else
205 0 : hsomtx2 = hsomtx
206 : end if !have_disentangled
207 :
208 0 : allocate(hwann(dim1,dim2,num_wann,num_wann,num_kpts))
209 0 : hwann=cmplx(0.0,0.0)
210 0 : wann_shift=0
211 0 : do k=1,num_kpts
212 0 : print*,"k=",k
213 0 : do m=1,num_wann
214 0 : do mp=1,num_wann
215 0 : do i=1,num_wann
216 0 : do j=1,num_wann
217 0 : do dir2=1,dim2
218 0 : do dir=1,dim1
219 : hwann(dir,dir2,mp,m,k)=hwann(dir,dir2,mp,m,k)+
220 : * conjg(u_matrix(j,mp,k))*
221 : * hsomtx2(dir,dir2,j,i,k)*
222 0 : * u_matrix(i,m,k)
223 : enddo
224 : enddo
225 : enddo !j
226 : enddo !i
227 : enddo !mp
228 : enddo !m
229 : enddo !k
230 :
231 : c************************************************************
232 : c Calculate matrix elements in real space.
233 : c***********************************************************
234 0 : write(oUnit,*)"transform to rs"
235 0 : allocate(hreal(dim1,dim2,num_wann,num_wann,rvecnum))
236 0 : hreal=cmplx(0.0,0.0)
237 0 : do rvecind=1,rvecnum
238 0 : do k=1,nkpts
239 : rdotk=tpi*( kpoints(1,k)*rvec(1,rvecind)+
240 : + kpoints(2,k)*rvec(2,rvecind)+
241 0 : + kpoints(3,k)*rvec(3,rvecind) )
242 0 : fac=cmplx(cos(rdotk),-sin(rdotk))
243 :
244 0 : do m2=1,num_wann
245 0 : do m1=1,num_wann
246 0 : do dir2=1,dim2
247 0 : do dir=1,dim1
248 :
249 : hreal(dir,dir2,m1,m2,rvecind)=
250 : & hreal(dir,dir2,m1,m2,rvecind)+
251 0 : & fac*hwann(dir,dir2,m1,m2,k)
252 : enddo !dir
253 : enddo !dir2
254 : enddo !m1
255 : enddo !m2
256 :
257 : enddo !k
258 : enddo !rvecind
259 0 : hreal=hreal/cmplx(real(nkpts),0.0)
260 :
261 0 : open(321,file=title,form='formatted')
262 0 : do rvecind=1,rvecnum
263 0 : r3=rvec(3,rvecind)
264 0 : r2=rvec(2,rvecind)
265 0 : r1=rvec(1,rvecind)
266 :
267 0 : do j=1,num_wann
268 0 : do i=1,num_wann
269 0 : do dir2=1,dim2
270 0 : do dir=1,dim1
271 : write(321,'(i3,1x,i3,1x,i3,1x,i3,1x,i3,
272 : & 1x,i3,1x,i3,1x,f20.8,1x,f20.8)')
273 0 : & r1,r2,r3,i,j,dir,dir2,
274 0 : & hreal(dir,dir2,i,j,rvecind)
275 : enddo !dir
276 : enddo !dir2
277 : enddo!i
278 : enddo !j
279 :
280 : enddo !rvecnum
281 0 : close(321)
282 :
283 : c deallocate(lwindow,u_matrix_opt,ndimwin)
284 : c deallocate(u_matrix,hwann,hreal)
285 :
286 0 : call timestop("wann_fft5")
287 0 : end subroutine wann_fft5
288 : end module m_wann_fft5
|