TABLE OF CONTENTS


ABINIT/littlegroup_q [ Functions ]

[ Top ] [ Functions ]

NAME

 littlegroup_q

FUNCTION

 Determines the symmetry operations by which reciprocal vector q is preserved,
 modulo a primitive reciprocal lattice vector, and the time-reversal symmetry.

COPYRIGHT

 Copyright (C) 1999-2018 ABINIT group (GMR)
 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

 nsym=number of space group symmetries
 qpt(3)= vector in reciprocal space
 symrec(3,3,nsym)=3x3 matrices of the group symmetries (reciprocal space)
 [prtvol]=integer flag defining the verbosity of output. =0 if no output is provided.
 prtgkk= integer flag. If 1 provide output of electron-phonon "gkk" matrix elements, for further
     treatment by mrggkk utility or anaddb utility. If 0 no output is provided.

OUTPUT

 symq(4,2,nsym)= 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
 timrev=1 if the time-reversal symmetry preserves the wavevector, modulo a reciprocal lattice vector (in principle, see below).

NOTES

 The condition is:
    $q =  O  S(q) - G$
 with O being either the identity or the time reversal symmetry (= inversion in reciprocal space)
 and G being a primitive vector of the reciprocal lattice.
 If the time-reversal (alone) also preserves q, modulo a lattice vector, then timrev is set to 1, otherwise 0.

TODO

 timrev is put to 1 only for Gamma.
 Better handling should be provided in further version.

PARENTS

      get_npert_rbz,m_bz_mesh,m_ddb,m_dvdb,m_dynmat,m_esymm,m_gkk,m_phgamma
      m_sigmaph,memory_eval,read_gkk,respfn

CHILDREN

      wrap2_pmhalf,wrtout

SOURCE

 51 #if defined HAVE_CONFIG_H
 52 #include "config.h"
 53 #endif
 54 
 55 #include "abi_common.h"
 56 
 57 
 58 subroutine littlegroup_q(nsym,qpt,symq,symrec,symafm,timrev,prtvol,use_sym)
 59 
 60  use defs_basis
 61 
 62  use m_numeric_tools,   only : wrap2_pmhalf
 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 'littlegroup_q'
 68  use interfaces_14_hidewrite
 69 !End of the abilint section
 70 
 71  implicit none
 72 
 73 !Arguments -------------------------------
 74 !scalars
 75  integer,intent(in) :: nsym
 76  integer,intent(in),optional :: prtvol,use_sym
 77  integer,intent(out) :: timrev
 78 !arrays
 79  integer,intent(in) :: symrec(3,3,nsym)
 80  integer,intent(in) :: symafm(nsym)
 81  integer,intent(out) :: symq(4,2,nsym)
 82  real(dp),intent(in) :: qpt(3)
 83 
 84 !Local variables -------------------------
 85 !scalars
 86  integer :: ii,isign,isym,itirev,my_prtvol
 87  real(dp),parameter :: tol=2.d-8
 88  real(dp) :: reduce
 89  character(len=500) :: message
 90 !arrays
 91  real(dp) :: difq(3),qsym(3),shift(3)
 92 
 93 ! *********************************************************************
 94 
 95  my_prtvol=0 ; if (PRESENT(prtvol)) my_prtvol=prtvol
 96 
 97 ! Initialise the array symq
 98  symq = 0
 99 
100  isym = symafm(1) ! just to fool abirules and use symafm for the moment
101 
102  do isym=1,nsym
103 !   if (symafm(isym) /= 1) cycle ! skip afm symops
104 ! TODO: check how much of the afm syms are coded in the rf part of the code. cf
105 ! test v3 / 12
106    do itirev=1,2
107      isign=3-2*itirev  ! isign is 1 without time-reversal, -1 with time-reversal
108 
109 !    Get the symmetric of the vector
110      do ii=1,3
111        qsym(ii)=qpt(1)*isign*symrec(ii,1,isym)&
112 &       +qpt(2)*isign*symrec(ii,2,isym)&
113 &       +qpt(3)*isign*symrec(ii,3,isym)
114      end do
115 
116 !    Get the difference between the symmetric and the original vector
117 
118      symq(4,itirev,isym)=1
119      do ii=1,3
120        difq(ii)=qsym(ii)-qpt(ii)
121 !      Project modulo 1 in the interval ]-1/2,1/2] such that difq = reduce + shift
122        call wrap2_pmhalf(difq(ii),reduce,shift(ii))
123        if(abs(reduce)>tol)symq(4,itirev,isym)=0
124      end do
125 
126 !    SP: When prtgkk is asked (GKK matrix element will be output), one has to
127 !    disable symmetries. There is otherwise a jauge problem with the unperturbed
128 !    and the perturbed wavefunctions. This leads to a +- 5% increase in computational
129 !    cost but provide the correct GKKs (i.e. the same as without the use of
130 !    symmerties.)
131 
132      if (PRESENT(use_sym)) then
133        if (use_sym == 0) then
134          symq(4,itirev,isym)=0
135          symq(4,itirev,1)=1
136        end if
137      end if
138 
139 !    If the operation succeded, change shift from real(dp) to integer, then exit loop
140      if(symq(4,itirev,isym)/=0)then
141        if (my_prtvol>0) then
142          if(itirev==1)write(message,'(a,i4,a)')' littlegroup_q : found symmetry',isym,' preserves q '
143          if(itirev==2)write(message,'(a,i4,a)')' littlegroup_q : found symmetry ',isym,' + TimeReversal preserves q '
144          call wrtout(std_out,message,'COLL')
145        end if
146 !      Uses the mathematical function NINT = nearest integer
147        do ii=1,3
148          symq(ii,itirev,isym)=nint(shift(ii))
149        end do
150      end if
151 
152    end do !itirev
153  end do !isym
154 
155 !Test time-reversal symmetry
156  timrev=1
157  do ii=1,3
158 !  Unfortunately, this version does not work yet ...
159 !  call wrap2_pmhalf(2*qpt(ii),reduce,shift(ii))
160 !  if(abs(reduce)>tol)timrev=0
161 !  So, this is left ...
162    if(abs(qpt(ii))>tol)timrev=0
163  end do
164 
165  if(timrev==1.and.my_prtvol>0)then
166    write(message, '(a,a,a)' )&
167 &   ' littlegroup_q : able to use time-reversal symmetry. ',ch10,&
168 &   '  (except for gamma, not yet able to use time-reversal symmetry)'
169    call wrtout(std_out,message,'COLL')
170  end if
171 
172 end subroutine littlegroup_q