TABLE OF CONTENTS
ABINIT/littlegroup_pert [ Functions ]
NAME
littlegroup_pert
FUNCTION
If syuse==0 and rfmeth==2, determines the set of symmetries that leaves a perturbation invariant. (Actually, all symmetries that leaves a q-wavevector invariant should be used to reduce the number of k-points for all perturbations. Unfortunately, one has to take into account the sign reversal of the perturbation under the symmetry operations, which makes GS routines not usable for the respfn code. The intermediate choice was to select only those that keep also the perturbation invariant. Note that the wavevector of the perturbation must also be invariant, a translation vector in real space is NOT allowed ).
COPYRIGHT
Copyright (C) 1999-2018 ABINIT group (XG, DRH) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
gprimd(3,3)=dimensional primitive translations for reciprocal space (bohr**-1) idir=direction of the perturbation indsym(4,nsym,natom)=indirect indexing of atom labels--see subroutine symatm for definition (if nsym>1) iout=if non-zero, output on unit iout ipert=characteristics of the perturbation natom= number of atoms nsym=number of space group symmetries rfmeth = 1 if non-stationary block 2 if stationary block 3 if third order derivatives symq(4,2,nsym)= Table computed by littlegroup_q. three first numbers define the G vector; fourth number is zero if the q-vector is not preserved, is 1 otherwise second index is one without time-reversal symmetry, two with time-reversal symmetry symafm(nsym)=(anti)ferromagnetic part of the symmetry operations symrec(3,3,nsym)=3x3 matrices of the group symmetries (reciprocal space) symrel(3,3,nsym)=3x3 matrices of the group symmetries (real space) syuse= flag to use the symmetries or not. If 0 usei it, if 1 do not use it. tnons(3,nsym)=nonsymmorphic translations of space group in terms of real space primitive translations (may be 0) [unit]=By default the routine writes to std_out and this is very annoying if we are inside a big loop. Use unit=dev_null or a negative integer to disable writing.
OUTPUT
nsym1 =number of space group symmetries that leaves the perturbation invariant symaf1(nsym1)=(anti)ferromagnetic part of the corresponding symmetry operations symrl1(3,3,nsym1)=corresponding 3x3 matrices of the group symmetries (real space) tnons1(3,nsym1)=corresponding nonsymmorphic translations of space group in terms of real space primitive translations (may be 0)!!
PARENTS
dfpt_looppert,get_npert_rbz,m_dvdb,read_gkk
CHILDREN
stresssym,wrtout
SOURCE
63 #if defined HAVE_CONFIG_H 64 #include "config.h" 65 #endif 66 67 #include "abi_common.h" 68 69 subroutine littlegroup_pert(gprimd,idir,indsym,iout,ipert,natom,nsym,nsym1, & 70 & rfmeth,symafm,symaf1,symq,symrec,symrel,symrl1,syuse,tnons,tnons1, & 71 & unit) ! Optional 72 73 use defs_basis 74 use m_profiling_abi 75 use m_errors 76 77 !This section has been created automatically by the script Abilint (TD). 78 !Do not modify the following lines by hand. 79 #undef ABI_FUNC 80 #define ABI_FUNC 'littlegroup_pert' 81 use interfaces_14_hidewrite 82 use interfaces_41_geometry, except_this_one => littlegroup_pert 83 !End of the abilint section 84 85 implicit none 86 87 !Arguments ------------------------------- 88 !scalars 89 integer,intent(in) :: idir,iout,ipert,natom,nsym,rfmeth,syuse 90 integer,intent(in),optional :: unit 91 integer,intent(out) :: nsym1 92 !arrays 93 integer,intent(in) :: indsym(4,nsym,natom),symafm(nsym),symq(4,2,nsym) 94 integer,intent(in) :: symrec(3,3,nsym),symrel(3,3,nsym) 95 integer,intent(out) :: symaf1(nsym),symrl1(3,3,nsym) 96 real(dp),intent(in) :: gprimd(3,3),tnons(3,nsym) 97 real(dp),intent(out) :: tnons1(3,nsym) 98 99 !Local variables ------------------------- 100 !scalars 101 integer :: idir1,ii,istr,isym,jj,nsym_test,tok,ount 102 character(len=500) :: msg 103 !arrays 104 integer :: sym_test(3,3,2) 105 real(dp) :: str_test(6) 106 107 ! ********************************************************************* 108 109 ount = std_out; if (present(unit)) ount = unit 110 111 nsym1=0 112 if((ipert==natom+3 .or. ipert==natom+4) .and. syuse==0 .and. rfmeth==2) then 113 ! Strain perturbation section 114 ! Use ground state routine which symmetrizes cartesian stress as a quick 115 ! and dirty test for the invariance of the strain (ipert,idir) under 116 ! each candidate symmetry 117 ! I am presently assuming that translations are acceptable because I dont 118 ! see why not. 119 120 istr=3*(ipert-natom-3)+idir 121 nsym_test=2 122 ! Store identity as first element for test 123 sym_test(:,:,1)=0 124 sym_test(1,1,1)=1; sym_test(2,2,1)=1; sym_test(3,3,1)=1 125 do isym=1,nsym 126 sym_test(:,:,2)=symrec(:,:,isym) 127 str_test(:)=0.0_dp 128 str_test(istr)=1.0_dp 129 call stresssym(gprimd,nsym_test,str_test,sym_test) 130 if(abs(str_test(istr)-1.0_dp)<tol8)then 131 ! The test has been successful ! 132 nsym1=nsym1+1 133 symaf1(nsym1)=symafm(isym) 134 do ii=1,3 135 tnons1(ii,nsym1)=tnons(ii,isym) 136 do jj=1,3 137 symrl1(ii,jj,nsym1)=symrel(ii,jj,isym) 138 end do 139 end do 140 end if 141 end do 142 143 else if(ipert>natom .or. syuse/=0 .or. rfmeth/=2)then 144 145 ! Not yet coded for d/dk or electric field perturbations 146 nsym1=1 147 do ii=1,3 148 tnons1(ii,1)=0._dp 149 symaf1(1)=1 150 do jj=1,3 151 symrl1(ii,jj,1)=0 152 if(ii==jj)symrl1(ii,jj,1)=1 153 end do 154 end do 155 156 else 157 158 do isym=1,nsym 159 ! Check that the symmetry operation preserves the wavevector 160 ! (a translation is NOT allowed) 161 if(symq(4,1,isym)==1 .and.& 162 & symq(1,1,isym)==0 .and.& 163 & symq(2,1,isym)==0 .and.& 164 & symq(3,1,isym)==0 )then 165 ! Check that the symmetry operation preserves the atom 166 if(ipert==indsym(4,isym,ipert))then 167 ! Check if the direction is preserved 168 tok=1 169 do idir1=1,3 170 if((idir1==idir.and.symrec(idir,idir1,isym)/=1) .or.& 171 & (idir1/=idir.and.symrec(idir,idir1,isym)/=0))then 172 tok=0 173 end if 174 end do 175 if(tok==1)then 176 ! All the tests have been successful ! 177 nsym1=nsym1+1 178 symaf1(nsym1)=symafm(isym) 179 do ii=1,3 180 tnons1(ii,nsym1)=tnons(ii,isym) 181 do jj=1,3 182 symrl1(ii,jj,nsym1)=symrel(ii,jj,isym) 183 end do 184 end do 185 end if 186 187 end if 188 end if 189 end do 190 end if 191 192 if (nsym1<1) then 193 write(msg,'(a,i0,a)')& 194 & ' The number of selected symmetries should be > 0, while it is nsym=',nsym1,'.' 195 MSG_BUG(msg) 196 end if 197 198 if (nsym1 /= 1) then 199 if (iout /= ount .and. iout > 0) then 200 write(msg,'(a,i5,a)')' Found ',nsym1,' symmetries that leave the perturbation invariant.' 201 call wrtout(iout,msg,'COLL') 202 end if 203 write(msg,'(a,i5,a)')' littlegroup_pert : found ',nsym1,' symmetries that leave the perturbation invariant :' 204 call wrtout(ount,msg,'COLL') 205 else 206 if (iout /= ount .and. iout > 0) then 207 write(msg,'(a,a)')' The set of symmetries contains',' only one element for this perturbation.' 208 call wrtout(iout,msg,'COLL') 209 end if 210 write(msg,'(a)')' littlegroup_pert : only one element in the set of symmetries for this perturbation :' 211 call wrtout(ount,msg,'COLL') 212 end if 213 214 if (ount > 0) then 215 do isym=1,nsym1 216 write(msg, '(9i4)' )((symrl1(ii,jj,isym),ii=1,3),jj=1,3) 217 call wrtout(ount,msg,'COLL') 218 end do 219 end if 220 221 end subroutine littlegroup_pert