From 24bdff49ca7ff19e46764ae019f82d9d661fbafd Mon Sep 17 00:00:00 2001 From: Jerome Jackson Date: Thu, 27 Jun 2024 09:51:26 +0100 Subject: [PATCH] use matrix names more consistently --- src/io.F90 | 2 +- src/library_extra.F90 | 110 ++++++++++++++------------------------ src/library_interface.F90 | 50 ++++++++--------- src/readwrite.F90 | 3 ++ src/wannier_prog.F90 | 30 ++++++----- 5 files changed, 82 insertions(+), 113 deletions(-) diff --git a/src/io.F90 b/src/io.F90 index 0200e952..713674b5 100644 --- a/src/io.F90 +++ b/src/io.F90 @@ -413,7 +413,7 @@ subroutine prterr(error, ie, istdout, istderr, comm) mesg = 'not set' if (mpirank(comm) == 0) then - ! fixme, report all failing ranks instead of lowest failing rank (current stand) + ! currently this printout will list only the lowest failing rank, not all failing ranks do j = mpisize(comm) - 1, 1, -1 call comms_no_sync_recv(je, 1, j, le, comm) diff --git a/src/library_extra.F90 b/src/library_extra.F90 index f37737c7..16cb5974 100644 --- a/src/library_extra.F90 +++ b/src/library_extra.F90 @@ -118,33 +118,38 @@ subroutine input_reader_special(common_data, seedname, istdout, istderr, ierr) if (common_data%num_bands > common_data%num_wann) then allocate (common_data%dis_manifold%ndimwin(common_data%num_kpts), stat=ierr) if (ierr /= 0) then - call set_error_alloc(error, 'Error allocating ndimwin in input_reader_special() call', common_data%comm) + call set_error_alloc(error, & + 'Error allocating ndimwin in input_reader_special() call', common_data%comm) call prterr(error, ierr, istdout, istderr, common_data%comm) return endif allocate (common_data%dis_manifold%nfirstwin(common_data%num_kpts), stat=ierr) if (ierr /= 0) then - call set_error_alloc(error, 'Error allocating nfirstwin in input_reader_special() call', common_data%comm) + call set_error_alloc(error, & + 'Error allocating nfirstwin in input_reader_special() call', common_data%comm) call prterr(error, ierr, istdout, istderr, common_data%comm) return endif allocate (common_data%dis_manifold%lwindow(common_data%num_bands, common_data%num_kpts), stat=ierr) if (ierr /= 0) then - call set_error_alloc(error, 'Error allocating lwindow in input_reader_special() call', common_data%comm) + call set_error_alloc(error, & + 'Error allocating lwindow in input_reader_special() call', common_data%comm) call prterr(error, ierr, istdout, istderr, common_data%comm) return endif endif allocate (common_data%wannier_data%centres(3, common_data%num_wann), stat=ierr) if (ierr /= 0) then - call set_error_alloc(error, 'Error allocating wannier_centres in input_reader_special() call', common_data%comm) + call set_error_alloc(error, & + 'Error allocating wannier_centres in input_reader_special() call', common_data%comm) call prterr(error, ierr, istdout, istderr, common_data%comm) return endif common_data%wannier_data%centres = 0.0_dp allocate (common_data%wannier_data%spreads(common_data%num_wann), stat=ierr) if (ierr /= 0) then - call set_error_alloc(error, 'Error in allocating wannier_spreads in input_reader_special() call', common_data%comm) + call set_error_alloc(error, & + 'Error in allocating wannier_spreads in input_reader_special() call', common_data%comm) call prterr(error, ierr, istdout, istderr, common_data%comm) return endif @@ -179,8 +184,10 @@ subroutine write_kmesh(common_data, istdout, istderr, ierr) ierr = 0 - if (mpirank(common_data%comm) == 0) then - if (.not. allocated(common_data%kmesh_info%nnlist)) call w90_create_kmesh(common_data, istdout, istderr, ierr) + if (mpirank(common_data%comm) == 0) then ! root only + if (.not. allocated(common_data%kmesh_info%nnlist)) then + call w90_create_kmesh(common_data, istdout, istderr, ierr) + endif call kmesh_write(common_data%exclude_bands, common_data%kmesh_info, & common_data%select_proj%auto_projections, common_data%proj_input, & @@ -214,14 +221,13 @@ subroutine overlaps(common_data, istdout, istderr, ierr) ierr = 0 if (.not. common_data%setup_complete) then - call set_error_fatal(error, 'kmesh is not setup before calling overlaps read routine', common_data%comm) + call set_error_fatal(error, 'Error: kmesh is not setup before reading overlap matrix', common_data%comm) call prterr(error, ierr, istdout, istderr, common_data%comm) return - if (ierr > 0) return endif ! projections are stored in u_opt - call overlap_read(common_data%kmesh_info, common_data%select_proj, common_data%u_opt, & + call overlap_read(common_data%kmesh_info, common_data%select_proj, common_data%u_matrix_opt, & common_data%m_matrix_local, common_data%num_bands, common_data%num_kpts, & common_data%num_proj, common_data%num_wann, common_data%print_output, & common_data%print_output%timing_level, cp_pp, common_data%use_bloch_phases, & @@ -231,33 +237,6 @@ subroutine overlaps(common_data, istdout, istderr, ierr) call prterr(error, ierr, istdout, istderr, common_data%comm) return endif - - ! Check Mmn(k,b) is symmetric in m and n for gamma_only case - ! if (gamma_only) call overlap_check_m_symmetry() - ! - ! If we don't need to disentangle we can now convert from A to U - ! And rotate M accordingly - ! Jan 2023, Jerome Jackson, moved overlap_project outside of overlap_read - ! if ((.not. disentanglement) .and. (.not. cp_pp) .and. (.not. use_bloch_phases)) then - ! if (.not. gamma_only) then - ! call overlap_project(sitesym, m_matrix_local, au_matrix, kmesh_info%nnlist, & - ! kmesh_info%nntot, num_bands, num_kpts, num_wann, timing_level, & - ! lsitesymmetry, stdout, timer, dist_k, error, comm) - ! else - ! call overlap_project_gamma(m_matrix_local, au_matrix, kmesh_info%nntot, num_wann, & - ! timing_level, stdout, timer, error, comm) - ! endif - ! if (allocated(error)) return - ! endif - ! - !~[aam] - !~ if( gamma_only .and. use_bloch_phases ) then - !~ write(stdout,'(1x,"+",76("-"),"+")') - !~ write(stdout,'(3x,a)') 'WARNING: gamma_only and use_bloch_phases ' - !~ write(stdout,'(3x,a)') ' M must be calculated from *real* Bloch functions' - !~ write(stdout,'(1x,"+",76("-"),"+")') - !~ end if - ![ysl-e] end subroutine overlaps subroutine read_eigvals(common_data, eigval, istdout, istderr, ierr) @@ -279,11 +258,13 @@ subroutine read_eigvals(common_data, eigval, istdout, istderr, ierr) ierr = 0 if (size(eigval, 1) /= common_data%num_bands) then - call set_error_fatal(error, 'eigval not dimensioned correctly (num_bands,num_kpts) in read_eigvals', common_data%comm) + call set_error_fatal(error, & + 'Error: eigval not dimensioned correctly (num_bands,num_kpts) in read_eigvals', common_data%comm) call prterr(error, ierr, istdout, istderr, common_data%comm) return elseif (size(eigval, 2) /= common_data%num_kpts) then - call set_error_fatal(error, 'eigval not dimensioned correctly (num_bands,num_kpts) in read_eigvals', common_data%comm) + call set_error_fatal(error, & + 'Error: eigval not dimensioned correctly (num_bands,num_kpts) in read_eigvals', common_data%comm) call prterr(error, ierr, istdout, istderr, common_data%comm) return endif @@ -294,7 +275,8 @@ subroutine read_eigvals(common_data, eigval, istdout, istderr, ierr) call prterr(error, ierr, istdout, istderr, common_data%comm) return else if (.not. eig_found) then - call set_error_fatal(error, 'failed to read eigenvalues file in read_eigvals call', common_data%comm) + call set_error_fatal(error, & + 'Error: failed to read eigenvalues file in read_eigvals', common_data%comm) call prterr(error, ierr, istdout, istderr, common_data%comm) return endif @@ -315,7 +297,7 @@ subroutine write_chkpt(common_data, label, istdout, istderr, ierr) type(lib_common_type), target, intent(in) :: common_data ! local variables - complex(kind=dp), allocatable :: u(:, :, :), uopt(:, :, :), m(:, :, :, :) + complex(kind=dp), allocatable :: m(:, :, :, :) integer, allocatable :: global_k(:) integer, pointer :: nw, nb, nk, nn integer :: rank, nkrank, ikg, ikl, istat @@ -323,32 +305,37 @@ subroutine write_chkpt(common_data, label, istdout, istderr, ierr) ierr = 0 rank = mpirank(common_data%comm) + nkrank = count(common_data%dist_kpoints == rank) nb => common_data%num_bands nk => common_data%num_kpts nn => common_data%kmesh_info%nntot nw => common_data%num_wann - if (.not. associated(common_data%u_opt)) then - call set_error_fatal(error, 'u_opt not set for write_chkpt call', common_data%comm) + if (.not. associated(common_data%u_matrix_opt)) then + call set_error_fatal(error, & + 'Error: u_matrix_opt not associated for write_chkpt call', common_data%comm) else if (.not. associated(common_data%u_matrix)) then - call set_error_fatal(error, 'u_matrix not set for write_chkpt call', common_data%comm) + call set_error_fatal(error, & + 'Error: u_matrix not associated for write_chkpt call', common_data%comm) else if (.not. associated(common_data%m_matrix_local)) then - call set_error_fatal(error, 'm_matrix_local not set for write_chkpt call', common_data%comm) + call set_error_fatal(error, & + 'Error: m_matrix_local not set for write_chkpt call', common_data%comm) endif if (allocated(error)) then call prterr(error, ierr, istdout, istderr, common_data%comm) return endif - nkrank = count(common_data%dist_kpoints == rank) allocate (global_k(nkrank), stat=istat) if (istat /= 0) then call set_error_alloc(error, 'Error allocating global_k in write_chkpt', common_data%comm) call prterr(error, ierr, istdout, istderr, common_data%comm) return endif - global_k = huge(1); ikl = 1 + + global_k = huge(1) + ikl = 1 do ikg = 1, nk if (rank == common_data%dist_kpoints(ikg)) then global_k(ikl) = ikg @@ -356,30 +343,23 @@ subroutine write_chkpt(common_data, label, istdout, istderr, ierr) endif enddo + ! reassemble full m matrix by MPI reduction + ! ! allocating and partially assigning the full matrix on all ranks and reducing is a terrible idea ! alternatively, allocate on root and use point-to-point ! or, if required only for checkpoint file writing, then use mpi-io (but needs to be ordered io, alas) ! or, even better, use parallel hdf5. JJ Nov 22 - allocate (u(nw, nw, nk), stat=istat) ! all kpts - if (istat /= 0) call set_error_alloc(error, 'Error allocating u in write_chkpt', common_data%comm) - allocate (uopt(nb, nw, nk), stat=istat) ! all kpts - if (istat /= 0) call set_error_alloc(error, 'Error allocating uopt in write_chkpt', common_data%comm) allocate (m(nw, nw, nn, nk), stat=istat) ! all kpts if (istat /= 0) call set_error_alloc(error, 'Error allocating m in write_chkpt', common_data%comm) if (allocated(error)) then call prterr(error, ierr, istdout, istderr, common_data%comm) return endif - - u(:, :, :) = common_data%u_matrix - uopt(:, :, :) = common_data%u_opt m(:, :, :, :) = 0.d0 - do ikl = 1, nkrank ikg = global_k(ikl) m(:, :, :, ikg) = common_data%m_matrix_local(1:nw, 1:nw, :, ikl) enddo - call comms_reduce(m(1, 1, 1, 1), nw*nw*nn*nk, 'SUM', error, common_data%comm) if (allocated(error)) then call prterr(error, ierr, istdout, istderr, common_data%comm) @@ -390,7 +370,8 @@ subroutine write_chkpt(common_data, label, istdout, istderr, ierr) call w90_wannier90_readwrite_write_chkpt(label, common_data%exclude_bands, & common_data%wannier_data, common_data%kmesh_info, & common_data%kpt_latt, nk, common_data%dis_manifold, & - nb, nw, u, uopt, m, common_data%mp_grid, & + nb, nw, common_data%u_matrix, & + common_data%u_matrix_opt, m, common_data%mp_grid, & common_data%real_lattice, & common_data%omega%invariant, & common_data%have_disentangled, & @@ -398,19 +379,9 @@ subroutine write_chkpt(common_data, label, istdout, istderr, ierr) common_data%seedname) endif - deallocate (u, stat=istat) - if (istat /= 0) then - call set_error_dealloc(error, 'Error deallocating u in write_chkpt', common_data%comm) - endif - deallocate (uopt, stat=istat) - if (istat /= 0) then - call set_error_dealloc(error, 'Error deallocating uopt in write_chkpt', common_data%comm) - endif deallocate (m, stat=istat) if (istat /= 0) then call set_error_dealloc(error, 'Error deallocating m in write_chkpt', common_data%comm) - endif - if (allocated(error)) then call prterr(error, ierr, istdout, istderr, common_data%comm) return endif @@ -449,7 +420,6 @@ subroutine read_chkpt(common_data, checkpoint, istdout, istderr, ierr) ! alternatively, allocate on root and use point-to-point ! or, if required only for checkpoint file writing, then use mpi-io (but needs to be ordered io, alas) ! or, even better, use parallel hdf5 - ! fixme, check allocation status of u, uopt? allocate (m(nw, nw, nn, nk), stat=istat) ! all kpts if (istat /= 0) then call set_error_alloc(error, 'Error allocating m in read_chkpt', common_data%comm) @@ -463,7 +433,7 @@ subroutine read_chkpt(common_data, checkpoint, istdout, istderr, ierr) call w90_readwrite_read_chkpt(common_data%dis_manifold, common_data%exclude_bands, & common_data%kmesh_info, common_data%kpt_latt, & common_data%wannier_data, m, common_data%u_matrix, & - common_data%u_opt, common_data%real_lattice, & + common_data%u_matrix_opt, common_data%real_lattice, & common_data%omega%invariant, common_data%mp_grid, nb, & nexclude, nk, nw, checkpoint, common_data%have_disentangled, & ispostw90, common_data%seedname, istdout, error, & @@ -476,7 +446,7 @@ subroutine read_chkpt(common_data, checkpoint, istdout, istderr, ierr) ! scatter from m_matrix to m_matrix_local (cf overlap_read) call w90_readwrite_chkpt_dist(common_data%dis_manifold, common_data%wannier_data, & - common_data%u_matrix, common_data%u_opt, m, & + common_data%u_matrix, common_data%u_matrix_opt, m, & common_data%m_matrix_local, common_data%omega%invariant, & nb, nk, nw, nn, checkpoint, common_data%have_disentangled, & common_data%dist_kpoints, error, common_data%comm) diff --git a/src/library_interface.F90 b/src/library_interface.F90 index 948781cb..030fd227 100644 --- a/src/library_interface.F90 +++ b/src/library_interface.F90 @@ -87,12 +87,10 @@ module w90_library !! pointer to matrix "a" (projections), dimensioned (number of bands, number WF, number of k-points) complex(kind=dp), pointer :: m_matrix_local(:, :, :, :) => null() !! pointer to rank local part of "m", dimensioned (number of bands, number of bands, FD neighbours, number of k-points on this rank) - complex(kind=dp), pointer :: m_matrix(:, :, :, :) => null() - !! pointer to matrix "m" (overlaps), dimensioned (number of bands, number of bands, FD neighbours, number of k-points) complex(kind=dp), pointer :: u_matrix(:, :, :) => null() !! pointer to matrix "u", dimensionsed (number of WF, number of WF, number of k-points) - complex(kind=dp), pointer :: u_opt(:, :, :) => null() - !! pointer to matrix "u_opt", dimensioned (number of bands, number of WF, number of k-points) + complex(kind=dp), pointer :: u_matrix_opt(:, :, :) => null() + !! pointer to matrix "u_matrix_opt", dimensioned (number of bands, number of WF, number of k-points) real(kind=dp), pointer :: eigval(:, :) => null() !! pointer to eigenvalue array dimensioned (number of bands, number of k-points) @@ -531,8 +529,8 @@ subroutine w90_disentangle(common_data, istdout, istderr, ierr) call set_error_fatal(error, 'Error: u_matrix not associated for w90_disentangle() call', common_data%comm) call prterr(error, ierr, istdout, istderr, common_data%comm) return - else if (.not. associated(common_data%u_opt)) then - call set_error_fatal(error, 'Error: u_opt not associated for w90_disentangle() call', common_data%comm) + else if (.not. associated(common_data%u_matrix_opt)) then + call set_error_fatal(error, 'Error: u_matrix_opt not associated for w90_disentangle() call', common_data%comm) call prterr(error, ierr, istdout, istderr, common_data%comm) return else if (.not. associated(common_data%eigval)) then @@ -545,16 +543,9 @@ subroutine w90_disentangle(common_data, istdout, istderr, ierr) if (common_data%dis_manifold%win_min == -huge(0.0_dp)) common_data%dis_manifold%win_min = minval(common_data%eigval) if (common_data%dis_manifold%win_max == huge(0.0_dp)) common_data%dis_manifold%win_max = maxval(common_data%eigval) if (common_data%dis_manifold%frozen_states) then - !write(*,*) shape(common_data%eigval) - !write(*,*) shape(common_data%dis_manifold%lwindow) - if (common_data%dis_manifold%froz_min == -huge(0.0_dp) .and. allocated(common_data%exclude_bands)) then - ioff = maxval(common_data%exclude_bands) - ioff = 0 !fixme - common_data%dis_manifold%froz_min = minval(common_data%eigval(ioff + 1, :)) + if (common_data%dis_manifold%froz_min == -huge(0.0_dp)) then + common_data%dis_manifold%froz_min = minval(common_data%eigval(:, :)) endif - !common_data%dis_manifold%froz_min = minval(common_data%eigval, mask=common_data%dis_manifold%lwindow) - !write(*,*)common_data%dis_manifold%froz_min,allocated(common_data%dis_manifold%lwindow),& - ! count(common_data%dis_manifold%lwindow) endif ! condition for disentanglement is number of bands > number of WF @@ -562,7 +553,7 @@ subroutine w90_disentangle(common_data, istdout, istderr, ierr) call dis_main(common_data%dis_control, common_data%dis_spheres, common_data%dis_manifold, & common_data%kmesh_info, common_data%kpt_latt, common_data%sitesym, & common_data%print_output, common_data%m_matrix_local, common_data%u_matrix, & - common_data%u_opt, common_data%eigval, common_data%real_lattice, & + common_data%u_matrix_opt, common_data%eigval, common_data%real_lattice, & common_data%omega%invariant, common_data%num_bands, common_data%num_kpts, & common_data%num_wann, common_data%gamma_only, common_data%lsitesymmetry, & istdout, common_data%timer, common_data%dist_kpoints, error, common_data%comm) @@ -571,6 +562,7 @@ subroutine w90_disentangle(common_data, istdout, istderr, ierr) return endif + !fixme, aliasing of input and output m_matrix_local here (resizing nb,nb -> nw,nw) call setup_m_loc(common_data%kmesh_info, common_data%print_output, common_data%m_matrix_local, & common_data%m_matrix_local, common_data%u_matrix, common_data%num_bands, & common_data%num_kpts, common_data%num_wann, common_data%timer, & @@ -606,8 +598,8 @@ subroutine w90_project_overlap(common_data, istdout, istderr, ierr) call set_error_fatal(error, 'm_matrix_local not set for w90_project_overlap call', common_data%comm) call prterr(error, ierr, istdout, istderr, common_data%comm) return - else if (.not. associated(common_data%u_opt)) then - call set_error_fatal(error, 'u_opt not set for w90_project_overlap call', common_data%comm) + else if (.not. associated(common_data%u_matrix_opt)) then + call set_error_fatal(error, 'u_matrix_opt not set for w90_project_overlap call', common_data%comm) call prterr(error, ierr, istdout, istderr, common_data%comm) return else if (.not. associated(common_data%u_matrix)) then @@ -625,11 +617,11 @@ subroutine w90_project_overlap(common_data, istdout, istderr, ierr) endif ! fixme, document! - common_data%u_matrix(:, :, :) = common_data%u_opt(:, :, :) ! u_opt contains initial projections - common_data%u_opt(:, :, :) = 0.d0 + common_data%u_matrix(:, :, :) = common_data%u_matrix_opt(:, :, :) ! u_matrix_opt contains initial projections + common_data%u_matrix_opt(:, :, :) = 0.d0 do ik = 1, common_data%num_kpts do iw = 1, common_data%num_wann - common_data%u_opt(iw, iw, ik) = 1.d0 + common_data%u_matrix_opt(iw, iw, ik) = 1.d0 enddo enddo @@ -746,7 +738,7 @@ subroutine w90_plot(common_data, istdout, istderr, ierr) common_data%kpoint_path, common_data%print_output, common_data%wannier_data, & common_data%wann_plot, common_data%ws_region, common_data%w90_calculation, & common_data%ham_k, common_data%ham_r, common_data%m_matrix_local, & - common_data%u_matrix, common_data%u_opt, common_data%eigval, & + common_data%u_matrix, common_data%u_matrix_opt, common_data%eigval, & common_data%real_lattice, common_data%wannier_centres_translated, & common_data%physics%bohr, common_data%irvec, common_data%mp_grid, & common_data%ndegen, common_data%shift_vec, common_data%nrpts, & @@ -788,7 +780,7 @@ subroutine w90_transport(common_data, istdout, istderr, ierr) common_data%output_file, common_data%real_space_ham, common_data%tran, & common_data%print_output, common_data%wannier_data, common_data%ws_region, & common_data%w90_calculation, common_data%ham_k, common_data%ham_r, & - common_data%u_matrix, common_data%u_opt, common_data%eigval, & + common_data%u_matrix, common_data%u_matrix_opt, common_data%eigval, & common_data%real_lattice, common_data%wannier_centres_translated, & common_data%irvec, common_data%mp_grid, common_data%ndegen, & common_data%shift_vec, common_data%nrpts, common_data%num_bands, & @@ -804,13 +796,13 @@ subroutine w90_transport(common_data, istdout, istderr, ierr) endif end subroutine w90_transport - subroutine w90_set_m_local(common_data, m_orig) ! m_matrix_local_orig + subroutine w90_set_m_local(common_data, m_matrix_local) ! m_matrix_local_orig implicit none type(lib_common_type), intent(inout) :: common_data - complex(kind=dp), intent(inout), target :: m_orig(:, :, :, :) + complex(kind=dp), intent(inout), target :: m_matrix_local(:, :, :, :) - common_data%m_matrix_local => m_orig + common_data%m_matrix_local => m_matrix_local end subroutine w90_set_m_local subroutine w90_set_u_matrix(common_data, u_matrix) @@ -822,13 +814,13 @@ subroutine w90_set_u_matrix(common_data, u_matrix) common_data%u_matrix => u_matrix end subroutine w90_set_u_matrix - subroutine w90_set_u_opt(common_data, u_opt) + subroutine w90_set_u_opt(common_data, u_matrix_opt) implicit none type(lib_common_type), intent(inout) :: common_data - complex(kind=dp), intent(inout), target :: u_opt(:, :, :) + complex(kind=dp), intent(inout), target :: u_matrix_opt(:, :, :) - common_data%u_opt => u_opt + common_data%u_matrix_opt => u_matrix_opt end subroutine w90_set_u_opt subroutine w90_set_eigval(common_data, eigval) diff --git a/src/readwrite.F90 b/src/readwrite.F90 index 7bbf1f9d..faefa2c6 100644 --- a/src/readwrite.F90 +++ b/src/readwrite.F90 @@ -2308,6 +2308,9 @@ subroutine w90_readwrite_read_chkpt_matrices(dis_manifold, kmesh_info, wannier_d ! U_matrix_opt read (chk_unit, err=124) (((u_matrix_opt(i, j, nkp), i=1, num_bands), j=1, num_wann), nkp=1, num_kpts) + else + ! if not read, u_matrix_opt must be explicitly zeroed + u_matrix_opt(:, :, :) = 0 endif ! U_matrix diff --git a/src/wannier_prog.F90 b/src/wannier_prog.F90 index deb5ab20..26081c04 100644 --- a/src/wannier_prog.F90 +++ b/src/wannier_prog.F90 @@ -73,9 +73,9 @@ program wannier character(len=:), allocatable :: seedname, progname, cpstatus character(len=:), pointer :: restart - complex(kind=dp), allocatable :: mloc(:, :, :, :) - complex(kind=dp), allocatable :: u(:, :, :) - complex(kind=dp), allocatable :: uopt(:, :, :) + complex(kind=dp), allocatable :: m_matrix_loc(:, :, :, :) + complex(kind=dp), allocatable :: u_matrix(:, :, :) + complex(kind=dp), allocatable :: u_matrix_opt(:, :, :) real(kind=dp), allocatable :: eigval(:, :) integer, allocatable :: distk(:) integer :: i, ctr @@ -160,9 +160,9 @@ program wannier if (ierr /= 0) stop ! special branch for writing nnkp file + ! exit immediately after writing the nnkp file if (pp) then - ! please only invoke on rank 0 - call write_kmesh(common_data, stdout, stderr, ierr) + call write_kmesh(common_data, stdout, stderr, ierr) ! only active on rank 0 if (ierr /= 0) stop if (rank == 0) close (unit=stderr, status='delete') #ifdef MPI @@ -184,15 +184,13 @@ program wannier nkl = count(distk == rank) ! number of kpoints this rank !write (*, *) 'rank, nw, nb, nk, nn, nk(rank): ', rank, nw, nb, nk, nn, nkl - allocate (mloc(nb, nb, nn, nkl)) - allocate (u(nw, nw, nk)) - allocate (uopt(nb, nw, nk)) + allocate (m_matrix_loc(nb, nb, nn, nkl)) + allocate (u_matrix(nw, nw, nk)) + allocate (u_matrix_opt(nb, nw, nk)) - call w90_set_m_local(common_data, mloc) ! we don't need global m - call w90_set_u_matrix(common_data, u) - call w90_set_u_opt(common_data, uopt) - - uopt = 0.d0 ! required when not disentangling + call w90_set_m_local(common_data, m_matrix_loc) ! we don't need global m + call w90_set_u_matrix(common_data, u_matrix) + call w90_set_u_opt(common_data, u_matrix_opt) ! restart system lovlp = .true. @@ -286,6 +284,12 @@ program wannier if (ierr /= 0) stop endif + if (need_eigvals) deallocate (eigval) + deallocate (distk) + deallocate (m_matrix_loc) + deallocate (u_matrix) + deallocate (u_matrix_opt) + call print_times(common_data, stdout) if (rank == 0) close (unit=stderr, status='delete') #ifdef MPI