TABLE OF CONTENTS
ABINIT/mk_irredpert [ Functions ]
NAME
mk_irredpert
FUNCTION
This routine finds the symop combinations needed to calculate other perturbations based on the present one.
COPYRIGHT
Copyright (C) 2004-2018 ABINIT group (MVer) This file is distributed under the terms of the GNU General Public Licence, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
indsym = indirect indexing array for symmetries iqptfull = qpoint index in full grid natom = number of atoms nbranch = 3*natom = number of phonon branches nqpt = total number of qpoints nsym = number of symops qpt = qpoint in reduced coordinates qtimrev = 1 or 0, use time reversal symmetry (for Gamma only) or do not symq = 1 if symmetry preserves present qpoint. From littlegroup_q symrel = symmetry operations in recip space (or real)
OUTPUT
irredpert(7,nbranch,nqpt) = indices for reconstructing perturbations explained in NOTES. If irredpert < 0 then there is no information for reconstructing the element, and it should be read in (it is irreducible).
NOTES
irredpert (1,,,) = iatsy1 irredpert (2,,,) = iatsy2 irredpert (3,,,) = isym irredpert (4,,,) = itim irredpert (5,,,) = factor for idir2 = 1 irredpert (6,,,) = factor for idir2 = 2 irredpert (7,,,) = factor for idir2 = 3
PARENTS
CHILDREN
SOURCE
51 #if defined HAVE_CONFIG_H 52 #include "config.h" 53 #endif 54 55 #include "abi_common.h" 56 57 58 subroutine mk_irredpert(indsym,iqptfull,irredpert,& 59 & natom,nbranch,nqpt,nsym,qpt,qtimrev,symq,symrel) 60 61 use defs_basis 62 use m_profiling_abi 63 64 !This section has been created automatically by the script Abilint (TD). 65 !Do not modify the following lines by hand. 66 #undef ABI_FUNC 67 #define ABI_FUNC 'mk_irredpert' 68 !End of the abilint section 69 70 implicit none 71 72 !Arguments ------------------------------------ 73 !scalars 74 integer,intent(in) :: iqptfull,natom,nbranch,nqpt,nsym,qtimrev 75 !arrays 76 integer,intent(in) :: indsym(4,nsym,natom),qpt(3),symq(4,2,nsym) 77 integer,intent(in) :: symrel(3,3,nsym) 78 integer,intent(out) :: irredpert(7,nbranch,nbranch,nqpt) 79 80 !Local variables ------------------------- 81 !scalars 82 integer :: found,iatom1,iatom2,iatsy1,iatsy2,idir1,idir2,idisy1,idisy2 83 integer :: ipert1,ipert2,ipertsy1,ipertsy2,isign,isym,ithree,itim 84 integer :: noccur,quit,quit1 85 real(dp) :: arg1,arg2,im,re 86 !arrays 87 integer :: indsyminv(4,nsym,natom),sym1(3,3),sym2(3,3) 88 89 ! ************************************************************************* 90 91 !------------------------------------------------------- 92 !for present symmetries, find equivalent perturbations 93 !------------------------------------------------------- 94 !irredpert (1,,,) = ipert1 95 !irredpert (2,,,) = ipert2 96 !irredpert (3,,,) = isym 97 !irredpert (4,,,) = itim 98 !irredpert (5,,,) = factor for idir2 = 1 99 !irredpert (6,,,) = factor for idir2 = 2 100 !irredpert (7,,,) = factor for idir2 = 3 101 ! 102 103 !Invert the indsym array 104 do iatom1=1,natom 105 do isym=1,nsym 106 indsyminv(1,nsym,indsym(4,isym,iatom1)) = -indsym(1,isym,iatom1) 107 indsyminv(2,nsym,indsym(4,isym,iatom1)) = -indsym(2,isym,iatom1) 108 indsyminv(3,nsym,indsym(4,isym,iatom1)) = -indsym(3,isym,iatom1) 109 indsyminv(4,isym,indsym(4,isym,iatom1)) = iatom1 110 end do 111 end do 112 113 !********************************************************************* 114 115 !DEBUG 116 !write(std_out,*)' mk_irredpert : enter mk_irredpert ' 117 !write(std_out,*)' mk_irredpert : qtimrev=',qtimrev 118 !ENDDEBUG 119 120 121 !Big Big Loop : symmetrize three times, because 122 !of some cases in which one element is not yet available 123 !at the first pass, and even at the second one ! 124 do ithree=1,3 125 126 ! Big loop on all elements 127 do iatom1=1,natom+2 128 do idir1=1,3 129 ipert1=(iatom1-1)*3+idir1 130 do iatom2=1,natom+2 131 do idir2=1,3 132 ipert2=(iatom2-1)*3+idir2 133 134 ! Will try to eliminate element (idir1,iatom1,idir2,iatom2) 135 ! so this element should not have been eliminated yet... 136 if(irredpert (1,ipert1,ipert2,iqptfull) > 0) cycle 137 ! write(std_out,*) 'trying to eliminate ',idir1,iatom1,idir2,iatom2 138 139 ! Loop on all symmetries, including time-reversal 140 quit1=0 141 do isym=1,nsym 142 do itim=1,2 143 isign=3-2*itim 144 145 if(symq(4,itim,isym)==0) cycle 146 found=1 147 148 ! Here select the symmetric of iatom1 149 if(iatom1<=natom)then 150 iatsy1=indsyminv(4,isym,iatom1) 151 sym1(:,:)=symrel(:,:,isym) 152 else 153 found=0 154 end if 155 156 ! Here select the symmetric of iatom2 157 if(iatom2<=natom)then 158 iatsy2=indsyminv(4,isym,iatom2) 159 ! try this: invert all relations symrec -> symrel 160 sym2(:,:)=symrel(:,:,isym) 161 else 162 found=0 163 end if 164 165 ! Now that a symmetric perturbation has been obtained, 166 ! including the expression of the symmetry matrix, see 167 ! if the symmetric values are available 168 ! first if found==1 169 if( found==1 ) then 170 171 noccur=0 172 quit=0 173 do idisy1=1,3 174 ipertsy1 = (iatsy1-1)*3+idisy1 175 do idisy2=1,3 176 ipertsy2 = (iatsy2-1)*3+idisy2 177 ! 178 ! NOTE : May be transpose on the next line 179 ! 180 if(sym1(idir1,idisy1)/=0 .and. sym2(idir2,idisy2)/=0 )then 181 if(irredpert (1,ipertsy1,ipertsy2,iqptfull) < 0)then 182 ! nothing, ipertsy1,ipertsy2 has not been eliminated by symmetry, 183 ! just leave found == 1 184 185 186 else 187 ! Not found: ipert1, ipert2 can not be reconstituted from 188 ! ipertsy1, ipertsy2 using isym 189 ! write(std_out,*) 'Not found: ',ipert1,ipert2,& 190 ! & ' can not be reconstituted from ', ipertsy1,ipertsy2,& 191 ! & ' using ', isym 192 found=0 193 quit=1 194 exit 195 end if 196 197 ! Here, in case the symmetric of the element 198 ! is the element, or the symmetric with 199 ! respect to permutation of perturbations 200 ! (some more conditions on the time-reversal 201 ! symmetry must be fulfilled although) 202 if( idisy1==idir1 .and. iatsy1==iatom1& 203 & .and. idisy2==idir2 .and. iatsy2==iatom2& 204 & .and.(isign==1 .or. qtimrev==1 & 205 & .or. (idir1==idir2 .and. iatom1==iatom2)))& 206 & then 207 noccur=noccur+sym1(idir1,idisy1)*sym2(idir2,idisy2) 208 else if( idisy1==idir2 .and. iatsy1==iatom2& 209 & .and. idisy2==idir1 .and. iatsy2==iatom1& 210 & .and.(isign==-1 .or. qtimrev==1& 211 & .or. (idir1==idir2 .and. iatom1==iatom2)))& 212 & then 213 noccur=noccur+sym1(idir1,idisy1)*sym2(idir2,idisy2) 214 end if 215 216 end if 217 end do 218 if(quit==1)exit 219 end do 220 end if 221 ! End first if found==1 222 223 ! SHOULD THIS BE CHECKED IN THE ELPHON CASE TOO? 224 ! second if found == 1 225 if(found==1)then 226 ! In case of phonons, need to take into account the 227 ! time-reversal symmetry, and the shift back to the unit cell 228 ! 229 if(ipert1<=natom .and. ipert2<=natom)then 230 ! 2) Shift the atoms back to the unit cell. 231 arg1=two_pi*( qpt(1)*indsyminv(1,isym,iatom1)& 232 & +qpt(2)*indsyminv(2,isym,iatom1)& 233 & +qpt(3)*indsyminv(3,isym,iatom1) ) 234 arg2=two_pi*( qpt(1)*indsyminv(1,isym,iatom2)& 235 & +qpt(2)*indsyminv(2,isym,iatom2)& 236 & +qpt(3)*indsyminv(3,isym,iatom2) ) 237 re=cos(arg1)*cos(arg2)+sin(arg1)*sin(arg2) 238 ! XG010117 Must use isign 239 im=isign*(cos(arg2)*sin(arg1)-cos(arg1)*sin(arg2)) 240 else 241 re=1.0_dp 242 im=0.0_dp 243 end if 244 245 ! Final check, could still fail if the 246 ! element was its own symmetric 247 if( abs(1.0_dp-re*noccur) < 1.0d-6& 248 & .and. abs(im*noccur) < 1.0d-6 )then 249 found=0 250 251 ! DEBUG 252 ! write(std_out,*)' element is its own symmetric ...' 253 ! ENDDEBUG 254 255 end if 256 ! End if element was its own symmetric 257 258 end if 259 ! End second if found == 1 260 261 ! third if found == 1 262 if(found==1)then 263 264 ! DEBUG 265 ! write(std_out,*)' all found ! isym, isign= ',isym,isign 266 ! write(std_out,'(9i4)' )((sym1(ii,jj),ii=1,3),jj=1,3) 267 ! write(std_out,'(9i4)' )((sym2(ii,jj),ii=1,3),jj=1,3) 268 ! ENDDEBUG 269 270 ! The element has been constructed ! 271 ! ipertsy1,ipertsy2 do not correspond to anything: all 272 ! pertsy's were scanned above, but all the needed ones 273 ! are present. 274 irredpert (1,ipert1,ipert2,iqptfull) = iatsy1 275 irredpert (2,ipert1,ipert2,iqptfull) = iatsy2 276 irredpert (3,ipert1,ipert2,iqptfull) = isym 277 irredpert (4,ipert1,ipert2,iqptfull) = itim 278 irredpert (5,ipert1,ipert2,iqptfull) = 0 279 irredpert (6,ipert1,ipert2,iqptfull) = 0 280 irredpert (7,ipert1,ipert2,iqptfull) = 0 281 282 ! Exit loop on symmetry operations 283 quit1=1 284 exit 285 286 end if 287 ! end third if found == 1 288 289 ! End loop on all symmetries + time-reversal 290 end do 291 ! end isym do 292 if(quit1==1)exit 293 end do 294 ! end itim do 295 296 ! End big loop on all elements 297 298 end do 299 end do 300 ! end ipert2 do 301 end do 302 end do 303 ! end ipert1 do 304 305 ! End Big Big Loop 306 end do 307 308 !DEBUG 309 !write(std_out,*) ' mk_irredpert : irredpert for qpt ',iqptfull,' = ' 310 !do ibranch=1,nbranch 311 !do jbranch=1,nbranch 312 !write(std_out,'(2i4,2x,7i6)') ibranch,jbranch,irredpert(:,ibranch,jbranch,iqptfull) 313 !end do 314 !end do 315 !ENDDEBUG 316 317 end subroutine mk_irredpert