Read checkpoint file IMPORTANT! If you change the chkpt format, adapt accordingly also the w90chk2chk.x utility!
Note on parallelization: this function should be called from the root node only!
This function should be called
subroutine param_read_chkpt()
!=================================================!
!! Read checkpoint file
!! IMPORTANT! If you change the chkpt format, adapt
!! accordingly also the w90chk2chk.x utility!
!!
!! Note on parallelization: this function should be called
!! from the root node only!
!!
!! This function should be called
!=================================================!
use w90_constants, only: eps6
use w90_io, only: io_error, io_file_unit, stdout, seedname
implicit none
integer :: chk_unit, nkp, i, j, k, l, ntmp, ierr
character(len=33) :: header
real(kind=dp) :: tmp_latt(3, 3), tmp_kpt_latt(3, num_kpts)
integer :: tmp_excl_bands(1:num_exclude_bands), tmp_mp_grid(1:3)
write (stdout, '(1x,3a)') 'Reading restart information from file ', trim(seedname), '.chk :'
chk_unit = io_file_unit()
open (unit=chk_unit, file=trim(seedname)//'.chk', status='old', form='unformatted', err=121)
! Read comment line
read (chk_unit) header
write (stdout, '(1x,a)', advance='no') trim(header)
! Consistency checks
read (chk_unit) ntmp ! Number of bands
if (ntmp .ne. num_bands) call io_error('param_read_chk: Mismatch in num_bands')
read (chk_unit) ntmp ! Number of excluded bands
if (ntmp .ne. num_exclude_bands) &
call io_error('param_read_chk: Mismatch in num_exclude_bands')
read (chk_unit) (tmp_excl_bands(i), i=1, num_exclude_bands) ! Excluded bands
do i = 1, num_exclude_bands
if (tmp_excl_bands(i) .ne. exclude_bands(i)) &
call io_error('param_read_chk: Mismatch in exclude_bands')
enddo
read (chk_unit) ((tmp_latt(i, j), i=1, 3), j=1, 3) ! Real lattice
do j = 1, 3
do i = 1, 3
if (abs(tmp_latt(i, j) - real_lattice(i, j)) .gt. eps6) &
call io_error('param_read_chk: Mismatch in real_lattice')
enddo
enddo
read (chk_unit) ((tmp_latt(i, j), i=1, 3), j=1, 3) ! Reciprocal lattice
do j = 1, 3
do i = 1, 3
if (abs(tmp_latt(i, j) - recip_lattice(i, j)) .gt. eps6) &
call io_error('param_read_chk: Mismatch in recip_lattice')
enddo
enddo
read (chk_unit) ntmp ! K-points
if (ntmp .ne. num_kpts) &
call io_error('param_read_chk: Mismatch in num_kpts')
read (chk_unit) (tmp_mp_grid(i), i=1, 3) ! M-P grid
do i = 1, 3
if (tmp_mp_grid(i) .ne. mp_grid(i)) &
call io_error('param_read_chk: Mismatch in mp_grid')
enddo
read (chk_unit) ((tmp_kpt_latt(i, nkp), i=1, 3), nkp=1, num_kpts)
do nkp = 1, num_kpts
do i = 1, 3
if (abs(tmp_kpt_latt(i, nkp) - kpt_latt(i, nkp)) .gt. eps6) &
call io_error('param_read_chk: Mismatch in kpt_latt')
enddo
enddo
read (chk_unit) ntmp ! nntot
if (ntmp .ne. nntot) &
call io_error('param_read_chk: Mismatch in nntot')
read (chk_unit) ntmp ! num_wann
if (ntmp .ne. num_wann) &
call io_error('param_read_chk: Mismatch in num_wann')
! End of consistency checks
read (chk_unit) checkpoint ! checkpoint
checkpoint = adjustl(trim(checkpoint))
read (chk_unit) have_disentangled ! whether a disentanglement has been performed
if (have_disentangled) then
read (chk_unit) omega_invariant ! omega invariant
! lwindow
if (.not. allocated(lwindow)) then
allocate (lwindow(num_bands, num_kpts), stat=ierr)
if (ierr /= 0) call io_error('Error allocating lwindow in param_read_chkpt')
endif
read (chk_unit, err=122) ((lwindow(i, nkp), i=1, num_bands), nkp=1, num_kpts)
! ndimwin
if (.not. allocated(ndimwin)) then
allocate (ndimwin(num_kpts), stat=ierr)
if (ierr /= 0) call io_error('Error allocating ndimwin in param_read_chkpt')
endif
read (chk_unit, err=123) (ndimwin(nkp), nkp=1, num_kpts)
! U_matrix_opt
if (.not. allocated(u_matrix_opt)) then
allocate (u_matrix_opt(num_bands, num_wann, num_kpts), stat=ierr)
if (ierr /= 0) call io_error('Error allocating u_matrix_opt in param_read_chkpt')
endif
read (chk_unit, err=124) (((u_matrix_opt(i, j, nkp), i=1, num_bands), j=1, num_wann), nkp=1, num_kpts)
endif
! U_matrix
if (.not. allocated(u_matrix)) then
allocate (u_matrix(num_wann, num_wann, num_kpts), stat=ierr)
if (ierr /= 0) call io_error('Error allocating u_matrix in param_read_chkpt')
endif
read (chk_unit, err=125) (((u_matrix(i, j, k), i=1, num_wann), j=1, num_wann), k=1, num_kpts)
! M_matrix
if (.not. allocated(m_matrix)) then
allocate (m_matrix(num_wann, num_wann, nntot, num_kpts), stat=ierr)
if (ierr /= 0) call io_error('Error allocating m_matrix in param_read_chkpt')
endif
read (chk_unit, err=126) ((((m_matrix(i, j, k, l), i=1, num_wann), j=1, num_wann), k=1, nntot), l=1, num_kpts)
! wannier_centres
read (chk_unit, err=127) ((wannier_centres(i, j), i=1, 3), j=1, num_wann)
! wannier spreads
read (chk_unit, err=128) (wannier_spreads(i), i=1, num_wann)
close (chk_unit)
write (stdout, '(a/)') ' ... done'
return
121 if (ispostw90) then
call io_error('Error opening '//trim(seedname)//'.chk in param_read_chkpt: did you run wannier90.x first?')
else
call io_error('Error opening '//trim(seedname)//'.chk in param_read_chkpt')
end if
122 call io_error('Error reading lwindow from '//trim(seedname)//'.chk in param_read_chkpt')
123 call io_error('Error reading ndimwin from '//trim(seedname)//'.chk in param_read_chkpt')
124 call io_error('Error reading u_matrix_opt from '//trim(seedname)//'.chk in param_read_chkpt')
125 call io_error('Error reading u_matrix from '//trim(seedname)//'.chk in param_read_chkpt')
126 call io_error('Error reading m_matrix from '//trim(seedname)//'.chk in param_read_chkpt')
127 call io_error('Error reading wannier_centres from '//trim(seedname)//'.chk in param_read_chkpt')
128 call io_error('Error reading wannier_spreads from '//trim(seedname)//'.chk in param_read_chkpt')
end subroutine param_read_chkpt