diff --git a/docs/docs/user_guide/wannier90/parameters.md b/docs/docs/user_guide/wannier90/parameters.md index 04a70459..a389cae1 100644 --- a/docs/docs/user_guide/wannier90/parameters.md +++ b/docs/docs/user_guide/wannier90/parameters.md @@ -151,7 +151,7 @@ physical value, L for a logical value and S for a text string. | num_print_cycles | I | Control frequency of printing | | write_r2mn | L | Write matrix elements of $r^2$ between WF to file | | guiding_centres | L | Use guiding centres | -| mv_functional | L | If `true`, use Marzari-Vanderbilt functional and if `false`, use Stengel-Spaldin functional | +| use_ss_functional | L | If `true`, use Stengel-Spaldin functional and if `false`, use Marzari-Vanderbilt functional | | num_guide_cycles | I | Frequency of guiding centres | | num_no_guide_iter | I | The number of iterations after which guiding centres are used | | trial_step \* | R | The trial step length for the parabolic line search during the minimisation of $\Omega$ | @@ -947,15 +947,18 @@ not possible if an explicit projection block is not defined). The default value is `false`. -### `logical :: mv_functional` +### `logical :: use_ss_functional` -If `true`, use Marzari-Vanderbilt functional and if `false`, use Stengel-Spaldin functional. +If `true`, use Stengel-Spaldin spread functional and if `false`, +use Marzari-Vanderbilt spread functional. -Both functionals converge to the same behavior if infinitely fine grid is used. Both of them are translationally invariant, but only Stengel-Spaldin functional is size consistent. +Both functionals converge to the same behavior if infinitely fine grid is used. +Both of them are translationally invariant, +but only Stengel-Spaldin functional is size consistent. For more information, refer to Phys. Rev. B 73, 075121. -The default value is `true`. +The default value is `false`. ### `integer :: num_guide_cycles` diff --git a/src/readwrite.F90 b/src/readwrite.F90 index 29e39be8..4a4b7bdc 100644 --- a/src/readwrite.F90 +++ b/src/readwrite.F90 @@ -965,7 +965,7 @@ subroutine w90_readwrite_clear_keywords(settings, comm) call w90_readwrite_get_keyword(settings, 'num_shells', found, error, comm) call w90_readwrite_get_keyword(settings, 'num_valence_bands', found, error, comm) call w90_readwrite_get_keyword(settings, 'num_wann', found, error, comm) - call w90_readwrite_get_keyword(settings, 'mv_functional', found, error, comm) + call w90_readwrite_get_keyword(settings, 'use_ss_functional', found, error, comm) call w90_readwrite_get_keyword(settings, 'one_dim_axis', found, error, comm) call w90_readwrite_get_keyword(settings, 'optimisation', found, error, comm) call w90_readwrite_get_keyword(settings, 'postproc_setup', found, error, comm) diff --git a/src/wannier90_readwrite.F90 b/src/wannier90_readwrite.F90 index c0f513d9..c6d27c06 100644 --- a/src/wannier90_readwrite.F90 +++ b/src/wannier90_readwrite.F90 @@ -713,8 +713,8 @@ subroutine w90_wannier90_readwrite_read_wannierise(settings, wann_control, num_w l_value=wann_control%guiding_centres%enable) if (allocated(error)) return - call w90_readwrite_get_keyword(settings, 'mv_functional', found, error, comm, & - l_value=wann_control%mv_functional) + call w90_readwrite_get_keyword(settings, 'use_ss_functional', found, error, comm, & + l_value=wann_control%use_ss_functional) if (allocated(error)) return call w90_readwrite_get_keyword(settings, 'num_guide_cycles', found, error, comm, & diff --git a/src/wannier90_types.F90 b/src/wannier90_types.F90 index aa09dfb7..3e367a51 100644 --- a/src/wannier90_types.F90 +++ b/src/wannier90_types.F90 @@ -182,7 +182,7 @@ module w90_wannier90_types real(kind=dp) :: conv_tol = 1.0e-10_dp integer :: conv_window type(guiding_centres_type) :: guiding_centres - logical :: mv_functional = .true. + logical :: use_ss_functional = .false. real(kind=dp) :: fixed_step = -999.0_dp real(kind=dp) :: trial_step = 2.0_dp logical :: precond = .false. diff --git a/src/wannierise.F90 b/src/wannierise.F90 index 4e0c791e..e53f51f3 100644 --- a/src/wannierise.F90 +++ b/src/wannierise.F90 @@ -435,7 +435,7 @@ subroutine wann_main(ham_logical, kmesh_info, kpt_latt, wann_control, omega, sit irguide = 0 if (wann_control%guiding_centres%enable .and. (wann_control%guiding_centres%num_no_guide_iter .le. 0)) then call wann_phases(csheet, sheet, rguide, irguide, num_wann, kmesh_info, num_kpts, & - wann_control%mv_functional, m_matrix_loc, rnkb, print_output%timing_level, & + wann_control%use_ss_functional, m_matrix_loc, rnkb, print_output%timing_level, & print_output%iprint, timer, nkrank, global_k, error, comm) if (allocated(error)) return @@ -450,7 +450,7 @@ subroutine wann_main(ham_logical, kmesh_info, kpt_latt, wann_control, omega, sit ! calculate initial centers and spread call wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, kmesh_info, & - num_kpts, print_output, wann_control%mv_functional, wann_control%constrain, & + num_kpts, print_output, wann_control%use_ss_functional, wann_control%constrain, & omega%invariant, ln_tmp_loc, m_matrix_loc, lambda_loc, first_pass, timer, & nkrank, global_k, error, comm) if (allocated(error)) return @@ -548,7 +548,7 @@ subroutine wann_main(ham_logical, kmesh_info, kpt_latt, wann_control, omega, sit (iter .gt. wann_control%guiding_centres%num_no_guide_iter) & .and. (mod(iter, wann_control%guiding_centres%num_guide_cycles) .eq. 0)) then call wann_phases(csheet, sheet, rguide, irguide, num_wann, kmesh_info, num_kpts, & - wann_control%mv_functional, m_matrix_loc, rnkb, print_output%timing_level, & + wann_control%use_ss_functional, m_matrix_loc, rnkb, print_output%timing_level, & print_output%iprint, timer, nkrank, global_k, error, comm) if (allocated(error)) return @@ -558,7 +558,7 @@ subroutine wann_main(ham_logical, kmesh_info, kpt_latt, wann_control, omega, sit ! calculate gradient of omega if (lsitesymmetry .or. wann_control%precond) then call wann_domega(dist_k, csheet, sheet, rave, num_wann, kmesh_info, num_kpts, & - wann_control%constrain, wann_control%mv_functional, lsitesymmetry, & + wann_control%constrain, wann_control%use_ss_functional, lsitesymmetry, & ln_tmp_loc, m_matrix_loc, & rnkb_loc, cdodq_loc, lambda_loc, print_output%timing_level, sitesym, & timer, nkrank, global_k, error, comm, print_output%iprint, cdodq) @@ -566,7 +566,7 @@ subroutine wann_main(ham_logical, kmesh_info, kpt_latt, wann_control, omega, sit else call wann_domega(dist_k, csheet, sheet, rave, num_wann, kmesh_info, num_kpts, & - wann_control%constrain, wann_control%mv_functional, lsitesymmetry, & + wann_control%constrain, wann_control%use_ss_functional, lsitesymmetry, & ln_tmp_loc, m_matrix_loc, & rnkb_loc, cdodq_loc, lambda_loc, print_output%timing_level, sitesym, & timer, nkrank, global_k, error, comm, print_output%iprint) @@ -627,7 +627,7 @@ subroutine wann_main(ham_logical, kmesh_info, kpt_latt, wann_control, omega, sit ! calculate spread at trial step call wann_omega(csheet, sheet, rave, r2ave, rave2, trial_spread, num_wann, kmesh_info, & - num_kpts, print_output, wann_control%mv_functional, wann_control%constrain, & + num_kpts, print_output, wann_control%use_ss_functional, wann_control%constrain, & omega%invariant, ln_tmp_loc, m_matrix_loc, lambda_loc, first_pass, timer, & nkrank, global_k, error, comm) if (allocated(error)) return @@ -690,7 +690,7 @@ subroutine wann_main(ham_logical, kmesh_info, kpt_latt, wann_control, omega, sit ! calculate the new centers and spread call wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, kmesh_info, & - num_kpts, print_output, wann_control%mv_functional ,wann_control%constrain, & + num_kpts, print_output, wann_control%use_ss_functional, wann_control%constrain, & omega%invariant, ln_tmp_loc, m_matrix_loc, lambda_loc, first_pass, timer, & nkrank, global_k, error, comm) if (allocated(error)) return @@ -884,7 +884,7 @@ subroutine wann_main(ham_logical, kmesh_info, kpt_latt, wann_control, omega, sit if (wann_control%guiding_centres%enable) then call wann_phases(csheet, sheet, rguide, irguide, num_wann, kmesh_info, num_kpts, & - wann_control%mv_functional, m_matrix_loc, rnkb, print_output%timing_level, & + wann_control%use_ss_functional, m_matrix_loc, rnkb, print_output%timing_level, & print_output%iprint, timer, nkrank, global_k, error, comm) if (allocated(error)) return endif @@ -1756,7 +1756,7 @@ end subroutine wann_main !================================================! - subroutine wann_phases(csheet, sheet, rguide, irguide, num_wann, kmesh_info, num_kpts, mv_functional, & + subroutine wann_phases(csheet, sheet, rguide, irguide, num_wann, kmesh_info, num_kpts, use_ss_functional, & m_matrix_loc, rnkb, timing_level, iprint, timer, nkrank, global_k, error, & comm, m_w) !================================================! @@ -1765,7 +1765,7 @@ subroutine wann_phases(csheet, sheet, rguide, irguide, num_wann, kmesh_info, num ! !================================================ - use w90_constants, only: eps6, cmplx_0, cmplx_i + use w90_constants, only: eps6, cmplx_0, cmplx_i, twopi use w90_io, only: io_stopwatch_start, io_stopwatch_stop use w90_utility, only: utility_inv3 use w90_comms, only: comms_allreduce, w90_comm_type, mpirank @@ -1794,7 +1794,7 @@ subroutine wann_phases(csheet, sheet, rguide, irguide, num_wann, kmesh_info, num complex(kind=dp), intent(out) :: csheet(:, :, :) !! Choice of phase complex(kind=dp), intent(in) :: m_matrix_loc(:, :, :, :) - logical, intent(in) :: mv_functional + logical, intent(in) :: use_ss_functional !local complex(kind=dp) :: csum(kmesh_info%nnh) @@ -1928,24 +1928,32 @@ subroutine wann_phases(csheet, sheet, rguide, irguide, num_wann, kmesh_info, num ! obtain branch cut choice guided by rguid sheet = 0.0_dp - do nkp = 1, num_kpts + if (use_ss_functional) then do nn = 1, kmesh_info%nntot do loop_wann = 1, num_wann ! sheet (loop_wann, nn, nkp) = 0.d0 do j = 1, 3 - if (.NOT. mv_functional) then - sheet(loop_wann, nn, 1) = sheet(loop_wann, nn, 1) & - + kmesh_info%bk(j, nn, nkp)*rguide(j, loop_wann) - else - sheet(loop_wann, nn, nkp) = sheet(loop_wann, nn, nkp) & - + kmesh_info%bk(j, nn, nkp)*rguide(j, loop_wann) - endif + sheet(loop_wann, nn, 1) = kmesh_info%bk(j, nn, 1)*rguide(j, loop_wann) enddo ! csheet (loop_wann, nn, nkp) = exp (ci * sheet (loop_wann, nn, nkp) ) enddo enddo - enddo - csheet = exp(cmplx_i*sheet) + csheet = exp(cmplx_i*sheet) + else + do nkp = 1, num_kpts + do nn = 1, kmesh_info%nntot + do loop_wann = 1, num_wann + ! sheet (loop_wann, nn, nkp) = 0.d0 + do j = 1, 3 + sheet(loop_wann, nn, nkp) = sheet(loop_wann, nn, nkp) & + + kmesh_info%bk(j, nn, nkp)*rguide(j, loop_wann) + enddo + ! csheet (loop_wann, nn, nkp) = exp (ci * sheet (loop_wann, nn, nkp) ) + enddo + enddo + enddo + csheet = exp(cmplx_i*sheet) + endif ! now check that we picked the proper sheet for the log ! of m_matrix. criterion: q_n^{k,b}=Im(ln(M_nn^{k,b})) + b \cdot r_n are @@ -1986,7 +1994,7 @@ end subroutine wann_phases !================================================! subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, kmesh_info, & - num_kpts, print_output, mv_functional, wann_slwf, omega_invariant, ln_tmp_loc, & + num_kpts, print_output, use_ss_functional, wann_slwf, omega_invariant, ln_tmp_loc, & m_matrix_loc, lambda_loc, first_pass, timer, nkrank, global_k, error, comm) !================================================! ! @@ -2022,7 +2030,7 @@ subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, real(kind=dp), intent(in) :: lambda_loc real(kind=dp), intent(in) :: omega_invariant - logical, intent(in) :: mv_functional + logical, intent(in) :: use_ss_functional real(kind=dp), intent(inout) :: ln_tmp_loc(:, :, :) real(kind=dp), intent(in) :: sheet(:, :, :) real(kind=dp), intent(out) :: r2ave(:) @@ -2044,20 +2052,20 @@ subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, if (print_output%timing_level > 1 .and. print_output%iprint > 0) call io_stopwatch_start('wann: omega', timer) - if (.NOT. mv_functional) then + if (use_ss_functional) then allocate (sum_mnn(num_wann, kmesh_info%nntot), stat=ierr) if (ierr /= 0) then call set_error_alloc(error, 'Error in allocating sum_mnn in wann_omega', comm) return endif - + sum_mnn = 0.0_dp do nn = 1, kmesh_info%nntot do n = 1, num_wann do nkp_loc = 1, nkrank nkp = global_k(nkp_loc) ! Note that this ln_tmp is defined differently wrt the one in wann_domega - sum_mnn(n, nn) = sum_mnn(n, nn) + csheet(n, nn, 1) * m_matrix_loc(n, n, nn, nkp_loc) + sum_mnn(n, nn) = sum_mnn(n, nn) + csheet(n, nn, 1)*m_matrix_loc(n, n, nn, nkp_loc) enddo enddo enddo @@ -2065,10 +2073,10 @@ subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, call comms_allreduce(sum_mnn(1, 1), num_wann*kmesh_info%nntot, 'SUM', error, comm) if (allocated(error)) return - sum_mnn = sum_mnn / real(num_kpts, dp) - + sum_mnn = sum_mnn/real(num_kpts, dp) + do nkp_loc = 1, nkrank - ln_tmp_loc(:, :, nkp_loc) = aimag(log(sum_mnn(:, :))) - sheet(:, :, 1) + ln_tmp_loc(:, :, nkp_loc) = aimag(log(sum_mnn(:, :))) - sheet(:, :, 1) enddo rave = 0.0_dp @@ -2087,13 +2095,13 @@ subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, call comms_allreduce(rave(1, 1), 3*num_wann, 'SUM', error, comm) if (allocated(error)) return - rave = -rave / real(num_kpts, dp) + rave = -rave/real(num_kpts, dp) rave2 = 0.0_dp do iw = 1, num_wann rave2(iw) = sum(rave(:, iw)*rave(:, iw)) enddo - + r2ave = 0.0_dp do nn = 1, kmesh_info%nntot r2ave(:) = r2ave(:) + kmesh_info%wb(nn)*(1.0_dp - abs(sum_mnn(:, nn))**2) @@ -2106,7 +2114,7 @@ subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, call set_error_dealloc(error, 'Error in deallocating sum_mnn in wann_omega', comm) return endif - else + else do nkp_loc = 1, nkrank nkp = global_k(nkp_loc) do nn = 1, kmesh_info%nntot @@ -2117,7 +2125,7 @@ subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, enddo enddo enddo - + rave = 0.0_dp do iw = 1, num_wann do ind = 1, 3 @@ -2141,7 +2149,7 @@ subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, rave2(iw) = sum(rave(:, iw)*rave(:, iw)) enddo - ! aam: is this useful? + ! aam: is this useful? !~ rtot=0.0_dp !~ do ind = 1, 3 !~ do loop_wann = 1, num_wann @@ -2154,7 +2162,7 @@ subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, do nkp_loc = 1, nkrank do nn = 1, kmesh_info%nntot mnn2 = real(m_matrix_loc(iw, iw, nn, nkp_loc) & - *conjg(m_matrix_loc(iw, iw, nn, nkp_loc)), kind=dp) + *conjg(m_matrix_loc(iw, iw, nn, nkp_loc)), kind=dp) r2ave(iw) = r2ave(iw) + kmesh_info%wb(nn)*(1.0_dp - mnn2 + ln_tmp_loc(iw, nn, nkp_loc)**2) enddo enddo @@ -2334,7 +2342,7 @@ subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, wann_spread%om_od = wann_spread%om_od/real(num_kpts, dp) - if (.NOT. mv_functional) then + if (use_ss_functional) then wann_spread%om_d = 0.0_dp do nn = 1, kmesh_info%nntot do n = 1, num_wann @@ -2346,7 +2354,7 @@ subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, call comms_allreduce(summ, 1, 'SUM', error, comm) if (allocated(error)) return - summ = summ / real(num_kpts, dp) + summ = summ/real(num_kpts, dp) wann_spread%om_d = wann_spread%om_d - kmesh_info%wb(nn)*abs(summ)**2 @@ -2358,7 +2366,7 @@ subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, call comms_allreduce(summ, 1, 'SUM', error, comm) if (allocated(error)) return - summ = summ / real(num_kpts, dp) + summ = summ/real(num_kpts, dp) wann_spread%om_d = wann_spread%om_d + kmesh_info%wb(nn)*summ enddo @@ -2375,10 +2383,10 @@ subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, enddo enddo enddo - + call comms_allreduce(wann_spread%om_d, 1, 'SUM', error, comm) if (allocated(error)) return - + wann_spread%om_d = wann_spread%om_d/real(num_kpts, dp) endif @@ -2392,7 +2400,7 @@ subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread, num_wann, end subroutine wann_omega !================================================! - subroutine wann_domega(dist_k, csheet, sheet, rave, num_wann, kmesh_info, num_kpts, wann_slwf, mv_functional, & + subroutine wann_domega(dist_k, csheet, sheet, rave, num_wann, kmesh_info, num_kpts, wann_slwf, use_ss_functional, & lsitesymmetry, ln_tmp_loc, m_matrix_loc, rnkb_loc, cdodq_loc, & lambda_loc, timing_level, sitesym, timer, nkrank, global_k, error, comm, & iprint, cdodq) @@ -2428,7 +2436,7 @@ subroutine wann_domega(dist_k, csheet, sheet, rave, num_wann, kmesh_info, num_kp integer, intent(in) :: timing_level, iprint integer, intent(in) :: nkrank integer, intent(in) :: global_k(:) - logical, intent(in) :: mv_functional + logical, intent(in) :: use_ss_functional real(kind=dp), intent(in) :: sheet(:, :, :) real(kind=dp), intent(out) :: rave(:, :) @@ -2484,7 +2492,7 @@ subroutine wann_domega(dist_k, csheet, sheet, rave, num_wann, kmesh_info, num_kp endif end if - if (.NOT. mv_functional) then + if (use_ss_functional) then allocate (sum_mnn(num_wann, kmesh_info%nntot), stat=ierr) if (ierr /= 0) then call set_error_alloc(error, 'Error in allocating sum_mnn in wann_omega', comm) @@ -2497,7 +2505,7 @@ subroutine wann_domega(dist_k, csheet, sheet, rave, num_wann, kmesh_info, num_kp do nkp_loc = 1, nkrank nkp = global_k(nkp_loc) ! Note that this ln_tmp is defined differently wrt the one in wann_domega - sum_mnn(n, nn) = sum_mnn(n, nn) + csheet(n, nn, 1) * m_matrix_loc(n, n, nn, nkp_loc) + sum_mnn(n, nn) = sum_mnn(n, nn) + csheet(n, nn, 1)*m_matrix_loc(n, n, nn, nkp_loc) enddo enddo enddo @@ -2505,10 +2513,10 @@ subroutine wann_domega(dist_k, csheet, sheet, rave, num_wann, kmesh_info, num_kp call comms_allreduce(sum_mnn(1, 1), num_wann*kmesh_info%nntot, 'SUM', error, comm) if (allocated(error)) return - sum_mnn = sum_mnn / real(num_kpts, dp) + sum_mnn = sum_mnn/real(num_kpts, dp) do nkp_loc = 1, nkrank - ln_tmp_loc(:, :, nkp_loc) = kmesh_info%wb(nn)*(aimag(log(sum_mnn(:, :))) - sheet(:, :, 1)) + ln_tmp_loc(:, :, nkp_loc) = kmesh_info%wb(nn)*(aimag(log(sum_mnn(:, :))) - sheet(:, :, 1)) enddo rave = 0.0_dp @@ -2517,7 +2525,7 @@ subroutine wann_domega(dist_k, csheet, sheet, rave, num_wann, kmesh_info, num_kp do nkp_loc = 1, nkrank nkp = global_k(nkp_loc) do nn = 1, kmesh_info%nntot - rave(ind, iw) = rave(ind, iw) + kmesh_info%bk(ind, nn, nkp) * ln_tmp_loc(iw, nn, nkp_loc) + rave(ind, iw) = rave(ind, iw) + kmesh_info%bk(ind, nn, nkp)*ln_tmp_loc(iw, nn, nkp_loc) enddo enddo enddo @@ -2538,31 +2546,45 @@ subroutine wann_domega(dist_k, csheet, sheet, rave, num_wann, kmesh_info, num_kp enddo enddo + sum_mnn = 0.0_dp + do n = 1, num_wann + do nn = 1, kmesh_info%nntot + do nkp_loc = 1, nkrank + nkp = global_k(nkp_loc) + ! Note that this ln_tmp is defined differently wrt the one in wann_domega + sum_mnn(n, nn) = sum_mnn(n, nn) + m_matrix_loc(n, n, nn, nkp_loc) + enddo + enddo + enddo + + call comms_allreduce(sum_mnn(1, 1), num_wann*kmesh_info%nntot, 'SUM', error, comm) + if (allocated(error)) return + cdodq_loc = cmplx_0 do nkp_loc = 1, nkrank do n = 1, num_wann do m = 1, num_wann do nn = 1, kmesh_info%nntot - nn2 = modulo(nn-1+kmesh_info%nntot/2, kmesh_info%nntot)+1 + nn2 = modulo(nn - 1 + kmesh_info%nntot/2, kmesh_info%nntot) + 1 cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) + & - kmesh_info%wb(nn) * m_matrix_loc(m, n, nn, nkp_loc) * & - conjg(sum_mnn(n, nn)) + kmesh_info%wb(nn)*m_matrix_loc(m, n, nn, nkp_loc)* & + conjg(sum_mnn(n, nn)) cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) - & - kmesh_info%wb(nn) * conjg(m_matrix_loc(n, m, nn2, nkp_loc)) * & - conjg(sum_mnn(m, nn)) + kmesh_info%wb(nn)*conjg(m_matrix_loc(n, m, nn2, nkp_loc))* & + conjg(sum_mnn(m, nn)) cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) - & - kmesh_info%wb(nn) * conjg(m_matrix_loc(n, m, nn, nkp_loc)) * & - sum_mnn(m, nn) + kmesh_info%wb(nn)*conjg(m_matrix_loc(n, m, nn, nkp_loc))* & + sum_mnn(m, nn) cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) + & - kmesh_info%wb(nn) * m_matrix_loc(m, n, nn2, nkp_loc) * & - sum_mnn(n, nn) + kmesh_info%wb(nn)*m_matrix_loc(m, n, nn2, nkp_loc)* & + sum_mnn(n, nn) enddo enddo enddo enddo - cdodq_loc = cdodq_loc / real(num_kpts, dp) - + cdodq_loc = cdodq_loc/real(num_kpts, dp) + if (present(cdodq)) then ! each process communicates its result to other processes cdodq(:, :, :) = 0.0_dp @@ -2611,7 +2633,7 @@ subroutine wann_domega(dist_k, csheet, sheet, rave, num_wann, kmesh_info, num_kp do nn = 1, kmesh_info%nntot do n = 1, num_wann r0kb(n, nn, nkp_loc) = sum(kmesh_info%bk(:, nn, nkp) & - *wann_slwf%centres(n, :)) + *wann_slwf%centres(n, :)) enddo enddo enddo @@ -2646,72 +2668,72 @@ subroutine wann_domega(dist_k, csheet, sheet, rave, num_wann, kmesh_info, num_kp if (n <= wann_slwf%slwf_num) then ! A[R^{k,b}]=(R-Rdag)/2 cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) & - + kmesh_info%wb(nn)*0.5_dp*(cr(m, n) - conjg(cr(n, m))) + + kmesh_info%wb(nn)*0.5_dp*(cr(m, n) - conjg(cr(n, m))) cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) & - - (crt(m, n)*ln_tmp_loc(n, nn, nkp_loc) & + - (crt(m, n)*ln_tmp_loc(n, nn, nkp_loc) & + conjg(crt(n, m)*ln_tmp_loc(m, nn, nkp_loc))) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) + *cmplx(0.0_dp, -0.5_dp, kind=dp) cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) & - - (crt(m, n)*rnkb_loc(n, nn, nkp_loc) & + - (crt(m, n)*rnkb_loc(n, nn, nkp_loc) & + conjg(crt(n, m)*rnkb_loc(m, nn, nkp_loc))) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) + *cmplx(0.0_dp, -0.5_dp, kind=dp) if (wann_slwf%constrain) then cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) + lambda_loc & - *(crt(m, n)*ln_tmp_loc(n, nn, nkp_loc) & - + conjg(crt(n, m)*ln_tmp_loc(m, nn, nkp_loc))) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) + *(crt(m, n)*ln_tmp_loc(n, nn, nkp_loc) & + + conjg(crt(n, m)*ln_tmp_loc(m, nn, nkp_loc))) & + *cmplx(0.0_dp, -0.5_dp, kind=dp) cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) & - + kmesh_info%wb(nn)*lambda_loc & - *(crt(m, n)*rnkb_loc(n, nn, nkp_loc) & - + conjg(crt(n, m)*rnkb_loc(m, nn, nkp_loc))) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) + + kmesh_info%wb(nn)*lambda_loc & + *(crt(m, n)*rnkb_loc(n, nn, nkp_loc) & + + conjg(crt(n, m)*rnkb_loc(m, nn, nkp_loc))) & + *cmplx(0.0_dp, -0.5_dp, kind=dp) cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) - lambda_loc & - *(crt(m, n)*ln_tmp_loc(n, nn, nkp_loc) & - + conjg(crt(n, m))*ln_tmp_loc(m, nn, nkp_loc)) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) + *(crt(m, n)*ln_tmp_loc(n, nn, nkp_loc) & + + conjg(crt(n, m))*ln_tmp_loc(m, nn, nkp_loc)) & + *cmplx(0.0_dp, -0.5_dp, kind=dp) cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) & - - kmesh_info%wb(nn)*lambda_loc & - *(r0kb(n, nn, nkp_loc)*crt(m, n) & - + r0kb(m, nn, nkp_loc)*conjg(crt(n, m))) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) + - kmesh_info%wb(nn)*lambda_loc & + *(r0kb(n, nn, nkp_loc)*crt(m, n) & + + r0kb(m, nn, nkp_loc)*conjg(crt(n, m))) & + *cmplx(0.0_dp, -0.5_dp, kind=dp) end if else cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) - kmesh_info%wb(nn) & - *0.5_dp*conjg(cr(n, m)) & - - conjg(crt(n, m)*(ln_tmp_loc(m, nn, nkp_loc) & + *0.5_dp*conjg(cr(n, m)) & + - conjg(crt(n, m)*(ln_tmp_loc(m, nn, nkp_loc) & + kmesh_info%wb(nn)*rnkb_loc(m, nn, nkp_loc))) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) + *cmplx(0.0_dp, -0.5_dp, kind=dp) if (wann_slwf%constrain) then cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) + lambda_loc & - *conjg(crt(n, m)*(ln_tmp_loc(m, nn, nkp_loc) & - + kmesh_info%wb(nn)*rnkb_loc(m, nn, nkp_loc))) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) & - - lambda_loc*(conjg(crt(n, m)) & - *ln_tmp_loc(m, nn, nkp_loc)) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) + *conjg(crt(n, m)*(ln_tmp_loc(m, nn, nkp_loc) & + + kmesh_info%wb(nn)*rnkb_loc(m, nn, nkp_loc))) & + *cmplx(0.0_dp, -0.5_dp, kind=dp) & + - lambda_loc*(conjg(crt(n, m)) & + *ln_tmp_loc(m, nn, nkp_loc)) & + *cmplx(0.0_dp, -0.5_dp, kind=dp) cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) & - - kmesh_info%wb(nn)*lambda_loc & - *r0kb(m, nn, nkp_loc)*conjg(crt(n, m)) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) + - kmesh_info%wb(nn)*lambda_loc & + *r0kb(m, nn, nkp_loc)*conjg(crt(n, m)) & + *cmplx(0.0_dp, -0.5_dp, kind=dp) end if end if else if (n <= wann_slwf%slwf_num) then cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) & - + kmesh_info%wb(nn)*cr(m, n)*0.5_dp & - - crt(m, n)*(ln_tmp_loc(n, nn, nkp_loc) & + + kmesh_info%wb(nn)*cr(m, n)*0.5_dp & + - crt(m, n)*(ln_tmp_loc(n, nn, nkp_loc) & + kmesh_info%wb(nn)*rnkb_loc(n, nn, nkp_loc)) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) + *cmplx(0.0_dp, -0.5_dp, kind=dp) if (wann_slwf%constrain) then cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) + lambda_loc & - *crt(m, n)*(ln_tmp_loc(n, nn, nkp_loc) & - + kmesh_info%wb(nn)*rnkb_loc(n, nn, nkp_loc)) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) & - - lambda_loc*crt(m, n)*ln_tmp_loc(n, nn, nkp_loc) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) + *crt(m, n)*(ln_tmp_loc(n, nn, nkp_loc) & + + kmesh_info%wb(nn)*rnkb_loc(n, nn, nkp_loc)) & + *cmplx(0.0_dp, -0.5_dp, kind=dp) & + - lambda_loc*crt(m, n)*ln_tmp_loc(n, nn, nkp_loc) & + *cmplx(0.0_dp, -0.5_dp, kind=dp) cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) - kmesh_info%wb(nn) & - *lambda_loc & - *r0kb(n, nn, nkp_loc)*crt(m, n) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) + *lambda_loc & + *r0kb(n, nn, nkp_loc)*crt(m, n) & + *cmplx(0.0_dp, -0.5_dp, kind=dp) end if else cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) @@ -2723,17 +2745,17 @@ subroutine wann_domega(dist_k, csheet, sheet, rave, num_wann, kmesh_info, num_kp do m = 1, num_wann ! A[R^{k,b}]=(R-Rdag)/2 cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) & - + kmesh_info%wb(nn)*0.5_dp & - *(cr(m, n) - conjg(cr(n, m))) + + kmesh_info%wb(nn)*0.5_dp & + *(cr(m, n) - conjg(cr(n, m))) ! -S[T^{k,b}]=-(T+Tdag)/2i ; T_mn = Rt_mn q_n cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) - & - (crt(m, n)*ln_tmp_loc(n, nn, nkp_loc) & + (crt(m, n)*ln_tmp_loc(n, nn, nkp_loc) & + conjg(crt(n, m)*ln_tmp_loc(m, nn, nkp_loc))) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) + *cmplx(0.0_dp, -0.5_dp, kind=dp) cdodq_loc(m, n, nkp_loc) = cdodq_loc(m, n, nkp_loc) - kmesh_info%wb(nn) & - *(crt(m, n)*rnkb_loc(n, nn, nkp_loc) & - + conjg(crt(n, m)*rnkb_loc(m, nn, nkp_loc))) & - *cmplx(0.0_dp, -0.5_dp, kind=dp) + *(crt(m, n)*rnkb_loc(n, nn, nkp_loc) & + + conjg(crt(n, m)*rnkb_loc(m, nn, nkp_loc))) & + *cmplx(0.0_dp, -0.5_dp, kind=dp) enddo enddo endif @@ -3083,7 +3105,7 @@ subroutine wann_main_gamma(kmesh_info, wann_control, omega, print_output, wannie irguide = 0 if (wann_control%guiding_centres%enable .and. (wann_control%guiding_centres%num_no_guide_iter .le. 0)) then call wann_phases(csheet, sheet, rguide, irguide, num_wann, kmesh_info, num_kpts, & - wann_control%mv_functional, m_matrix, rnkb, print_output%timing_level, & + wann_control%use_ss_functional, m_matrix, rnkb, print_output%timing_level, & print_output%iprint, timer, num_kpts, global_k, error, comm) ! no plellisation so num_kpts_local = num_kpts if (allocated(error)) return irguide = 1 @@ -3175,7 +3197,7 @@ subroutine wann_main_gamma(kmesh_info, wann_control, omega, print_output, wannie (iter .gt. wann_control%guiding_centres%num_no_guide_iter) & .and. (mod(iter, wann_control%guiding_centres%num_guide_cycles) .eq. 0)) then call wann_phases(csheet, sheet, rguide, irguide, num_wann, kmesh_info, num_kpts, & - wann_control%mv_functional, m_matrix, rnkb, print_output%timing_level, & + wann_control%use_ss_functional, m_matrix, rnkb, print_output%timing_level, & print_output%iprint, timer, num_kpts, global_k, error, comm, m_w) ! num_kpts_loc == num_kpts here if (allocated(error)) return irguide = 1 @@ -3286,7 +3308,7 @@ subroutine wann_main_gamma(kmesh_info, wann_control, omega, print_output, wannie if (wann_control%guiding_centres%enable) then call wann_phases(csheet, sheet, rguide, irguide, num_wann, kmesh_info, num_kpts, & - wann_control%mv_functional, m_matrix, rnkb, print_output%timing_level, & + wann_control%use_ss_functional, m_matrix, rnkb, print_output%timing_level, & print_output%iprint, timer, num_kpts, global_k, error, comm) ! num_kpts_loc == num_kpts here if (allocated(error)) return endif diff --git a/test-suite/tests/testw90_example34_2/silicon.win b/test-suite/tests/testw90_example34_2/silicon.win index 2e354ac8..3b49775d 100644 --- a/test-suite/tests/testw90_example34_2/silicon.win +++ b/test-suite/tests/testw90_example34_2/silicon.win @@ -11,7 +11,7 @@ num_print_cycles = 10 length_unit = bohr -mv_functional = .false. +use_ss_functional = .true. begin projections !! !! Bond-centred s-orbitals diff --git a/tutorials/tutorial34/Marzari-Vanderbilt/silicon.win b/tutorials/tutorial34/Marzari-Vanderbilt/silicon.win index 5c8dea9d..55e9e9dc 100644 --- a/tutorials/tutorial34/Marzari-Vanderbilt/silicon.win +++ b/tutorials/tutorial34/Marzari-Vanderbilt/silicon.win @@ -9,7 +9,7 @@ num_print_cycles = 10 length_unit = bohr -mv_functional = .true. +use_ss_functional = .false. begin projections !! !! Bond-centred s-orbitals diff --git a/tutorials/tutorial34/Stengel-Spaldin/silicon.win b/tutorials/tutorial34/Stengel-Spaldin/silicon.win index 253361b0..ffc6931c 100644 --- a/tutorials/tutorial34/Stengel-Spaldin/silicon.win +++ b/tutorials/tutorial34/Stengel-Spaldin/silicon.win @@ -9,7 +9,7 @@ num_print_cycles = 10 length_unit = bohr -mv_functional = .false. +use_ss_functional = .true. begin projections !! !! Bond-centred s-orbitals