Branch data Line data Source code
1 : : !!****m* ABINIT/m_wannier_io
2 : : !! NAME
3 : : !! m_wannier_io
4 : : !!
5 : : !! FUNCTION
6 : : !! subroutines for writting Wannier90 related files.
7 : : !!
8 : : !! COPYRIGHT
9 : : !! Copyright (C) 2005-2026 ABINIT group (BAmadon, CEspejo, FJollet, TRangel, DRH, hexu)
10 : : !! This file is distributed under the terms of the
11 : : !! GNU General Public License, see ~abinit/COPYING
12 : : !! or http://www.gnu.org/copyleft/gpl.txt .
13 : : !!
14 : : !! SOURCE
15 : :
16 : : #if defined HAVE_CONFIG_H
17 : : #include "config.h"
18 : : #endif
19 : :
20 : : #include "abi_common.h"
21 : :
22 : : module m_wannier_io
23 : :
24 : : use defs_basis
25 : : use defs_wannier90
26 : : use m_abicore
27 : : use m_errors
28 : : use m_atomdata
29 : : use m_xmpi
30 : : use m_sort
31 : : use m_dtfil
32 : : #ifdef FC_NAG
33 : : use f90_unix_dir
34 : : #endif
35 : :
36 : : use m_dtset, only : dataset_type
37 : : use defs_abitypes, only : MPI_type
38 : : use m_fft, only: fourwf
39 : : use m_io_tools, only : open_file, get_unit
40 : : use m_fftcore, only : sphereboundary
41 : : use m_abstract_wf, only: abstract_wf, cg_cprj, wfd_wf
42 : : use m_ebands, only : ebands_t
43 : :
44 : : implicit none
45 : :
46 : : private
47 : :
48 : : public :: write_Amn
49 : : public :: compute_and_write_unk
50 : : public :: write_eigenvalues
51 : : public :: write_Mmn
52 : : public :: read_chkunit
53 : : !!***
54 : :
55 : : contains
56 : :
57 : : ! Write amn file
58 [ - + ]: 19 : subroutine write_Amn(A_matrix, fname, nsppol, mband, nkpt, num_bands, nwan, band_in, test_matrix_output)
59 : : ! TODO use the A_matrix sizes instead of the nsppol, mband, nkpt, nwan
60 : : complex(dp),pointer :: A_matrix(:,:,:,:)
61 : : !type(dataset_type),intent(in) :: dtset
62 : : logical, intent(in) :: band_in(:, :), test_matrix_output
63 : : integer, intent(in) :: nsppol, num_bands(nsppol), nwan(nsppol), mband, nkpt
64 : : character(len=fnlen), intent(in) :: fname(nsppol)
65 : :
66 : : !Local variables-------------------------------
67 : : !scalars
68 : : character(len=1000) :: msg
69 : 19 : integer :: iun(nsppol), units(2)
70 : : integer :: isppol, ikpt, iband, iwan, jband, ii, jj
71 : :
72 [ + + ]: 57 : units = [std_out, ab_out]
73 : :
74 : : ! below is copied/modified from m_mlwfovlp.F90
75 [ + + ]: 39 : do isppol=1,nsppol
76 : : ! TODO : relpace this.
77 [ + - - + ]: 20 : if (open_file(trim(fname(isppol)), msg, newunit=iun(isppol), form="formatted", status="unknown", action="write") /= 0) then
78 : 0 : ABI_ERROR(msg)
79 : : end if
80 : 20 : write(iun(isppol),*) 'Projections from Abinit : mband,nkpt,nwan. indices: iband1,iwan,ikpt'
81 : 39 : write(iun(isppol),*) num_bands(isppol),nkpt,nwan(isppol)
82 : : end do
83 : :
84 [ + + ]: 39 : do isppol=1,nsppol
85 [ + + ]: 743 : do ikpt=1,nkpt
86 [ + + ]: 5364 : do iwan=1,nwan(isppol)
87 : 4640 : jband=0
88 [ + + ]: 72608 : do iband=1,mband
89 [ + + ]: 71904 : if(band_in(iband,isppol)) then
90 : 59648 : jband=jband+1
91 : 59648 : write(iun(isppol),'(3i6,13x,3x,2f18.14)')jband,iwan,ikpt,A_matrix(jband,iwan,ikpt,isppol)
92 : : end if !band_in
93 : : end do !iband
94 : : end do !iwan
95 : : end do !ikpt
96 : : end do !isppol
97 : : !
98 [ + + ]: 39 : do isppol=1,nsppol
99 : 20 : close(iun(isppol))
100 [ + - ]: 20 : write(msg, '(3a)' )' ',trim(fname(isppol)),' written'
101 : 39 : call wrtout(std_out,msg)
102 : : end do
103 : :
104 : : ! Individual matrix elements are gauge-dependent for spinor wavefunctions,
105 : : ! especially inside degenerate subspaces. The caller disables this test-only
106 : : ! output for nspinor == 2, while the complete .amn file is still written.
107 [ + + ]: 19 : if (.not. test_matrix_output) return
108 : : !
109 : : ! Write down part of the matrix to the output file
110 : : ! This is for the automatic tests
111 : : !
112 : 17 : write(msg, '(4a)' ) ch10,&
113 : 17 : ' Writing top of the initial projections matrix: A_mn(ik)',ch10,&
114 : 34 : ' m=1:3, n=1:3, ik=1'
115 : 17 : call wrtout(units, msg)
116 : : !
117 : : ! just write down the first 3 elements
118 : : !
119 [ + + ]: 35 : do isppol=1,nsppol
120 : 18 : write(msg, '( " " )')
121 [ + + ]: 18 : if (nsppol>1 ) then
122 [ + + ]: 2 : if (isppol==1) write(msg,'(2a)')trim(msg),' spin up:'
123 [ + + ]: 2 : if (isppol==2) write(msg,'(2a)')trim(msg),' spin down:'
124 : : end if
125 [ + + ]: 72 : do ii=1,3
126 [ - + ]: 54 : if(ii>num_bands(isppol)) cycle
127 : 54 : write(msg,'(3a)') trim(msg),ch10,'; ( '
128 [ + + ]: 216 : do jj=1,3
129 [ - + ]: 162 : if(jj>nwan(isppol))cycle
130 : 162 : write(msg, '(a,2f11.6,a)') trim(msg),&
131 : 378 : & A_matrix(ii,jj,1,isppol),' , '
132 : : end do
133 : 72 : write(msg,'(2a)') trim(msg),' ) '
134 : : end do
135 : 35 : call wrtout(units, msg)
136 : : end do
137 : : !
138 : : ! Now write down bottom of the matrix
139 : : !
140 : 17 : write(msg, '(4a)' ) ch10,&
141 : 17 : & ' Writing bottom of the initial projections matrix: A_mn(ik)',ch10,&
142 : 34 : & ' m=num_bands-2:num_bands, n=nwan-2:nwan, ik=nkpt'
143 : 17 : call wrtout(units, msg)
144 : :
145 [ + + ]: 35 : do isppol=1,nsppol
146 : 18 : write(msg, '( " " )')
147 [ + + ]: 18 : if (nsppol>1 ) then
148 [ + + ]: 2 : if (isppol==1) write(msg,'(2a)')trim(msg),' spin up:'
149 [ + + ]: 2 : if (isppol==2) write(msg,'(2a)')trim(msg),' spin down:'
150 : : end if
151 [ + + ]: 72 : do ii=num_bands(isppol)-2,num_bands(isppol)
152 [ - + ]: 54 : if(ii<1) cycle
153 : 54 : write(msg,'(3a)') trim(msg),ch10,'; ( '
154 [ + + ]: 216 : do jj=nwan(isppol)-2,nwan(isppol)
155 [ - + ]: 162 : if(jj<1)cycle
156 : 162 : write(msg, '(a,2f11.6,a)') trim(msg),&
157 : 378 : & A_matrix(ii,jj,nkpt,isppol),' , '
158 : : end do
159 : 72 : write(msg,'(2a)') trim(msg),' ) '
160 : : end do
161 : 35 : call wrtout(units, msg)
162 : : end do !isppol
163 : :
164 : : end subroutine write_Amn
165 : :
166 : :
167 : : !-----------------------------------------------------------------------
168 : : !@brief Write down the unk matrix
169 : : !-----------------------------------------------------------------------
170 : 4 : subroutine compute_and_write_unk(wfnname, usepaw, w90prtunk, &
171 : : & mpi_enreg, ngfft, nsppol, nspinor, &
172 [ - + ]: 4 : & nkpt, mband, mpw, mgfftc, mkmem, nprocs, rank, npwarr, &
173 [ - + ]: 4 : & band_in, dtset, kg, mywfc)
174 : : !TODO split the calculation of unk with writting.
175 : : character(len=fnlen), intent(inout) :: wfnname
176 : : integer ,intent(in) :: usepaw
177 : : integer, intent(in) ::w90prtunk
178 : : type(MPI_type),intent(in) :: mpi_enreg
179 : : integer, intent(in) :: ngfft(18), nsppol, nspinor, nkpt, mpw, mgfftc, mband, mkmem
180 : : integer, intent(in) :: nprocs, rank, npwarr(:)
181 : : logical, intent(in) :: band_in(:, :)
182 : : type(dataset_type),intent(in) :: dtset
183 : :
184 : : !integer, intent(in) :: iwav(:, :,:,:),kg(3,mpw*mkmem)
185 : : integer, intent(in) :: kg(3,mpw*mkmem)
186 : : !real(dp),intent(in) :: cg(2,mpw*nspinor*mband*mkmem*nsppol)
187 : : class(abstract_wf), intent(inout) :: mywfc
188 : :
189 : : !Local variables-------------------------------
190 : : !scalars
191 : : integer :: iun_plot
192 : : integer :: isppol, ikpt, ikg, iband, ig
193 : : integer :: n1, n2, n3, n4, n5, n6, cplex, mgfft, npw_k
194 : : integer :: tim_fourwf
195 : : real(dp) :: weight
196 : 8 : integer :: spacing, nband_inc(nsppol)
197 : 4 : integer,allocatable:: kg_k(:,:)
198 : 4 : real(dp),allocatable :: denpot(:,:,:), cwavef(:,:), fofgout(:,:),fofr(:,:,:,:)
199 : 4 : integer,allocatable :: gbound(:,:)
200 : : integer :: n1tmp, n2tmp, n3tmp, jj1, jj2, jj3, ipw, ispinor
201 : : character(len=1000) :: msg
202 : :
203 [ + + ]: 4 : if(usepaw==1) then
204 : 2 : write(msg, '( a,a,a,a,a,a,a,a,a)')ch10,&
205 : 2 : & " WARNING: The UNK matrices will not contain the correct wavefunctions ",ch10,&
206 : 2 : & " since we are just writing the plane wave contribution.",ch10,&
207 : 2 : & " The contribution from inside the spheres is missing. ",ch10,&
208 : 4 : & " However, these files can be used for plotting purposes",ch10
209 : 2 : call wrtout(std_out, msg)
210 : : end if
211 : : !
212 : 4 : spacing = w90prtunk
213 : 4 : write(msg, '( 8a,i3,2a)')ch10,&
214 : 4 : & " UNK files will be written.",ch10,&
215 : 4 : & " According to the chosen value of w90prtunk",ch10,&
216 : 4 : & " the wavefunctions are to be written ",ch10, &
217 : 8 : & " at every ", spacing," records.",ch10
218 : 4 : call wrtout(std_out, msg)
219 : : !
220 [ - + + - : 12 : ABI_MALLOC(kg_k,(3,mpw))
+ - + - -
+ - + ]
221 : 4 : n1=ngfft(1)
222 : 4 : n2=ngfft(2)
223 : 4 : n3=ngfft(3)
224 : 4 : n4=ngfft(4)
225 : 4 : n5=ngfft(5)
226 : 4 : n6=ngfft(6)
227 : 4 : cplex=1
228 : 4 : mgfft=mgfftc ! error
229 [ + + ]: 8 : do isppol=1,nsppol
230 : 4 : ikg=0
231 [ + + ]: 96 : do ikpt=1,nkpt
232 : : !
233 : : ! MPI:cycle over k-points not treated by this node
234 : : !
235 [ - + ]: 88 : if (nprocs>1 ) then !sometimes we can have just one processor
236 [ # # ]: 0 : if ( ABS(MPI_enreg%proc_distrb(ikpt,1,isppol)-rank) /=0) CYCLE
237 : : end if
238 : : !
239 : 88 : npw_k=npwarr(ikpt)
240 [ + + + + ]: 109456 : kg_k(:,1:npw_k)=kg(:,1+ikg:npw_k+ikg)
241 [ - + + - : 440 : ABI_MALLOC(denpot,(cplex*n4,n5,n6))
- + + - -
+ + - + -
+ - - + -
+ ]
242 [ - + + - : 264 : ABI_MALLOC(cwavef,(2,npw_k))
+ - + - -
+ - + ]
243 [ - + + - : 440 : ABI_MALLOC(fofr,(2,n4,n5,n6))
- + + - -
+ + - + -
+ - - + -
+ ]
244 [ + - - + : 352 : ABI_MALLOC(gbound,(2*mgfft+8,2))
+ - + - +
- + - - +
- + ]
245 [ - + + - : 176 : ABI_MALLOC(fofgout,(2,npw_k))
- + - + ]
246 : :
247 : : !iun_plot=1000+ikpt+ikpt*(isppol-1)
248 [ + - ]: 88 : if (nspinor == 1) then
249 : 88 : write(wfnname,'("UNK",I5.5,".",I1)') ikpt, isppol
250 [ # # ]: 0 : else if (nspinor == 2) then
251 : 0 : write(wfnname,'("UNK",I5.5,".NC")') ikpt
252 : : end if
253 [ - + ]: 88 : if (open_file(trim(wfnname), msg, newunit=iun_plot, &
254 : : form="unformatted", status="unknown", action="write") /= 0) then
255 : 0 : ABI_ERROR(msg)
256 : : endif
257 : :
258 : : ! open (unit=iun_plot, file=wfnname,form='formatted')
259 : : ! open(unit=iun_plot, file=wfnname,form='unformatted')
260 : : ! optimizing grid for UNK files
261 : 88 : n1tmp = n1/spacing
262 : 88 : n2tmp = n2/spacing
263 : 88 : n3tmp = n3/spacing
264 [ - + ]: 88 : if( mod(n1,spacing) /= 0) then
265 : 0 : n1tmp = n1tmp + 1
266 : : end if
267 [ - + ]: 88 : if( mod(n2,spacing) /= 0) then
268 : 0 : n2tmp = n2tmp + 1
269 : : end if
270 [ - + ]: 88 : if( mod(n3,spacing) /= 0) then
271 : 0 : n3tmp = n3tmp + 1
272 : : end if
273 : : ! write(iun_plot,*) n1tmp,n2tmp,n3tmp,ikpt,nband_inc
274 : 88 : write(iun_plot) n1tmp,n2tmp,n3tmp,ikpt,nband_inc(isppol)
275 : : ! gbound=zero
276 : 88 : call sphereboundary(gbound,mywfc%hdr%istwfk(ikpt),kg_k,mgfft,npw_k)
277 : 88 : write(std_out,*) " writes UNK file for ikpt, spin=",ikpt,isppol
278 [ + + + + : 2506536 : denpot(:,:,:)=zero
+ + ]
279 : 88 : weight = one
280 [ + + ]: 1080 : do iband=1,mband
281 [ + - ]: 1080 : if(band_in(iband,isppol)) then
282 : : ! TODO: check if this is the right order.
283 [ + + ]: 1984 : do ispinor = 1, nspinor
284 [ + + + + ]: 606806 : cwavef = zero
285 [ + + ]: 202930 : do ipw = 1, npw_k
286 : : !do ig=1,npw_k*dtset%nspinor
287 : : !cwavef(1,ig)=cg(1,ipw+iwav(ispinor, iband,ikpt,isppol))
288 : : !cwavef(2,ig)=cg(2,ipw+iwav(ispinor, iband,ikpt,isppol))
289 : : !ig = ipw + (ispinor-1)*npw_k
290 : 201938 : ig = ipw
291 : 201938 : cwavef(1,ig)=mywfc%cg_elem(1,ipw, ispinor, iband, ikpt, isppol)
292 : 202930 : cwavef(2,ig)=mywfc%cg_elem(2,ipw, ispinor, iband, ikpt, isppol)
293 : : end do
294 : 992 : tim_fourwf=0
295 : : call fourwf(cplex,denpot,cwavef,fofgout,fofr,&
296 : : & gbound,gbound,mywfc%hdr%istwfk(ikpt),kg_k,kg_k,mgfft,&
297 : : & mpi_enreg,1,ngfft,npw_k,npw_k,n4,n5,n6,0,&
298 : 992 : & tim_fourwf,weight,weight,gpu_option=dtset%gpu_option)
299 : : ! do jj3=1,n3,spacing
300 : : ! do jj2=1,n2,spacing
301 : : ! do jj1=1,n1,spacing
302 : : ! write(iun_plot,*) fofr(1,jj1,jj2,jj3), fofr(2,jj1,jj2,jj3)
303 : : ! end do !jj1
304 : : ! end do !jj2
305 : : ! end do !jj3
306 : : ! unformatted (must be one record)
307 : 27241344 : write(iun_plot) (((fofr(1,jj1,jj2,jj3),fofr(2,jj1,jj2,jj3),&
308 [ - + - - : 27273824 : jj1=1,n1,spacing),jj2=1,n2,spacing),jj3=1,n3,spacing)
+ - - + -
- + - - +
- - - + -
+ + + + -
+ + + - +
+ ]
309 : : end do !ispinor
310 : : end if !iband
311 : : end do ! iband
312 [ - + ]: 88 : ABI_FREE(cwavef)
313 [ - + ]: 88 : ABI_FREE(fofr)
314 [ - + ]: 88 : ABI_FREE(gbound)
315 [ - + ]: 88 : ABI_FREE(denpot)
316 [ - + ]: 88 : ABI_FREE(fofgout)
317 : 88 : ikg=ikg+npw_k
318 : 92 : close(iun_plot)
319 : : end do ! ikpt
320 : : end do ! nsppol
321 [ - + ]: 4 : ABI_FREE(kg_k)
322 : : !
323 : 4 : write(msg, '(4a)' )ch10, ' ','UNK files written',ch10
324 : 4 : call wrtout(std_out, msg)
325 : 4 : end subroutine compute_and_write_unk
326 : :
327 : : !-----------------------------------------------------------------------
328 : : ! This subroutine writes the eigenvalues in the w90 format
329 : : ! (see w90io.F90 in w90 package)
330 : : !-----------------------------------------------------------------------
331 [ - + - + ]: 19 : subroutine write_eigenvalues(filew90_eig,eigen, band_in, eigenvalues_w, &
332 : : & nsppol, nkpt, mband, dtset, rank, master )
333 : :
334 : : integer, intent(in) :: nsppol, nkpt, mband, rank, master
335 : : character(len=fnlen),intent(in) :: filew90_eig(nsppol)
336 : : logical, intent(in) :: band_in(:, :)
337 : : real(dp), intent(in) :: eigen(mband,nkpt,nsppol)
338 : : real(dp), intent(inout) :: eigenvalues_w(:, :, :)
339 : : type(dataset_type),intent(in) :: dtset
340 : :
341 : : !Local variables-------------------------------
342 : : !scalars
343 : 38 : integer :: iun(nsppol), isppol, band_index, iband, jband, nband_k, ikpt
344 : : character(len=1000) :: msg
345 : : ! Assign file unit numbers
346 [ + - ]: 19 : if(rank==master) then
347 [ + + ]: 39 : do isppol=1,nsppol
348 [ + - - + ]: 20 : if (open_file(trim(filew90_eig(isppol)), msg, newunit=iun(isppol), &
349 : 19 : form="formatted", status="unknown", action="write") /= 0) then
350 : 0 : ABI_ERROR(msg)
351 : : endif
352 : : end do
353 : : end if !rank==master
354 : : ! Loop to write eigenvalues
355 : : band_index=0
356 [ + + ]: 39 : do isppol=1,nsppol
357 [ + + ]: 743 : do ikpt=1,nkpt
358 : 704 : nband_k=dtset%nband(ikpt+(isppol-1)*nkpt)
359 : 704 : jband=0
360 [ + + ]: 9232 : do iband=1,mband
361 [ + + ]: 9232 : if(band_in(iband,isppol)) then
362 : 6656 : jband=jband+1
363 : : ! Writing data
364 [ + - ]: 6656 : if (rank==master) then
365 : 6656 : write(iun(isppol), '(2i6,4x,f10.5)' ) jband,ikpt,Ha_eV*eigen(iband, ikpt, isppol)
366 : : end if
367 : : !eigen(iband+band_index)
368 : : ! save eigenvalues
369 : 6656 : eigenvalues_w(jband,ikpt,isppol)=Ha_eV*eigen(iband, ikpt, isppol)
370 : : !eigen(iband+band_index)
371 : : end if
372 : : end do !iband
373 : 20 : band_index=band_index+nband_k
374 : : end do !ikpt
375 : : end do !nsppol
376 [ + - ]: 19 : if(rank==master) then
377 [ + + ]: 39 : do isppol=1,nsppol
378 : 39 : close(iun(isppol))
379 : : end do
380 : 19 : write(msg, '(a,a)' ) ch10,' mlwfovlp : eigenvalues written'
381 : 19 : call wrtout(std_out, msg)
382 : : end if !master
383 : :
384 : 19 : end subroutine write_eigenvalues
385 : :
386 : :
387 : : !-----------------------------------------------------------------------
388 : : ! This subroutine writes the overlap matrix in the w90 format
389 : : ! (see w90io.F90 in w90 package)
390 : : !-----------------------------------------------------------------------
391 : : ! TODO: this routine is doing much more than writing Mmn!! It also calculates etc. The io module should be segregated and cleaned
392 [ - + - + : 19 : subroutine write_Mmn(filew90_mmn, band_in, cm1, ovikp, g1, M_matrix, &
- + ]
393 [ - + - + ]: 19 : & nkpt, nsppol, nntot, mband, num_bands, msg, iam_master, test_matrix_output)
394 : :
395 : : integer, intent(in) :: nsppol, nntot, nkpt
396 : : integer, intent(in) :: mband, num_bands(:)
397 : : character(len=fnlen), intent(in) :: filew90_mmn(nsppol)
398 : : logical, intent(in) :: band_in(mband, nsppol)
399 : : integer,intent(in):: ovikp(:,:)
400 : : integer, intent(in) :: g1(:, :, :)
401 : : real(dp), intent(in) :: cm1(:,:,:,:,:,:)
402 : : logical, intent(in) :: iam_master, test_matrix_output
403 : : complex(dp),intent(inout) :: M_matrix(:,:,:,:,:)
404 : :
405 : : !Local variables-------------------------------
406 : : !scalars
407 : : integer :: isppol, ikpt1, intot,ii, jj,jband1, iband1, jband2, iband2
408 : 38 : integer :: iun(nsppol), units(2)
409 : : character(len=1000) :: msg
410 : :
411 [ + + ]: 57 : units = [std_out, ab_out]
412 : :
413 [ + + ]: 39 : do isppol=1,nsppol !we write separate output files for each isppol
414 : 20 : iun(isppol)=220+isppol
415 [ + - ]: 39 : if (iam_master) then
416 : 20 : open(unit=iun(isppol),file=filew90_mmn(isppol),form='formatted',status='unknown')
417 : 20 : write(iun(isppol),*) "nnkp version 90"
418 : 20 : write(iun(isppol),*) num_bands(isppol),nkpt,nntot
419 : : end if
420 : : end do
421 : :
422 [ + + ]: 39 : do isppol=1,nsppol
423 [ + + ]: 724 : do ikpt1=1,nkpt
424 [ + + ]: 7108 : do intot=1,nntot
425 [ + - ]: 6384 : if( iam_master) then
426 : 6384 : write(iun(isppol),'(2i6,3x,3x,3i5)') ikpt1,ovikp(ikpt1,intot),(g1(jj,ikpt1,intot),jj=1,3)
427 : : end if
428 : : jband2=0
429 [ + + ]: 87536 : do iband2=1,mband ! the first index is faster
430 [ + + ]: 86832 : if(band_in(iband2,isppol)) then
431 : 65472 : jband2=jband2+1
432 : 65472 : jband1=0
433 [ + + ]: 1061056 : do iband1=1,mband
434 [ + + ]: 1061056 : if(band_in(iband1,isppol)) then
435 : 907008 : jband1=jband1+1
436 [ + - ]: 907008 : if(iam_master) write(iun(isppol),*) cm1(1,iband1,iband2,intot,ikpt1,isppol),cm1(2,iband1,iband2,intot,ikpt1,isppol)
437 : : M_matrix(jband1,jband2,intot,ikpt1,isppol)=&
438 : 907008 : & cmplx(cm1(1,iband1,iband2,intot,ikpt1,isppol),cm1(2,iband1,iband2,intot,ikpt1,isppol), kind=dp )
439 : : ! write(2211,*) ikpt1,intot,iband1,iband2
440 : : ! write(2211,*) cm1(1,iband1,iband2,intot,ikpt1,isppol),cm1(2,iband1,iband2,intot,ikpt1,isppol)
441 : : end if ! band_in(iband1)
442 : : end do ! iband1
443 : : end if ! band_in(iband2)
444 : : end do ! iband2
445 : : end do !intot
446 : : end do !ikpt
447 [ + - ]: 39 : if( iam_master ) then
448 : 20 : close(iun(isppol))
449 [ + - ]: 20 : write(msg, '(3a)' ) ' ',trim(filew90_mmn(isppol)),' written'
450 : 20 : call wrtout(std_out, msg)
451 : : end if !rank==master
452 : : end do !isppol
453 : :
454 : :
455 : : ! Individual spinor overlap elements depend on the arbitrary gauge selected
456 : : ! in degenerate subspaces. The complete .mmn file above is always written.
457 [ + - + + ]: 19 : if(iam_master .and. test_matrix_output) then
458 : 17 : write(msg, '(4a)' ) ch10,&
459 : 17 : ' Writing top of the overlap matrix: M_mn(ikb,ik)',ch10,&
460 : 34 : ' m=n=1:3, ikb=1, ik=1'
461 : 17 : call wrtout(units, msg)
462 : :
463 : : ! just write down the first 3 elements
464 [ + + ]: 35 : do isppol=1,nsppol
465 : 18 : write(msg, '( " " )')
466 [ + + ]: 18 : if (nsppol>1 ) then
467 [ + + ]: 2 : if (isppol==1) write(msg,'(2a)')trim(msg),' spin up:'
468 [ + + ]: 2 : if (isppol==2) write(msg,'(2a)')trim(msg),' spin down:'
469 : : end if
470 [ + + ]: 72 : do ii=1,3
471 [ - + ]: 54 : if(ii>num_bands(isppol)) cycle
472 : 54 : write(msg,'(3a)') trim(msg),ch10,'; ( '
473 [ + + ]: 216 : do jj=1,3
474 [ - + ]: 162 : if(jj>num_bands(isppol))cycle
475 : 216 : write(msg, '(a,2f11.6,a)') trim(msg),M_matrix(ii,jj,1,1,isppol),' , '
476 : : end do
477 : 72 : write(msg,'(2a)') trim(msg),' ) '
478 : : end do
479 : 35 : call wrtout(units, msg)
480 : : end do
481 : :
482 : : ! Now write down bottom of the matrix
483 : 17 : write(msg, '(4a)' ) ch10,&
484 : 17 : ' Writing bottom of the overlap matrix: M_mn(ikb,ik)',ch10,&
485 : 34 : ' m=n=num_bands-2:num_bands, ikb=nntot, ik=nkpt'
486 : 17 : call wrtout(units, msg)
487 : :
488 [ + + ]: 35 : do isppol=1,nsppol
489 : 18 : write(msg, '( " " )')
490 [ + + ]: 18 : if (nsppol>1 ) then
491 [ + + ]: 2 : if (isppol==1) write(msg,'(2a)')trim(msg),' spin up:'
492 [ + + ]: 2 : if (isppol==2) write(msg,'(2a)')trim(msg),' spin down:'
493 : : end if
494 [ + + ]: 72 : do ii=num_bands(isppol)-2,num_bands(isppol)
495 [ - + ]: 54 : if(ii<1) cycle
496 : 54 : write(msg,'(3a)') trim(msg),ch10,'; ( '
497 [ + + ]: 216 : do jj=num_bands(isppol)-2,num_bands(isppol)
498 [ - + ]: 162 : if(jj<1)cycle
499 : 216 : write(msg, '(a,2f11.6,a)') trim(msg),M_matrix(ii,jj,nntot,nkpt,isppol),' , '
500 : : end do !j
501 : 72 : write(msg,'(2a)') trim(msg),' ) '
502 : : end do !ii
503 : 35 : call wrtout(units, msg)
504 : : end do !isppol
505 : : end if !rank==master
506 : :
507 : 19 : end subroutine write_Mmn
508 : :
509 : : !!****f* mlwfovlp/read_chkunit
510 : : !! NAME
511 : : !! read_chkunit
512 : : !!
513 : : !! FUNCTION
514 : : !! Function which reads the .chk file produced by Wannier90
515 : : !!
516 : : !! INPUTS
517 : : !!
518 : : !! OUTPUT
519 : : !!
520 : : !! SOURCE
521 : :
522 : 0 : subroutine read_chkunit(seed_name,nkpt,ndimwin,ierr)
523 : :
524 : : !Arguments ------------------------------------
525 : : !scalars
526 : : integer,intent(in) :: nkpt
527 : : character(len=*),intent(in) :: seed_name
528 : : integer,intent(out) :: ierr
529 : : !arrays
530 : : integer,intent(out) :: ndimwin(nkpt)
531 : :
532 : : !Local variables-------------------------------
533 : : !string
534 : : character(len=fnlen) :: fname
535 : : !scalars
536 : : integer :: chk_unit,ios,ikpt
537 : : logical :: have_disentangled
538 : :
539 : : !************************************************************************
540 : :
541 : 0 : chk_unit=get_unit()
542 [ # # ]: 0 : fname=TRIM(seed_name)//'.chk'
543 : 0 : open(unit=chk_unit,file=fname,form='unformatted',status='old',iostat=ios)
544 : :
545 : 0 : ierr=0
546 : 0 : read(chk_unit) ! header ! Date and time
547 : 0 : read(chk_unit) ! ((real_lattice(i,j),i=1,3),j=1,3) ! Real lattice
548 : 0 : read(chk_unit) ! ((recip_lattice(i,j),i=1,3),j=1,3) ! Reciprocal lattice
549 : 0 : read(chk_unit) ! num_kpts
550 : 0 : read(chk_unit) ! ((kpt_latt(i,nkp),i=1,3),nkp=1,num_kpts) ! K-points
551 : 0 : read(chk_unit) ! nntot ! Number of nearest k-point neighbours
552 : 0 : read(chk_unit) ! num_wann ! Number of wannier functions
553 : 0 : read(chk_unit) ! chkpt1 ! Position of checkpoint
554 : 0 : read(chk_unit) have_disentangled ! Whether a disentanglement has been performed
555 [ # # ]: 0 : if (have_disentangled) then
556 : : !read(chk_unit) ! omega_invariant ! Omega invariant
557 : : !read(chk_unit) ((lwindow(i,nkp),i=1,num_bands),nkp=1,num_kpts)
558 : 0 : read(chk_unit) (ndimwin(ikpt),ikpt=1,nkpt)
559 : : !read(chk_unit) (((u_matrix_opt(i,j,nkp),i=1,num_bands),j=1,num_wann),nkp=1,num_kpts)
560 : : else
561 : : ! this is not expected. we should have disentanglement. Report the error.
562 : 0 : ierr=-1
563 : : end if
564 : : !read(chk_unit) (((u_matrix(i,j,k),i=1,num_wann),j=1,num_wann),k=1,num_kpts) ! U_matrix
565 : : !read(chk_unit) ((((m_matrix(i,j,k,l,1),i=1,num_wann),j=1,num_wann),k=1,nntot),l=1,num_kpts) ! M_matrix
566 : : !read(chk_unit) ((wannier_centres(i,j),i=1,3),j=1,num_wann)
567 : 0 : close(chk_unit)
568 : :
569 : 0 : end subroutine read_chkunit
570 : : !!***
571 : :
572 : : end module m_wannier_io
|
