TABLE OF CONTENTS


ABINIT/dfpt_symph [ Functions ]

[ Top ] [ Functions ]

NAME

 dfpt_symph

FUNCTION

 Determine the symmetry character of the different phonon modes.

COPYRIGHT

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

 acell(3)=length scales of primitive translations (bohr)
 eigvec(2*3*natom*3*natom)=eigenvectors of the dynamical matrix
 indsym(4,nsym,natom)=indirect indexing array : for each
   isym,iatom, fourth element is label of atom into which iatom is sent by
   INVERSE of symmetry operation isym; first three elements are the primitive
   translations which must be subtracted after the transformation to get back
   to the original unit cell.
 iout=unit number to which output is written
 natom=number of atoms in unit cell
 nsym=number of space group symmetries
 phfrq(3*natom)=phonon frequencies (Hartree units)
 rprim(3,3)=dimensionless primitive translations in real space
 symrel(3,3,nsym)=matrices of the group symmetries (real space)

OUTPUT

PARENTS

      anaddb,m_phonons

CHILDREN

      matr3inv,mkrdim,wrtout

SOURCE

 41 #if defined HAVE_CONFIG_H
 42 #include "config.h"
 43 #endif
 44 
 45 #include "abi_common.h"
 46 
 47 
 48 subroutine dfpt_symph(iout,acell,eigvec,indsym,natom,nsym,phfrq,rprim,symrel)
 49 
 50  use defs_basis
 51  use m_profiling_abi
 52 
 53  use m_geometry,   only : mkrdim
 54 
 55 !This section has been created automatically by the script Abilint (TD).
 56 !Do not modify the following lines by hand.
 57 #undef ABI_FUNC
 58 #define ABI_FUNC 'dfpt_symph'
 59  use interfaces_14_hidewrite
 60  use interfaces_32_util
 61 !End of the abilint section
 62 
 63  implicit none
 64 
 65 !Arguments ------------------------------------
 66 !scalars
 67  integer,intent(in) :: iout,natom,nsym
 68 !arrays
 69  integer,intent(in) :: indsym(4,nsym,natom),symrel(3,3,nsym)
 70  real(dp),intent(in) :: acell(3),eigvec(2*3*natom*3*natom),phfrq(3*natom)
 71  real(dp),intent(in) :: rprim(3,3)
 72 
 73 !Local variables -------------------------
 74 !scalars
 75  integer :: iad1,iad2,iad3,iatom,idir,ii1,ii2,ii3,imode,isym,itol,jad,jatom,jj
 76  integer :: jmode,kk,ntol
 77  character(len=500) :: message
 78 !arrays
 79  integer,allocatable :: degeneracy(:),integer_characters(:),symind(:,:)
 80  real(dp) :: gprimd(3,3),rprimd(3,3)
 81  real(dp),allocatable :: eigvtr(:),redvec(:),redvtr(:),symph(:,:)
 82 
 83 !******************************************************************
 84 
 85 !Compute dimensional primitive translations rprimd and its inverse gprimd
 86  call mkrdim(acell,rprim,rprimd)
 87  call matr3inv(rprimd,gprimd)
 88 
 89 !Build the symmetry index (inverse of indsym(4,:,:))
 90  ABI_ALLOCATE(symind,(nsym,natom))
 91  do isym=1,nsym
 92    do iatom=1,natom
 93      symind(isym,indsym(4,isym,iatom))=iatom
 94    end do
 95  end do
 96 
 97  ABI_ALLOCATE(symph,(nsym,3*natom))
 98  ABI_ALLOCATE(redvec,(2*3*natom))
 99  ABI_ALLOCATE(redvtr,(2*3*natom))
100  ABI_ALLOCATE(eigvtr,(2*3*natom))
101 
102 !Loop over the vibration modes
103  do imode=1,3*natom
104 
105 !  Compute eigvec for this mode in reduced coordinates redvec
106    do iatom=1,natom
107      iad1=3*(iatom-1)+1
108      ii1=2*3*natom*(imode-1)+2*(iad1-1)+1
109      iad2=3*(iatom-1)+2
110      ii2=2*3*natom*(imode-1)+2*(iad2-1)+1
111      iad3=3*(iatom-1)+3
112      ii3=2*3*natom*(imode-1)+2*(iad3-1)+1
113      do idir=1,3
114        jad=3*(iatom-1)+idir
115        jj=2*(jad-1)+1
116        redvec(jj)=gprimd(1,idir)*eigvec(ii1)+&
117 &       gprimd(2,idir)*eigvec(ii2)+&
118 &       gprimd(3,idir)*eigvec(ii3)
119        redvec(jj+1)=gprimd(1,idir)*eigvec(ii1+1)+&
120 &       gprimd(2,idir)*eigvec(ii2+1)+&
121 &       gprimd(3,idir)*eigvec(ii3+1)
122      end do !idir
123    end do !iatom
124 
125 !  Apply each transformation to redvec and store at the correct location in redvtr (iatom -> jatom)
126    do isym=1,nsym
127      do iatom=1,natom
128        jatom=symind(isym,iatom)
129        iad1=3*(iatom-1)+1
130        ii1=2*(iad1-1)+1
131        iad2=3*(iatom-1)+2
132        ii2=2*(iad2-1)+1
133        iad3=3*(iatom-1)+3
134        ii3=2*(iad3-1)+1
135        do idir=1,3
136          jad=3*(jatom-1)+idir
137          jj=2*(jad-1)+1
138          redvtr(jj)=dble(symrel(idir,1,isym))*redvec(ii1)+&
139 &         dble(symrel(idir,2,isym))*redvec(ii2)+&
140 &         dble(symrel(idir,3,isym))*redvec(ii3)
141          redvtr(jj+1)=dble(symrel(idir,1,isym))*redvec(ii1+1)+&
142 &         dble(symrel(idir,2,isym))*redvec(ii2+1)+&
143 &         dble(symrel(idir,3,isym))*redvec(ii3+1)
144 
145        end do !idir
146      end do !iatom
147 
148 !    Compute redvtr in cartesian coordinates eigvtr
149      do iatom=1,natom
150        iad1=3*(iatom-1)+1
151        ii1=2*(iad1-1)+1
152        iad2=3*(iatom-1)+2
153        ii2=2*(iad2-1)+1
154        iad3=3*(iatom-1)+3
155        ii3=2*(iad3-1)+1
156        do idir=1,3
157          jad=3*(iatom-1)+idir
158          jj=2*(jad-1)+1
159          eigvtr(jj)=rprimd(idir,1)*redvtr(ii1)+&
160 &         rprimd(idir,2)*redvtr(ii2)+&
161 &         rprimd(idir,3)*redvtr(ii3)
162          eigvtr(jj+1)=rprimd(idir,1)*redvtr(ii1+1)+&
163 &         rprimd(idir,2)*redvtr(ii2+1)+&
164 &         rprimd(idir,3)*redvtr(ii3+1)
165        end do !idir
166      end do !iatom
167 
168 !    Compute scalar product...
169      symph(isym,imode)=0.0_dp
170      do jad=1,3*natom
171        jj=2*(jad-1)+1
172        kk=2*3*natom*(imode-1)+2*(jad-1)+1
173        symph(isym,imode)=symph(isym,imode)+eigvtr(jj)*eigvec(kk)+eigvtr(jj+1)*eigvec(kk+1)
174      end do
175 
176    end do !isym
177  end do !imode
178 
179 !Treat degeneracies (different tolerances will be tried)
180 !Compute the order of the degeneracy, and
181 !attribute it to the lowest of the degenerate modes
182 !Also attribute the characters to the lowest mode
183 !When all the characters are integers, consider that the
184 !mode is non-degenerate. The maximum difference in frequency
185 !that is tolerated is on the order of 4cm-1 (which is large...)
186  ABI_ALLOCATE(degeneracy,(3*natom))
187  ABI_ALLOCATE(integer_characters,(3*natom))
188  degeneracy(:)=1
189  integer_characters(:)=0
190  do itol=1,20
191    ntol=itol
192    do imode=3*natom,2,-1
193      if(integer_characters(imode)==0)then
194        do jmode=imode-1,1,-1
195          if(integer_characters(jmode)==0)then
196            if(abs(phfrq(imode)-phfrq(jmode))<itol*tol6)then
197              degeneracy(jmode)=degeneracy(jmode)+degeneracy(imode)
198              degeneracy(imode)=0
199              symph(:,jmode)=symph(:,jmode)+symph(:,imode)
200              symph(:,imode)=0.0_dp
201            end if
202          end if !integer_characters(jmode)==0
203        end do !jmode
204      end if !integer_characters(imode)==0
205    end do !imode
206    do imode=1,3*natom
207      if(maxval(abs( symph(:,imode)-nint(symph(:,imode)) ))<0.05_dp)then
208        integer_characters(imode)=1
209      end if
210    end do
211    if(sum(integer_characters(:))==3*natom)exit
212  end do !itol
213 
214 !DEBUG
215 !write(std_out,*)' dfpt_symph : degeneracy=',degeneracy(:)
216 !ENDDEBUG
217 
218  write(message,'(a,a,es8.2,a)')ch10,&
219 & ' Analysis of degeneracies and characters (maximum tolerance=',ntol*tol6,' a.u.)'
220  call wrtout(iout,message,'COLL')
221  call wrtout(std_out,message,'COLL')
222 
223  do imode=1,3*natom
224    if(degeneracy(imode)/=0)then
225      write(message,'(a,i4)') ' Symmetry characters of vibration mode #',imode
226      call wrtout(iout,message,'COLL')
227      call wrtout(std_out,message,'COLL')
228      if(degeneracy(imode)>=2)then
229        if(degeneracy(imode)==2) write(message,'(a,i4)') &
230 &       '        degenerate with vibration mode #',imode+1
231        if(degeneracy(imode)>=3) write(message,'(a,i4,a,i4)') &
232 &       '       degenerate with vibration modes #',imode+1,' to ',imode+degeneracy(imode)-1
233        call wrtout(iout,message,'COLL')
234        call wrtout(std_out,message,'COLL')
235      end if
236      do jj=1,(nsym-1)/16+1
237        write(message,'(16f5.1)') (symph(isym,imode),isym=(jj-1)*16+1,min(nsym,jj*16))
238        call wrtout(iout,message,'COLL')
239        call wrtout(std_out,message,'COLL')
240      end do
241    end if
242  end do !imode
243 
244  ABI_DEALLOCATE(degeneracy)
245  ABI_DEALLOCATE(integer_characters)
246  ABI_DEALLOCATE(eigvtr)
247  ABI_DEALLOCATE(redvtr)
248  ABI_DEALLOCATE(redvec)
249  ABI_DEALLOCATE(symph)
250  ABI_DEALLOCATE(symind)
251 
252 end subroutine dfpt_symph