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_pauli_rs
8 : USE m_juDFT
9 : contains
10 0 : subroutine wann_pauli_rs(
11 0 : > rvecnum,rvec,kpoints,
12 : > jspins_in,nkpts,l_bzsym,film,l_onedimens,
13 : > l_soc,band_min,band_max,neigd,
14 0 : > l_socmmn0,l_ndegen,ndegen,wan90version,
15 : > l_unformatted)
16 : c*************************************************
17 : c Calculate the matrix elements of the
18 : c Pauli matrix in real space from the
19 : c files WF1.chk (and WF1_um.dat) produced
20 : c by wannier90.
21 : c FF, January 2009
22 : c*************************************************
23 : use m_constants, only:pimach, ImagUnit
24 : use m_wann_read_umatrix
25 :
26 : implicit none
27 : integer, intent(in) :: rvecnum
28 : integer, intent(in) :: rvec(:,:)
29 : real, intent(in) :: kpoints(:,:)
30 : integer, intent(in) :: jspins_in
31 : integer, intent(in) :: nkpts
32 : logical,intent (in) :: l_bzsym,l_soc
33 : logical,intent(in) :: film
34 : integer,intent(in) :: band_min(2),band_max(2),neigd
35 : logical, intent(in) :: l_socmmn0
36 : logical, intent(in) :: l_ndegen
37 : integer, intent(in) :: ndegen(:)
38 : integer, intent(in) :: wan90version
39 : logical, intent(in) :: l_unformatted
40 :
41 : integer :: ikpt,jspins
42 : integer :: kpts
43 : logical :: l_file
44 : integer :: num_wann,num_kpts,num_nnmax,jspin
45 : integer :: kspin,kkspin
46 : integer :: wann_shift,num_wann2
47 : integer :: i,j,k,m,info,r1,r2,r3,dummy1
48 : integer :: dummy2,dummy3
49 : integer :: hopmin,hopmax,counter,m1,m2
50 : integer :: num_bands2
51 : integer,allocatable :: iwork(:)
52 : real,allocatable :: energy(:,:),ei(:)
53 : real,allocatable :: eigw(:,:),rwork(:)
54 : complex,allocatable :: work(:),vec(:,:)
55 0 : complex,allocatable :: u_matrix(:,:,:,:),hwann(:,:,:,:)
56 0 : complex,allocatable :: hreal(:,:,:,:)
57 : complex,allocatable :: hrealsoc(:,:,:,:,:,:,:)
58 : complex,allocatable :: hwannsoc(:,:,:,:,:)
59 0 : complex,allocatable :: paulimat(:,:,:,:)
60 0 : complex,allocatable :: paulimat_opt(:,:,:,:)
61 0 : complex,allocatable :: paulimat2(:,:,:,:)
62 : complex :: fac,eulav,eulav1
63 : real :: tmp_omi,rdotk,tpi,minenerg,maxenerg
64 : real, allocatable :: minieni(:),maxieni(:)
65 : character :: jobz,uplo
66 : integer :: kpt,band,lee,lwork,lrwork,liwork,n,lda
67 : complex :: value(4)
68 : logical :: um_format
69 : logical :: repro_eig
70 : logical :: l_chk,l_proj,l_onedimens
71 : logical :: have_disentangled
72 0 : integer,allocatable :: ndimwin(:,:)
73 0 : logical,allocatable :: lwindow(:,:,:)
74 : integer :: chk_unit,nkp,ntmp,ierr
75 : character(len=33) :: header
76 : character(len=20) :: checkpoint
77 : real :: tmp_latt(3,3), tmp_kpt_latt(3,nkpts)
78 : real :: omega_invariant
79 0 : complex,allocatable :: u_matrix_opt(:,:,:,:)
80 : integer :: num_bands,counter1,counter2
81 : logical :: l_umdat
82 : real,allocatable :: eigval2(:,:)
83 : real,allocatable :: eigval_opt(:,:)
84 : real :: scale,a,b
85 : character(len=2) :: spinspin12(0:2)
86 : character(len=3) :: spin12(2)
87 : character(len=6) :: filename
88 : integer :: jp,mp,kk
89 : integer :: dir,rvecind
90 : integer :: spin1,spin2
91 0 : complex,allocatable :: amn(:,:,:)
92 0 : complex,allocatable :: paulimatmix(:,:,:,:)
93 : data spinspin12/' ','.1' , '.2'/
94 : data spin12/'WF1','WF2'/
95 : integer :: nbnd,nwfs,fullnkpts
96 0 : complex,allocatable :: hwannmix(:,:,:,:)
97 0 : complex,allocatable :: perpmag(:,:,:,:)
98 :
99 0 : tpi=2*pimach()
100 :
101 0 : jspins=jspins_in
102 0 : if(l_soc)jspins=1
103 :
104 0 : write(6,*)"nkpts=",nkpts
105 :
106 : c*****************************************************
107 : c get num_bands and num_wann from the proj file
108 : c*****************************************************
109 0 : do j=jspins,0,-1
110 0 : inquire(file=trim('proj'//spinspin12(j)),exist=l_file)
111 0 : if(l_file)then
112 0 : filename='proj'//spinspin12(j)
113 0 : exit
114 : endif
115 : enddo
116 0 : if(l_file)then
117 0 : open (203,file=trim(filename),status='old')
118 0 : rewind (203)
119 : else
120 : CALL judft_error("no proj/proj.1/proj.2",calledby
121 0 : + ="wann_pauli_rs")
122 : endif
123 0 : read (203,*) num_wann,num_bands
124 0 : close (203)
125 0 : write(6,*)'According to proj there are ',num_bands,' bands'
126 0 : write(6,*)"and ",num_wann," wannier functions."
127 :
128 0 : allocate( u_matrix_opt(num_bands,num_wann,nkpts,jspins) )
129 0 : allocate( u_matrix(num_wann,num_wann,nkpts,jspins) )
130 0 : allocate( lwindow(num_bands,nkpts,jspins) )
131 0 : allocate( ndimwin(nkpts,jspins) )
132 :
133 0 : do jspin=1,jspins !spin loop
134 :
135 :
136 : c****************************************************************
137 : c read in chk
138 : c****************************************************************
139 0 : num_kpts=nkpts
140 :
141 : call wann_read_umatrix2(
142 : > nkpts,num_wann,num_bands,
143 : > um_format,jspin,wan90version,
144 : < have_disentangled,
145 : < lwindow(:,:,jspin),
146 : < ndimwin(:,jspin),
147 : < u_matrix_opt(:,:,:,jspin),
148 0 : < u_matrix(:,:,:,jspin))
149 :
150 :
151 : enddo
152 :
153 0 : jspin=1
154 :
155 0 : num_bands2=jspins*num_bands
156 0 : num_wann2=jspins*num_wann
157 :
158 0 : allocate( paulimat(3,num_bands2,num_bands2,nkpts) )
159 0 : paulimat=0.0
160 :
161 : c****************************************************
162 : c Read the file "WF1.socmmn0".
163 : c****************************************************
164 0 : if(l_unformatted)then
165 0 : if(l_soc)then
166 0 : open(304,file='WF1.socmmn0_unf',form='unformatted')
167 : else
168 0 : open(304,file='WF1.mmn0_unf',form='unformatted')
169 : endif
170 0 : read(304)nbnd,fullnkpts,nwfs
171 0 : allocate(amn(nbnd,nwfs,fullnkpts))
172 0 : read(304)amn(1:nbnd,1:nwfs,1:fullnkpts)
173 0 : do nkp=1,num_kpts
174 0 : do i=1,num_bands
175 0 : do j=1,num_bands
176 0 : paulimat(3,j,i,nkp)=amn(i,j,nkp)
177 : enddo !j
178 : enddo !i
179 : enddo !nkp
180 0 : close(304)
181 0 : deallocate(amn)
182 : else !l_unformatted
183 0 : if(l_soc)then
184 0 : open(304,file='WF1.socmmn0',form='formatted')
185 : else
186 0 : open(304,file='WF1.mmn0',form='formatted')
187 : endif
188 0 : read(304,*)
189 0 : read(304,*)
190 0 : do nkp=1,num_kpts
191 0 : do i=1,num_bands
192 0 : do j=1,num_bands
193 0 : read(304,*)dummy1,dummy2,dummy3,a,b
194 0 : paulimat(3,i,j,nkp)=cmplx(a,b)
195 : enddo !j
196 : enddo !i
197 : enddo !nkp
198 0 : close(304)
199 : endif
200 : c****************************************************
201 : c Read the file "WF2.socmmn0".
202 : c****************************************************
203 0 : if(l_unformatted)then
204 0 : if(l_soc)then
205 0 : open(304,file='WF2.socmmn0_unf',form='unformatted')
206 : else
207 0 : open(304,file='WF2.mmn0_unf',form='unformatted')
208 : endif
209 0 : read(304)nbnd,fullnkpts,nwfs
210 0 : allocate(amn(nbnd,nwfs,fullnkpts))
211 0 : read(304)amn(1:nbnd,1:nwfs,1:fullnkpts)
212 0 : do nkp=1,num_kpts
213 0 : do i=1+(jspins-1)*num_bands,num_bands2
214 0 : do j=1+(jspins-1)*num_bands,num_bands2
215 : paulimat(3,i,j,nkp)=
216 : & paulimat(3,i,j,nkp)-
217 : & amn(j-(jspins-1)*num_bands,
218 0 : & i-(jspins-1)*num_bands,nkp)
219 : enddo !j
220 : enddo !i
221 : enddo !nkp
222 0 : close(304)
223 0 : deallocate(amn)
224 : else !l_unformatted
225 0 : if(l_soc)then
226 0 : open(304,file='WF2.socmmn0',form='formatted')
227 : else
228 0 : open(304,file='WF2.mmn0',form='formatted')
229 : endif
230 0 : read(304,*)
231 0 : read(304,*)
232 0 : do nkp=1,num_kpts
233 0 : do i=1+(jspins-1)*num_bands,num_bands2
234 0 : do j=1+(jspins-1)*num_bands,num_bands2
235 0 : read(304,*)dummy1,dummy2,dummy3,a,b
236 : paulimat(3,i,j,nkp)=
237 0 : & paulimat(3,i,j,nkp)-cmplx(a,b)
238 : enddo !j
239 : enddo !i
240 : enddo !nkp
241 0 : close(304)
242 : endif !l_unformatted
243 :
244 : c****************************************************
245 : c Read the file "updown.mmn0".
246 : c****************************************************
247 0 : if(l_unformatted)then
248 0 : open(304,file='updown.mmn0_unf',form='unformatted')
249 0 : read(304)nbnd,fullnkpts,nwfs
250 0 : allocate(amn(nbnd,nwfs,fullnkpts))
251 0 : read(304)amn(1:nbnd,1:nwfs,1:fullnkpts)
252 0 : close(304)
253 0 : deallocate(amn)
254 0 : do nkp=1,fullnkpts
255 0 : do i=1,num_bands
256 0 : do j=1,num_bands
257 : paulimat(1,j,i+num_bands*(jspins-1),nkp)=
258 0 : & conjg(amn(j,i,nkp))
259 : enddo !j
260 : enddo !i
261 0 : paulimat(2,:,:,nkp)= -ImagUnit*paulimat(1,:,:,nkp)
262 : paulimat(1,:,:,nkp)=paulimat(1,:,:,nkp)+
263 0 : & transpose(conjg( paulimat(1,:,:,nkp) ))
264 : paulimat(2,:,:,nkp)=paulimat(2,:,:,nkp)+
265 0 : & transpose(conjg( paulimat(2,:,:,nkp) ))
266 : enddo !nkp
267 :
268 : else !l_unformatted
269 0 : open(304,file='updown.mmn0',form='formatted')
270 0 : read(304,*)
271 0 : read(304,*)
272 0 : do nkp=1,num_kpts
273 0 : do i=1,num_bands
274 0 : do j=1,num_bands
275 0 : read(304,*)dummy1,dummy2,dummy3,a,b
276 0 : paulimat(1,j,i+num_bands*(jspins-1),nkp)=cmplx(a,-b)
277 : enddo !j
278 : enddo !i
279 0 : paulimat(2,:,:,nkp)= -ImagUnit*paulimat(1,:,:,nkp)
280 : paulimat(1,:,:,nkp)=paulimat(1,:,:,nkp)+
281 0 : & transpose(conjg( paulimat(1,:,:,nkp) ))
282 : paulimat(2,:,:,nkp)=paulimat(2,:,:,nkp)+
283 0 : & transpose(conjg( paulimat(2,:,:,nkp) ))
284 : enddo !nkp
285 0 : close(304)
286 : endif !l_unformatted
287 :
288 : c****************************************************************
289 : c Calculate matrix elements of Pauli in the basis of
290 : c rotated Bloch functions.
291 : c****************************************************************
292 0 : allocate( paulimat2(3,num_wann2,num_wann2,nkpts) )
293 : write(6,*)"calculate matrix elements of spin operator
294 0 : & between wannier orbitals"
295 :
296 0 : if(have_disentangled) then
297 :
298 0 : allocate( paulimat_opt(3,num_bands2,num_bands2,nkpts) )
299 0 : allocate( paulimatmix(3,num_wann2,num_bands2,nkpts))
300 :
301 0 : do nkp=1,num_kpts
302 : counter1=0
303 0 : do m=1,num_bands2
304 0 : spin1=(m-1)/num_bands
305 0 : if(lwindow(m-spin1*num_bands,nkp,spin1+1))then
306 0 : counter1=counter1+1
307 0 : counter2=0
308 0 : do mp=1,num_bands2
309 0 : spin2=(mp-1)/num_bands
310 0 : if(lwindow(mp-spin2*num_bands,nkp,spin2+1))then
311 0 : counter2=counter2+1
312 0 : do k=1,3
313 : paulimat_opt(k,counter2,counter1,nkp)=
314 0 : & paulimat(k,mp,m,nkp)
315 : enddo
316 : endif
317 : enddo !mp
318 : endif
319 : enddo !m
320 : enddo
321 :
322 0 : paulimatmix=0.0
323 0 : do spin1=0,jspins-1
324 0 : do spin2=0,jspins-1
325 0 : do nkp=1,num_kpts
326 0 : do jp=1,num_wann
327 0 : do m=1,ndimwin(nkp,spin2+1)
328 0 : do mp=1,ndimwin(nkp,spin1+1)
329 0 : do k=1,3
330 :
331 : paulimatmix(k,jp+spin1*num_wann,m+spin2*ndimwin(nkp,1),nkp)=
332 : = paulimatmix(k,jp+spin1*num_wann,m+spin2*ndimwin(nkp,1),nkp)+
333 : & conjg(u_matrix_opt(mp,jp,nkp,spin1+1))*
334 : & paulimat_opt(k,mp+spin1*ndimwin(nkp,1),
335 0 : & m+spin2*ndimwin(nkp,1),nkp)
336 : enddo !k
337 : enddo !mp
338 : enddo !m
339 : enddo !jp
340 : enddo !nkp
341 : enddo !spin2
342 : enddo !spin1
343 :
344 0 : paulimat2=0.0
345 0 : do spin1=0,jspins-1
346 0 : do spin2=0,jspins-1
347 0 : do nkp=1,num_kpts
348 0 : do j=1,num_wann
349 0 : do jp=1,num_wann
350 0 : do m=1,ndimwin(nkp,spin2+1)
351 :
352 0 : do k=1,3
353 :
354 :
355 : paulimat2(k,jp+spin1*num_wann,j+spin2*num_wann,nkp)=
356 : = paulimat2(k,jp+spin1*num_wann,j+spin2*num_wann,nkp)+
357 : & paulimatmix(k,jp+spin1*num_wann,
358 :
359 : & m+spin2*ndimwin(nkp,1),nkp)*
360 0 : & u_matrix_opt(m,j,nkp,spin2+1)
361 : enddo !k
362 : enddo !m
363 : enddo !jp
364 : enddo !j
365 : enddo !nkp
366 : enddo !spin2
367 : enddo !spin1
368 :
369 :
370 0 : deallocate(paulimatmix)
371 0 : deallocate(paulimat_opt)
372 :
373 : else
374 0 : paulimat2=paulimat
375 : end if !have_disentangled
376 :
377 0 : allocate(hwann(3,num_wann2,num_wann2,num_kpts))
378 0 : hwann=cmplx(0.0,0.0)
379 0 : wann_shift=0
380 :
381 0 : allocate(hwannmix(3,num_wann2,num_wann2,num_kpts))
382 0 : hwannmix=cmplx(0.0,0.0)
383 :
384 0 : do spin1=0,jspins-1
385 0 : do spin2=0,jspins-1
386 0 : do k=1,num_kpts
387 0 : do i=1,num_wann
388 0 : do mp=1,num_wann
389 0 : do j=1,num_wann
390 :
391 0 : do kk=1,3
392 : hwannmix(kk,mp+spin1*num_wann,i+spin2*num_wann,k)=
393 : = hwannmix(kk,mp+spin1*num_wann,i+spin2*num_wann,k)+
394 : * conjg(u_matrix(j,mp,k,spin1+1))*
395 : * paulimat2(kk,j+spin1*num_wann,
396 0 : * i+spin2*num_wann,k)
397 :
398 : enddo !kk
399 : enddo !j
400 : enddo !mp
401 :
402 : enddo !i
403 : enddo !k
404 : enddo !spin2
405 : enddo !spin1
406 :
407 0 : do spin1=0,jspins-1
408 0 : do spin2=0,jspins-1
409 0 : do k=1,num_kpts
410 0 : do m=1,num_wann
411 0 : do i=1,num_wann
412 0 : do mp=1,num_wann
413 0 : do kk=1,3
414 : hwann(kk,mp+spin1*num_wann,m+spin2*num_wann,k)=
415 : = hwann(kk,mp+spin1*num_wann,m+spin2*num_wann,k)+
416 : * hwannmix(kk,mp+spin1*num_wann,
417 : * i+spin2*num_wann,k)*
418 0 : * u_matrix(i,m,k,spin2+1)
419 : enddo !kk
420 : enddo !mp
421 : enddo !i
422 : enddo !m
423 : enddo !k
424 : enddo !spin2
425 : enddo !spin1
426 :
427 0 : deallocate(hwannmix)
428 :
429 :
430 :
431 : c************************************************************
432 : c Calculate matrix elements in real space.
433 : c***********************************************************
434 0 : write(6,*)"calculate pauli-mat in rs"
435 0 : allocate(hreal(3,num_wann2,num_wann2,rvecnum))
436 0 : hreal = cmplx(0.0,0.0)
437 0 : if(l_ndegen)then
438 0 : do rvecind=1,rvecnum
439 0 : do k=1,nkpts
440 : rdotk=tpi*( kpoints(1,k)*rvec(1,rvecind)+
441 : & kpoints(2,k)*rvec(2,rvecind)+
442 0 : & kpoints(3,k)*rvec(3,rvecind) )
443 0 : fac=cmplx(cos(rdotk),-sin(rdotk))/real(ndegen(rvecind))
444 0 : do m2=1,num_wann2
445 0 : do m1=1,num_wann2
446 0 : do dir=1,3
447 : hreal(dir,m1,m2,rvecind)=
448 : & hreal(dir,m1,m2,rvecind)+
449 0 : & fac*hwann(dir,m1,m2,k)
450 : enddo !dir
451 : enddo !m1
452 : enddo !m2
453 : enddo !k
454 : enddo !rvecind
455 : else
456 0 : do rvecind=1,rvecnum
457 0 : do k=1,nkpts
458 : rdotk=tpi*( kpoints(1,k)*rvec(1,rvecind)+
459 : & kpoints(2,k)*rvec(2,rvecind)+
460 0 : & kpoints(3,k)*rvec(3,rvecind) )
461 0 : fac=cmplx(cos(rdotk),-sin(rdotk))
462 0 : do m2=1,num_wann2
463 0 : do m1=1,num_wann2
464 0 : do dir=1,3
465 : hreal(dir,m1,m2,rvecind)=
466 : & hreal(dir,m1,m2,rvecind)+
467 0 : & fac*hwann(dir,m1,m2,k)
468 : enddo !dir
469 : enddo !m1
470 : enddo !m2
471 : enddo !k
472 : enddo !rvecind
473 : endif !l_ndegen
474 0 : hreal=hreal/cmplx(real(nkpts),0.0)
475 :
476 :
477 0 : if(l_unformatted)then
478 0 : allocate(perpmag(num_wann2,num_wann2,3,rvecnum))
479 :
480 0 : if(l_ndegen)then
481 : open(321,file='rspauli_ndegen_unf'//spinspin12(1),
482 : & form='unformatted'
483 : #ifdef CPP_INTEL
484 : & ,convert='BIG_ENDIAN'
485 : #endif
486 0 : & )
487 : else
488 : open(321,file='rspauli_unf'//spinspin12(1),
489 : & form='unformatted'
490 : #ifdef CPP_INTEL
491 : & ,convert='BIG_ENDIAN'
492 : #endif
493 0 : & )
494 : endif
495 :
496 0 : do rvecind=1,rvecnum
497 0 : do j=1,num_wann2
498 0 : do i=1,num_wann2
499 0 : do kk=1,3
500 0 : perpmag(i,j,kk,rvecind)=hreal(kk,i,j,rvecind)
501 : enddo !kk
502 : enddo !i
503 : enddo !j
504 : enddo !rvecind
505 0 : write(321)perpmag
506 0 : deallocate(perpmag)
507 0 : close(321)
508 :
509 : else !l_unformatted
510 :
511 0 : if(l_ndegen)then
512 : open(321,file='rspauli_ndegen'//spinspin12(1),
513 0 : & form='formatted')
514 : else
515 : open(321,file='rspauli'//spinspin12(1),
516 0 : & form='formatted')
517 : endif
518 0 : do rvecind=1,rvecnum
519 0 : r3=rvec(3,rvecind)
520 0 : r2=rvec(2,rvecind)
521 0 : r1=rvec(1,rvecind)
522 0 : do j=1,num_wann2
523 0 : do i=1,num_wann2
524 0 : do kk=1,3
525 : write(321,'(i3,1x,i3,1x,i3,1x,i3,
526 : & 1x,i3,1x,i3,1x,f20.8,1x,f20.8)')
527 0 : & r1,r2,r3,i,j,kk,hreal(kk,i,j,rvecind)
528 : enddo !kk
529 : enddo !i
530 : enddo !j
531 : enddo !rvecnum
532 :
533 0 : close(321)
534 :
535 : endif
536 :
537 0 : deallocate(lwindow,u_matrix_opt,ndimwin)
538 0 : deallocate(u_matrix,hwann,hreal)
539 :
540 0 : end subroutine wann_pauli_rs
541 : end module m_wann_pauli_rs
|
