1! --- 2! Copyright (C) 1996-2016 The SIESTA group 3! This file is distributed under the terms of the 4! GNU General Public License: see COPYING in the top directory 5! or http://www.gnu.org/copyleft/gpl.txt . 6! See Docs/Contributors.txt for a list of contributors. 7! --- 8subroutine compute_pw_matrix( nncount, bvectorsfrac ) 9! 10! In this subroutine we compute the matrix elements of the plane wave, 11! for all the wave vectors that connect a given k-point to its nearest 12! neighbours 13! 14! ------------------------------ INPUT ----------------------------------------- 15! integer nncount : The number of nearest neighbours belonging to 16! each k-point of the Monkhorst-Pack mesh 17! real(dp) bvectorsfrac(:,:) : Vectors that connect each mesh k-point 18! to its nearest neighbours. 19! ------------------------------ OUTPUT ---------------------------------------- 20! Matrix elements of the plane wave for all the required wave vectors 21! connecting a given k-point to its nearest neighbours 22! This array (delkmatgen) is defined in the module m_planewavematrixvar, 23! file m_planewavematrixvar.F90 24! ------------------------------ BEHAVIOUR ------------------------------------- 25! First of all, we allocate the arrays where the matrix elements of the 26! plane waves will be stored 27! 28! Loop on all the neighbour k-points 29! 30! For each neighbour, compute the vector that connects each mesh point 31! to its nearest neighbour in bohr^-1. 32! Remember that bvectorsfrac are read in reduced units, so we have 33! to multiply then by the reciprocal lattice vector. 34! This is done in the function getkvector. 35! 36! Call to the subroutine that computes the matrix elements of the plane 37! wave in the grid (planewavematrix) 38! 39! Store the matrix elements of the plane wave for a given wave vector 40! (array delkmat) 41! in a permanent array, 42! (array delkmatgen) 43! So delkmat can be rewritten for the next wave vector. 44! 45! ------------------------------ UNITS ----------------------------------------- 46! bvectorsfrac are given in reduced units. 47! ------------------------------------------------------------------------------ 48 49! Used module procedures 50 use alloc, only: re_alloc ! Re-allocation routine 51 use m_planewavematrix, only: planewavematrix ! To compute the plane wave 52 ! matrix elements in the grid 53 54! Used module variables 55 use precision, only: dp ! Real douple precision type 56 use m_planewavematrixvar, only: delkmat ! Matrix elements of a plane wave 57 use m_planewavematrixvar, only: delkmatgen ! Matrix elements of a plane wave 58 ! (saved) 59 use sparse_matrices, only: maxnh ! Maximum number of orbitals 60 ! overlapping 61 use parallel, only: IOnode ! Input/Output node 62 63!! For debugging 64! use m_writedelk, only: write_delk 65! use atomlist, only: no_u, no_s, iaorb, iphorb 66! use atomlist, only: indxuo, qtot 67! use sparse_matrices, only: maxnh, listh, listhptr, numh, xijo, S 68! use m_spin, only: nspin 69! use siesta_options, only: temp 70! use sys, only: die 71!! End debugging 72 73 implicit none 74 75! The number of nearest neighbours belonging to 76! each k-point of the Monkhorst-Pack mesh 77 integer, intent(in) :: nncount 78! The vectors b that connect each mesh-point k to its nearest neighbours 79 real(dp), intent(in) :: bvectorsfrac(3,nncount) 80 81! 82! Internal variables 83! 84 integer :: inn 85 real(dp) :: bvector(3) 86 87!! For debugging 88! logical :: onlygamma 89!! End debugging 90 91 if( IOnode ) & 92 & write(6,'(/,a)') & 93 & 'compute_pw_matrix: Computing the matrix elements of a plane wave' 94 95! First of all, we allocate the arrays where the matrix elements of the 96! plane waves will be stored 97 if ( .not. associated(delkmatgen) ) then 98 nullify(delkmatgen) 99 call re_alloc(delkmatgen, 1, nncount, 1, maxnh, & 100 & name='delkmatgen', routine='compute_pw_matrix') 101 endif 102 if ( .not. associated(delkmat) ) then 103 nullify(delkmat) 104 call re_alloc(delkmat, 1, maxnh, & 105 & name='delkmat',routine='compute_pw_matrix') 106 endif 107 108! Loop on neighbour k-points 109 do inn = 1, nncount 110 111! For each neighbour, compute the vector that connects each mesh point 112! to its nearest neighbour in bohr^-1 (done in the subroutine getkvector). 113! Remember that bvectorsfrac are read in reduced units, so we have 114! to multiply then by the reciprocal lattice vector. 115 call getkvector( bvectorsfrac(:,inn), bvector ) 116 117!!$ if( IOnode ) then 118!!$ write(6,'(a)') & 119!!$ & 'compute_pw_matrix: Vector connecting a k-point with its neighbours' 120!!$ write(6,'(a,3f12.5)') & 121!!$ & 'compute_pw_matrix:', bvector(:) 122!!$ endif 123 124! Call to the subroutine that computes the matrix elements of the plane 125! wave in the grid (planewavematrix) 126 call planewavematrix( -1, bvector ) 127 128! Store the matrix elements of the plane wave for a given wave vector 129! (array delkmat) 130! in a permanent array, 131! (array delkmatgen) 132! So delkmat can be rewritten for the next wave vector. 133 delkmatgen(inn,:) = delkmat(:) 134 135!! For debugging 136! onlygamma = .false. 137! call write_delk( onlygamma, no_u, no_s, nspin, indxuo, maxnh, & 138! & numh, listhptr, listh, delkmat, S, qtot, temp, xijo) 139! call die() 140!! End debugging 141 142 enddo 143 144end subroutine compute_pw_matrix 145 146subroutine getkvector(frac, vector) 147! 148! Converts from fractional to Bohr^-1 coordinates in reciprocal space. 149! 150 151use m_siesta2wannier90, only : reclatvec ! Reciprocal lattice vectors 152use precision, only : dp ! Real double precision type 153 154integer :: i,j 155real(dp),intent(in) :: frac(3) 156real(dp),intent(out) :: vector(3) 157 158 159 do i = 1, 3 160 vector(i) = 0.0_dp 161 do j = 1, 3 162 vector(i) = frac(j)*reclatvec(i,j) + vector(i) 163 enddo 164 enddo 165 166end subroutine getkvector 167 168