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! --- 8module m_ts_kpoints 9! 10! Routines that are related to TS kpoint sampling 11! 12!============================================================================== 13! CONTAINS: 14! 1) setup_ts_scf_kscell 15! 2) setup_ts_kpoint_grid 16! 3) ts_write_k_points 17! 4) ts_iokp 18 19 use precision, only : dp 20 21 implicit none 22 23 private 24 save 25 26!===================== K-POINT RELATED VARIABLES ============================== 27!============================================================================== 28! Contains data structures and routines to deal with the kpoint-grid 29! for the self-consistent calculation with Green Functions. 30! Original verison modified so that the kpoint sampling is done for the 31! direction perpendicular to the transport directins 32! Version created to integrate TRANSIESTA in siesta2.3 33!============================================================================== 34! The kpoint mesh parameters that are public may be accessed in other 35! parts of the code. With respesct to the original names of the variables used 36! in the m_kpoint_grid module, a "ts_" was added. Also in the SIESTA module 37! kscell and kdispl are not public, but are made public (ts_kscell and 38! ts_kdispl) here. 39 40!============================================================================== 41! DETAILS: To obtain the kpoints for the GFs calculations it uses the same 42! scheme as for SIESTA but puts ts_kscell(3,3)=1 and ts_kdispl(3)=0.0 43! 44!============================================================================== 45 46 logical, public :: ts_scf_kgrid_first_time = .true. 47 logical, public :: ts_Gamma 48 integer, public :: ts_nkpnt ! Total number of k-points 49 real(dp), public :: ts_eff_kgrid_cutoff ! Effective kgrid_cutoff 50 51 real(dp), pointer, public :: ts_kweight(:) 52 real(dp), pointer, public :: ts_kpoint(:,:) 53 54 integer, public :: ts_kscell(3,3) = 0 55 real(dp), public :: ts_kdispl(3) = 0.0_dp 56 57 logical, public :: ts_firm_displ = .false. 58 logical, public :: ts_user_requested_mp = .false. 59 logical, public :: ts_user_requested_cutoff = .false. 60 61 public :: setup_ts_scf_kscell, setup_ts_kpoint_grid 62 public :: ts_write_k_points 63 64contains 65 66 subroutine setup_ts_scf_kscell( cell ) 67 68! ***************** INPUT ********************************************** 69! real*8 cell(3,3) : Unit cell vectors in real space cell(ixyz,ivec) 70! type(Elec) Elecs(:) : Electrodes 71 72! The relevant fdf labels are kgrid_cutoff and kgrid_Monkhorst_Pack. 73! If both are present, kgrid_Monkhorst_Pack has priority. If none is 74! present, the cutoff default is zero, producing only the gamma point. 75! Examples of fdf data specifications: 76! kgrid_cutoff 50. Bohr 77! %block kgrid_Monkhorst_Pack # Defines kscell and kdispl 78! 4 0 0 0.50 # (kscell(i,1),i=1,3), kdispl(1) 79! 0 4 0 0.50 # (kscell(i,2),i=1,3), kdispl(2) 80! 0 0 4 0.50 # (kscell(i,3),i=1,3), kdispl(3) 81! %endblock kgrid_Monkhorst_Pack 82! ********************************************************************** 83 84! Modules 85 86 use parallel, only : IONode 87 use m_minvec, only : minvec 88 use fdf 89 use sys, only : die 90#ifdef MPI 91 use mpi_siesta 92#endif 93 94 use m_ts_tdir, only: ts_tidx 95 use m_ts_global_vars, only : TSmode 96 97 implicit none 98 99 ! Passed variables 100 real(dp), intent(in) :: cell(3,3) 101 102 ! Internal variables 103 integer :: i, j, factor(3,3), expansion_factor 104 real(dp) :: scmin(3,3), vmod, cutoff 105 real(dp) :: ctransf(3,3) 106 107 type(block_fdf) :: bfdf 108 type(parsed_line), pointer :: pline 109 110 logical :: bool 111 real(dp), parameter :: defcut = 0.0_dp 112 113 bool = fdf_block('TS.kgrid_Monkhorst_Pack',bfdf) 114 if ( .not. bool ) & 115 bool = fdf_block('kgrid_Monkhorst_Pack',bfdf) 116 117 if ( bool ) then 118 ts_user_requested_mp = .true. 119 do i = 1,3 120 if (fdf_bline(bfdf,pline)) then 121 ts_kscell(1,i) = fdf_bintegers(pline,1) 122 ts_kscell(2,i) = fdf_bintegers(pline,2) 123 ts_kscell(3,i) = fdf_bintegers(pline,3) 124 if ( fdf_bnvalues(pline) > 3 ) then 125 ts_kdispl(i) = mod(fdf_bvalues(pline,4), 1._dp) 126 else 127 ts_kdispl(i) = 0._dp 128 end if 129 else 130 call die( 'setup_ts_scf_kscell: ERROR no data in' // & 131 'kgrid_Monkhorst_Pack block' ) 132 endif 133 enddo 134 call fdf_bclose(bfdf) 135 ts_firm_displ = .true. 136 137 else 138 139 if ( IONode .and. TSmode ) then 140 write(*,*) 'WARNING !!!' 141 write(*,*) 'TS kgrid determined first with 3D cell !!!' 142 write(*,*) 'Specifying only cutoff in Electrode AND Scattering calculations might lead to problems !!' 143 end if 144 145 cutoff = fdf_get('kgrid_cutoff',defcut,'Bohr') 146 ts_user_requested_cutoff = (cutoff /= defcut) 147 148 ts_kdispl(1:3) = 0.0_dp ! Might be changed later 149 ts_firm_displ = .false. 150 151 ! Find equivalent rounded unit-cell 152 call minvec( cell, scmin, ctransf ) 153 154 expansion_factor = 1 155 do j = 1,3 156 factor(j,1:3) = 0 157 vmod = sqrt(dot_product(scmin(1:3,j),scmin(1:3,j))) 158 factor(j,j) = int(2.0_dp*cutoff/vmod) + 1 159 expansion_factor = expansion_factor * factor(j,j) 160 enddo 161 ! Generate actual supercell skeleton 162 ts_kscell = matmul(ctransf, factor) 163 164 if ( expansion_factor == 1 ) then 165 ts_kscell = 0 166 ts_kscell(1,1) = 1 167 ts_kscell(2,2) = 1 168 ts_kscell(3,3) = 1 169 end if 170 171 end if 172 173 ! In case of TSmode we have a transiesta run. 174 ! This means that we truncate the k-points in the transport direction. 175 ! However, if we are dealing with an electrode calculation we simply allow it 176 ! to have the same cell (the check made will disregard the transport directions k-points) 177 if ( TSmode .and. ts_tidx > 0 ) then 178 i = ts_tidx 179 ts_kscell(:,i) = 0 180 ts_kscell(i,:) = 0 181 ts_kscell(i,i) = 1 182 ts_kdispl(i) = 0._dp 183 end if 184 185 end subroutine setup_ts_scf_kscell 186 187 subroutine setup_ts_kpoint_grid( ucell ) 188 189 ! SIESTA Modules 190 USE precision, only : dp 191 USE fdf, only : fdf_get 192 USE m_find_kgrid, only : find_kgrid 193 USE parallel, only : IONode 194 use m_ts_global_vars, only : TSmode 195 use kpoint_grid 196 197 ! Local Variables 198 real(dp), intent(in) :: ucell(3,3) 199 200 if ( .not. TSMode ) then 201 202 ! Same as SIESTA (mainly to be able to write the TSHS files). 203 ts_Gamma = gamma_SCF 204 ts_nkpnt = nkpnt 205 ts_eff_kgrid_cutoff = eff_kgrid_cutoff 206 ts_kweight => kweight 207 ts_kpoint => kpoint 208 ts_kscell = kscell 209 ts_kdispl = kdispl 210 211 ! Since this is not transiesta we immediately return 212 ! and then we also skip printing out transiesta information. 213 return 214 215 end if 216 217 if ( ts_scf_kgrid_first_time ) then 218 219 nullify(ts_kweight,ts_kpoint) 220 if ( fdf_get('SpinSpiral', .false.) ) then 221 call die('transiesta: Does not work with spin-spiral') 222 end if 223 224 call setup_ts_scf_kscell(ucell) 225 226 ts_scf_kgrid_first_time = .false. 227 228 else 229 if ( ts_user_requested_mp ) then 230 ! no need to set up the kscell again 231 else 232 call setup_ts_scf_kscell(ucell) 233 endif 234 endif 235 236 call find_kgrid(ucell,ts_kscell,ts_kdispl,ts_firm_displ, & 237 .true., & 238 ts_nkpnt,ts_kpoint,ts_kweight, ts_eff_kgrid_cutoff) 239 240 ts_Gamma = ts_nkpnt == 1 .and. & 241 dot_product(ts_kpoint(:,1),ts_kpoint(:,1)) < 1.0e-20_dp 242 243 if (IONode) call ts_write_k_points() 244 245 end subroutine setup_ts_kpoint_grid 246 247 subroutine ts_write_k_points() 248 USE siesta_options, only: writek 249 USE siesta_cml 250 251 implicit none 252 253 integer :: ik, ix, i 254 255 if ( writek ) then 256 write(*,'(/,a)') 'transiesta: ts_k-point coordinates (Bohr**-1) and weights:' 257 write(*,'(a,i4,3f12.6,3x,f12.6)') & 258 ('transiesta: ', ik, (ts_kpoint(ix,ik),ix=1,3), ts_kweight(ik), & 259 ik=1,ts_nkpnt) 260 endif 261 262 ! Always write the TranSIESTA k-points 263 call ts_iokp( ts_nkpnt, ts_kpoint, ts_kweight ) 264 265 write(*,'(/,a,i6)') 'transiesta: k-grid: Number of Green function k-points =', ts_nkpnt 266 write(*,'(a)') 'transiesta: k-grid: Supercell and displacements' 267 write(*,'(a,3i4,3x,f8.3)') 'transiesta: k-grid: ', & 268 (ts_kscell(i,1),i=1,3), ts_kdispl(1) 269 write(*,'(a,3i4,3x,f8.3)') 'transiesta: k-grid: ', & 270 (ts_kscell(i,2),i=1,3), ts_kdispl(2) 271 write(*,'(a,3i4,3x,f8.3)') 'transiesta: k-grid: ', & 272 (ts_kscell(i,3),i=1,3), ts_kdispl(3) 273 if (cml_p) then 274 call cmlStartPropertyList(xf=mainXML, title="Transiesta k-points", & 275 dictRef="siesta:ts_kpoints") 276 call cmlAddProperty(xf=mainXML, value=ts_nkpnt, dictref='siesta:ts_nkpnt', & 277 units="cmlUnits:countable") 278 do ik=1, ts_nkpnt 279 call cmlAddKPoint(xf=mainXML, coords=ts_kpoint(:,ik), weight=ts_kweight(ik)) 280 enddo 281 call cmlEndPropertyList(mainXML) 282 call cmlAddProperty(xf=mainXML, value=ts_kscell, dictref='siesta:ts_kscell', & 283 units="siestaUnits:Ang") 284 call cmlAddProperty(xf=mainXML, value=ts_kdispl, dictref='siesta:ts_kdispl', & 285 units="siestaUnits:Ang") 286 endif 287 end subroutine ts_write_k_points 288 289 subroutine ts_iokp( nk, points, weight ) 290! ******************************************************************* 291! Saves TranSIESTA k-points (only writing) Bohr^-1 292! Emilio Artacho, Feb. 1999 293! Modified by Nick Papior Andersen to not overwrite the SIESTA k-points 294! ********** INPUT ************************************************** 295! integer nk : Number of TS k-points 296! real*8 points(3,nk) : TS k-point coordinates 297! real*8 weight(3,nk) : TS k-point weight 298! ******************************************************************* 299 use fdf 300 use files, only : slabel, label_length 301 302 integer :: nk 303 real(dp) :: points(3,nk), weight(nk) 304 external :: io_assign, io_close 305 306! Internal 307 character(len=label_length+5) :: fname 308 integer :: iu, ik, ix 309! ------------------------------------------------------------------- 310 311 fname = trim( slabel ) // '.TSKP' 312 313 call io_assign( iu ) 314 open( iu, file=fname, form='formatted', status='unknown' ) 315 316 write(iu,'(i6)') nk 317 write(iu,'(i6,3f12.6,3x,f12.6)') & 318 (ik, (points(ix,ik),ix=1,3), weight(ik), ik=1,nk) 319 320 call io_close( iu ) 321 322 end subroutine ts_iokp 323 324end module m_ts_kpoints 325