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 Kpoint_grid 9! 10! Contains data structures and routines to deal with the kpoint-grid 11! for the self-consistent calculation 12! Other uses (bands, optical, polarization) have their own structures. 13! 14 USE precision, only : dp 15 16 implicit none 17 18 public :: setup_kpoint_grid, scf_kgrid_first_time, gamma_scf, & 19 nkpnt, kweight, kpoint, kscell, kdispl 20 public :: eff_kgrid_cutoff 21 22 private 23 24 logical :: scf_kgrid_first_time = .true. 25 logical :: gamma_scf 26 integer :: nkpnt ! Total number of k-points 27 real(dp) :: eff_kgrid_cutoff ! Effective kgrid_cutoff 28 29 real(dp), pointer :: kweight(:) => null() 30 real(dp), pointer :: kpoint(:,:) => null() 31 32 integer, dimension(3,3) :: kscell = 0 33 real(dp), dimension(3) :: kdispl = 0.0_dp 34 logical :: user_requested_mp = .false. 35 logical :: user_requested_cutoff = .false. 36 37 logical, save :: time_reversal_symmetry = .true. 38 logical :: firm_displ = .false. 39 40 CONTAINS 41 42 subroutine setup_Kpoint_grid( ucell ) 43 USE parallel, only : Node 44 USE fdf, only : fdf_defined, fdf_get 45 USE m_find_kgrid, only : find_kgrid 46 use m_spin, only: TrSym 47 48 implicit none 49 real(dp) :: ucell(3,3) 50 logical :: spiral 51 52 if (scf_kgrid_first_time) then 53 spiral = fdf_defined('SpinSpiral') 54 ! Allow the user to control the use of time-reversal-symmetry 55 ! By default, it is on, except for "spin-spiral" calculations 56 ! and/or non-collinear calculations 57 time_reversal_symmetry = fdf_get( & 58 "TimeReversalSymmetryForKpoints", & 59 (.not. spiral) .and. TrSym) 60 call setup_scf_kscell(ucell, firm_displ) 61 62 scf_kgrid_first_time = .false. 63 64 else 65 if ( user_requested_mp ) then 66 ! no need to set up the kscell again 67 else 68 ! This was wrong in the old code 69 call setup_scf_kscell(ucell, firm_displ) 70 endif 71 endif 72 73 call find_kgrid(ucell,kscell,kdispl,firm_displ, & 74 time_reversal_symmetry, & 75 nkpnt,kpoint,kweight, eff_kgrid_cutoff) 76 77 gamma_scf = (nkpnt == 1 .and. & 78 dot_product(kpoint(:,1),kpoint(:,1)) < 1.0e-20_dp) 79 80 if (Node .eq. 0) call siesta_write_k_points() 81 82 end subroutine setup_Kpoint_grid 83 84!-------------------------------------------------------------------- 85 subroutine setup_scf_kscell( cell, firm_displ ) 86 87! ***************** INPUT ********************************************** 88! real*8 cell(3,3) : Unit cell vectors in real space cell(ixyz,ivec) 89! ***************** OUTPUT ********************************************* 90! logical firm_displ : User-specified displacements (firm)? 91 92! The relevant fdf labels are kgrid_cutoff and kgrid_Monkhorst_Pack. 93! If both are present, kgrid_Monkhorst_Pack has priority. If none is 94! present, the cutoff default is zero, producing only the gamma point. 95! Examples of fdf data specifications: 96! kgrid_cutoff 50. Bohr 97! %block kgrid_Monkhorst_Pack # Defines kscell and kdispl 98! 4 0 0 0.50 # (kscell(i,1),i=1,3), kdispl(1) 99! 0 4 0 0.50 # (kscell(i,2),i=1,3), kdispl(2) 100! 0 0 4 0.50 # (kscell(i,3),i=1,3), kdispl(3) 101! %endblock kgrid_Monkhorst_Pack 102! ********************************************************************** 103 104! Modules 105 106 use precision, only : dp 107 use parallel, only : Node 108 use m_minvec, only : minvec 109 use sys, only : die 110 use fdf 111 112 implicit none 113 114! Passed variables 115 real(dp), intent(in) :: cell(3,3) 116 logical, intent(out) :: firm_displ 117 118! Internal variables 119 integer i, j, factor(3,3), expansion_factor 120 real(dp) scmin(3,3), vmod, cutoff 121 real(dp) ctransf(3,3) 122 logical mp_input 123 124 real(dp), parameter :: defcut = 0.0_dp 125 integer, dimension(3,3), parameter :: unit_matrix = & 126 reshape ((/1,0,0,0,1,0,0,0,1/), (/3,3/)) 127 128 type(block_fdf) :: bfdf 129 type(parsed_line), pointer :: pline 130 131 132 mp_input = fdf_block('kgrid_Monkhorst_Pack',bfdf) 133 if ( mp_input ) then 134 user_requested_mp = .true. 135 do i= 1, 3 136 if (.not. fdf_bline(bfdf,pline)) & 137 call die('setup_scf_kscell: ERROR in ' // & 138 'kgrid_Monkhorst_Pack block') 139 kscell(1,i) = fdf_bintegers(pline,1) 140 kscell(2,i) = fdf_bintegers(pline,2) 141 kscell(3,i) = fdf_bintegers(pline,3) 142 if ( fdf_bnvalues(pline) > 3 ) then 143 kdispl(i) = mod(fdf_bvalues(pline,4), 1._dp) 144 else 145 kdispl(i) = 0._dp 146 end if 147 enddo 148 call fdf_bclose(bfdf) 149 firm_displ = .true. 150 151 else 152 153 cutoff = fdf_physical('kgrid_cutoff',defcut,'Bohr') 154 if (cutoff /= defcut) then 155 !! write(6,"(a,f10.5)") "Kgrid cutoff input: ", cutoff 156 user_requested_cutoff = .true. 157 endif 158 159 kdispl(1:3) = 0.0_dp ! Might be changed later 160 firm_displ = .false. ! In future we might add new options 161 ! for user-specified displacements 162 163 ! Find equivalent rounded unit-cell 164 call minvec( cell, scmin, ctransf ) 165 166 expansion_factor = 1 167 do j = 1,3 168 factor(j,1:3) = 0 169 vmod = sqrt(dot_product(scmin(1:3,j),scmin(1:3,j))) 170 factor(j,j) = int(2.0_dp*cutoff/vmod) + 1 171 expansion_factor = expansion_factor * factor(j,j) 172 enddo 173 ! Generate actual supercell skeleton 174 kscell = matmul(ctransf, factor) 175 ! Avoid confusing permutations 176 if (expansion_factor == 1) then 177 kscell = unit_matrix 178 endif 179 endif 180 181 end subroutine setup_scf_kscell 182 183 subroutine siesta_write_k_points() 184 USE siesta_options, only: writek 185 USE units, only: Ang 186 USE siesta_cml 187 188 implicit none 189 190 integer :: ik, ix, i 191 external :: iokp 192 193 if ( writek ) then 194 write(6,'(/,a)') 'siesta: k-point coordinates (Bohr**-1) and weights:' 195 write(6,'(a,i4,3f12.6,3x,f12.6)') & 196 ('siesta: ', ik, (kpoint(ix,ik),ix=1,3), kweight(ik), & 197 ik=1,nkpnt) 198 endif 199 ! Always write KP file 200 call iokp( nkpnt, kpoint, kweight ) 201 write(6,'(/a,i6)') 'siesta: k-grid: Number of k-points =', nkpnt 202 write(6,'(a,f10.3,a)') 'siesta: k-grid: Cutoff (effective) =', & 203 eff_kgrid_cutoff/Ang, ' Ang' 204 write(6,'(a)') 'siesta: k-grid: Supercell and displacements' 205 write(6,'(a,3i4,3x,f8.3)') 'siesta: k-grid: ', & 206 (kscell(i,1),i=1,3), kdispl(1) 207 write(6,'(a,3i4,3x,f8.3)') 'siesta: k-grid: ', & 208 (kscell(i,2),i=1,3), kdispl(2) 209 write(6,'(a,3i4,3x,f8.3)') 'siesta: k-grid: ', & 210 (kscell(i,3),i=1,3), kdispl(3) 211 if (cml_p) then 212 call cmlStartPropertyList(xf=mainXML, title="k-points", & 213 dictRef="siesta:kpoints") 214 call cmlAddProperty(xf=mainXML, value=nkpnt, dictref='siesta:nkpnt', & 215 units="cmlUnits:countable") 216 do ik=1, nkpnt 217 call cmlAddKPoint(xf=mainXML, coords=kpoint(:,ik), weight=kweight(ik)) 218 enddo 219 call cmlAddProperty(xf=mainXML, value=eff_kgrid_cutoff/Ang, & 220 dictref='siesta:kcutof', units='siestaUnits:angstrom') 221 call cmlEndPropertyList(mainXML) 222 call cmlAddProperty(xf=mainXML, value=kscell, dictref='siesta:kscell', & 223 units="siestaUnits:Ang") 224 call cmlAddProperty(xf=mainXML, value=kdispl, dictref='siesta:kdispl', & 225 units="siestaUnits:Ang") 226 endif 227 end subroutine siesta_write_k_points 228 229END MODULE Kpoint_grid 230