Find the B1 weights for a set of shells specified by the user
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| type(kmesh_input_type), | intent(in) | :: | kmesh_input | |||
| type(print_output_type), | intent(in) | :: | print_output | |||
| real(kind=dp), | intent(out) | :: | bweight(max_shells) | |||
| real(kind=dp), | intent(in) | :: | dnn(kmesh_input%search_shells) | |||
| real(kind=dp), | intent(in) | :: | kpt_cart(:,:) | |||
| real(kind=dp), | intent(in) | :: | recip_lattice(3,3) | |||
| integer, | intent(in) | :: | lmn(:,:) | |||
| integer, | intent(in) | :: | multi(kmesh_input%search_shells) | |||
| integer, | intent(in) | :: | num_kpts | |||
| integer, | intent(in) | :: | stdout | |||
| type(timer_list_type), | intent(inout) | :: | timer | |||
| type(w90_error_type), | intent(out), | allocatable | :: | error | ||
| type(w90_comm_type), | intent(in) | :: | comm |
subroutine kmesh_shell_fixed(kmesh_input, print_output, bweight, dnn, kpt_cart, recip_lattice, & lmn, multi, num_kpts, stdout, timer, error, comm) !================================================ ! !! Find the B1 weights for a set of shells specified by the user ! !================================================ use w90_constants, only: eps7 use w90_io, only: io_stopwatch_start, io_stopwatch_stop use w90_types, only: kmesh_input_type, print_output_type, timer_list_type implicit none ! arguments type(print_output_type), intent(in) :: print_output type(kmesh_input_type), intent(in) :: kmesh_input type(timer_list_type), intent(inout) :: timer type(w90_error_type), allocatable, intent(out) :: error type(w90_comm_type), intent(in) :: comm integer, intent(in) :: num_kpts integer, intent(in) :: stdout integer, intent(in) :: lmn(:, :) integer, intent(in) :: multi(kmesh_input%search_shells) ! the number of kpoints in the shell real(kind=dp), intent(in) :: recip_lattice(3, 3) real(kind=dp), intent(in) ::kpt_cart(:, :) real(kind=dp), intent(in) :: dnn(kmesh_input%search_shells) ! the bvectors real(kind=dp), intent(out) :: bweight(max_shells) ! local variables real(kind=dp), allocatable :: bvector(:, :, :) real(kind=dp) :: singv(kmesh_input%num_shells) real(kind=dp) :: amat(max_shells, kmesh_input%num_shells) real(kind=dp) :: umat(max_shells, max_shells) real(kind=dp) :: vmat(kmesh_input%num_shells, kmesh_input%num_shells) real(kind=dp) :: smat(kmesh_input%num_shells, max_shells) integer, parameter :: lwork = max_shells*10 real(kind=dp) :: work(lwork) real(kind=dp), parameter :: target(6) = (/1.0_dp, 1.0_dp, 1.0_dp, 0.0_dp, 0.0_dp, 0.0_dp/) logical :: b1sat integer :: ierr, loop_i, loop_j, loop_b, loop_s, info real(kind=dp) :: delta integer :: loop, shell if (print_output%timing_level > 1) call io_stopwatch_start('kmesh: shell_fixed', timer) allocate (bvector(3, maxval(multi), kmesh_input%num_shells), stat=ierr) if (ierr /= 0) then call set_error_alloc(error, 'Error allocating bvector in kmesh_shell_fixed', comm) return end if bvector = 0.0_dp; bweight = 0.0_dp amat = 0.0_dp; umat = 0.0_dp; vmat = 0.0_dp; smat = 0.0_dp; singv = 0.0_dp if (print_output%iprint > 0) then write (stdout, '(1x,a)') '| The b-vectors are set in the win file |' end if do shell = 1, kmesh_input%num_shells ! get the b vectors for this shell call kmesh_get_bvectors(kmesh_input, print_output, & bvector(:, 1:multi(kmesh_input%shell_list(shell)), shell), kpt_cart, & recip_lattice, dnn(kmesh_input%shell_list(shell)), lmn, 1, & multi(kmesh_input%shell_list(shell)), num_kpts, timer, error, comm) if (allocated(error)) return end do if (print_output%iprint >= 3) then do shell = 1, kmesh_input%num_shells write (stdout, '(1x,a8,1x,I2,a14,1x,I2,49x,a)') '| Shell:', shell, ' Multiplicity:', & multi(kmesh_input%shell_list(shell)), '|' do loop = 1, multi(kmesh_input%shell_list(shell)) write (stdout, '(1x,a10,I2,1x,a1,4x,3f12.6,5x,a9,9x,a)') '| b-vector ', loop, ':', & bvector(:, loop, shell)/print_output%lenconfac, '('//trim(print_output%length_unit)//'^-1)', '|' end do end do end if do loop_s = 1, kmesh_input%num_shells do loop_b = 1, multi(kmesh_input%shell_list(loop_s)) amat(1, loop_s) = amat(1, loop_s) + bvector(1, loop_b, loop_s)*bvector(1, loop_b, loop_s) amat(2, loop_s) = amat(2, loop_s) + bvector(2, loop_b, loop_s)*bvector(2, loop_b, loop_s) amat(3, loop_s) = amat(3, loop_s) + bvector(3, loop_b, loop_s)*bvector(3, loop_b, loop_s) amat(4, loop_s) = amat(4, loop_s) + bvector(1, loop_b, loop_s)*bvector(2, loop_b, loop_s) amat(5, loop_s) = amat(5, loop_s) + bvector(2, loop_b, loop_s)*bvector(3, loop_b, loop_s) amat(6, loop_s) = amat(6, loop_s) + bvector(3, loop_b, loop_s)*bvector(1, loop_b, loop_s) end do end do info = 0 call dgesvd('A', 'A', max_shells, kmesh_input%num_shells, amat, max_shells, singv, umat, & max_shells, vmat, kmesh_input%num_shells, work, lwork, info) if (info < 0) then if (print_output%iprint > 0) then write (stdout, '(1x,a,1x,I1,1x,a)') 'kmesh_shell_fixed: Argument', abs(info), & 'of dgesvd is incorrect' end if call set_error_fatal(error, 'kmesh_shell_fixed: Problem with Singular Value Decomposition', comm) return else if (info > 0) then call set_error_fatal(error, 'kmesh_shell_fixed: Singular Value Decomposition did not converge', comm) return end if if (any(abs(singv) < eps7)) then call set_error_fatal(error, 'kmesh_shell_fixed: Singular Value Decomposition has found a very small singular value', comm) return end if smat = 0.0_dp do loop_s = 1, kmesh_input%num_shells smat(loop_s, loop_s) = 1/singv(loop_s) end do bweight(1:kmesh_input%num_shells) = matmul(transpose(vmat), matmul(smat, matmul(transpose(umat), target))) if (print_output%iprint >= 2) then do loop_s = 1, kmesh_input%num_shells write (stdout, '(1x,a,I2,a,f12.7,5x,a8,36x,a)') '| Shell: ', loop_s, & ' w_b ', bweight(loop_s)*print_output%lenconfac**2, '('//trim(print_output%length_unit)//'^2)', '|' end do end if !check b1 b1sat = .true. if (.not. kmesh_input%skip_B1_tests) then do loop_i = 1, 3 do loop_j = loop_i, 3 delta = 0.0_dp do loop_s = 1, kmesh_input%num_shells do loop_b = 1, multi(kmesh_input%shell_list(loop_s)) delta = delta + bweight(loop_s)*bvector(loop_i, loop_b, loop_s)*bvector(loop_j, loop_b, loop_s) end do end do if (loop_i == loop_j) then if (abs(delta - 1.0_dp) > kmesh_input%tol) b1sat = .false. end if if (loop_i /= loop_j) then if (abs(delta) > kmesh_input%tol) b1sat = .false. end if end do end do end if if (.not. b1sat) then call set_error_fatal(error, 'kmesh_shell_fixed: B1 condition not satisfied', comm) return end if if (print_output%timing_level > 1) call io_stopwatch_stop('kmesh: shell_fixed', timer) return end subroutine kmesh_shell_fixed