TABLE OF CONTENTS


ABINIT/m_sorth_ph [ Modules ]

[ Top ] [ Modules ]

NAME

FUNCTION

COPYRIGHT

  Copyright (C) 2008-2018 ABINIT group (MVer, FDortu, MVeithen)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

SOURCE

17 #if defined HAVE_CONFIG_H
18 #include "config.h"
19 #endif
20 
21 #include "abi_common.h"
22 
23 module m_sortph
24 
25  use defs_basis
26  use m_abicore
27  use m_errors
28 
29  use m_io_tools,   only : open_file
30 
31  implicit none
32 
33  private
34 
35  complex(dpc),save,allocatable :: eigvecLast(:,:)
36 
37  public :: end_sortph
38  public :: sortph
39 
40  ! Logical units used to write data.
41  integer,private,save :: udispl=-1,ufreq=-1

m_sortph/end_sortph [ Functions ]

[ Top ] [ Functions ]

NAME

 end_sortph

FUNCTION

 Deallocate array for sortph

INPUTS

OUTPUT

  Only deallocation

NOTES

PARENTS

      harmonic_thermo,m_phonons

CHILDREN

SOURCE

66 subroutine end_sortph()
67 
68 
69 !This section has been created automatically by the script Abilint (TD).
70 !Do not modify the following lines by hand.
71 #undef ABI_FUNC
72 #define ABI_FUNC 'end_sortph'
73 !End of the abilint section
74 
75  if (allocated(eigvecLast))  then
76    ABI_DEALLOCATE(eigvecLast)
77  end if
78 
79  if (ufreq /= -1) then
80    close(ufreq); ufreq = -1
81  end if
82  if (udispl /= -1) then
83    close(udispl); udispl = -1
84  end if
85 
86 end subroutine end_sortph

m_sortph/sortph [ Functions ]

[ Top ] [ Functions ]

NAME

 sortph

FUNCTION

 Sort the energies in order to have fine phonon dispersion curves
 It is best not to include the gamma point in the list

 MODIFIED
 Takeshi Nishimatsu

INPUTS

  eigvec(2*3*natom*3*natom)= contain
  the eigenvectors of the dynamical matrix.
  displ(2*3*natom*3*natom)= contain
   the displacements of atoms in cartesian coordinates.
   The first index means either the real or the imaginary part,
   The second index runs on the direction and the atoms displaced
   The third index runs on the modes.
  filnam=name of output files
   hacmm1,hartev,harthz,xkb= different conversion factors
  natom= number of atom
  phfrq(3*natom)= phonon frequencies in Hartree

OUTPUT

  (only writing ?)

NOTES

 Called by one processor only

PARENTS

      harmonic_thermo,m_phonons

CHILDREN

SOURCE

126 subroutine sortph(eigvec,displ,filnam, natom,phfrq)
127 
128 
129 !This section has been created automatically by the script Abilint (TD).
130 !Do not modify the following lines by hand.
131 #undef ABI_FUNC
132 #define ABI_FUNC 'sortph'
133 !End of the abilint section
134 
135 implicit none
136 
137 !Arguments -----------------------------------
138 !scalars
139 integer,intent(in) :: natom
140 character(len=*),intent(in) :: filnam
141 !arrays
142 real(dp),intent(in) :: eigvec(2,3,natom,3,natom)
143 real(dp),intent(in) :: displ(2*3*natom*3*natom)
144 real(dp),intent(in) :: phfrq(3*natom)
145 
146 !Local variables-------------------------------
147 !scalars
148 integer :: iatom,imode,j,idir1,idir2,ipert1,ipert2,i1,i2
149 character(len=fnlen) :: file_displ,file_freq
150 character(len=20) :: fmt_phfrq
151 character(len=500) :: msg
152 !arrays
153 logical     ::               mask(3*natom)
154 real(dp)    ::           phfrqNew(3*natom)
155 complex(dpc) ::           displIn(3*natom,3*natom)
156 complex(dpc) ::           displNew(3*natom,3*natom)
157 complex(dpc) ::           eigvecIn(3*natom,3*natom)
158 complex(dpc) ::          eigvecNew(3*natom,3*natom)
159 complex(dpc) ::   transpose_eigvec(3*natom,3*natom)
160 real(dp)    ::     abs_similarity(3*natom,3*natom)  !|<displNew|displLast>|
161 ! *********************************************************************
162 
163 do ipert2=1,natom
164   do idir2=1,3
165     i2=idir2+(ipert2-1)*3
166     do ipert1=1,natom
167       do idir1=1,3
168         i1=idir1+(ipert1-1)*3
169         eigvecIn(i1,i2)=cmplx(eigvec(1,idir1,ipert1,idir2,ipert2),eigvec(2,idir1,ipert1,idir2,ipert2))
170         displIn(i1,i2)=cmplx(displ(1+2*(i1-1)+2*3*natom*(i2-1)),displ(2+2*(i1-1)+2*3*natom*(i2-1)))
171       end do
172     end do
173   end do
174 end do
175 
176  if(.not.allocated(eigvecLast)) then
177    file_freq  = trim(filnam)//".freq" !---------------------------------------------------
178    write(std_out,'(a,a)' )' sortph : opening file ',trim(file_freq)
179    if (open_file(file_freq,msg,newunit=ufreq,STATUS='replace',ACTION='write') /= 0) then
180      MSG_ERROR(msg)
181    end if
182    file_displ = trim(filnam)//".displ" !--------------------------------------------------
183    write(std_out,'(a,a)' )' sortph : opening file ',trim(file_displ)
184    if (open_file(file_displ,msg,newunit=udispl,STATUS='replace',ACTION='write') /= 0) then
185      MSG_ERROR(msg)
186    end if
187    ABI_ALLOCATE(eigvecLast,(3*natom,3*natom))
188    phfrqNew(:)   =  phfrq(:)
189    displNew(:,:) =  displIn(:,:)
190    eigvecNew(:,:) = eigvecIn(:,:)
191  else
192 !  Avoid gfortran 4.2.1 bug, with which you CANNOT conjg(transpose(displ))
193    transpose_eigvec = transpose(eigvecIn)
194    abs_similarity = abs(matmul(conjg(transpose_eigvec),eigvecLast))
195    mask(:) = .true.
196    phfrqNew(:)   =  phfrq(:)
197    displNew(:,:) =  displIn(:,:)
198    eigvecNew(:,:) = eigvecIn(:,:)
199  end if
200 
201 
202 !Write frequencies in a file
203  write(fmt_phfrq,'(a,i3,a)') '(', 3*natom, 'e18.10)'
204  write(ufreq,fmt_phfrq) (phfrqNew(j),j=1,3*natom)
205 
206 !write displacements in a file
207 ! NB: sqrt still returns a complex number could be using modulus or something simpler
208    do imode=1,3*natom
209      do iatom=1,natom
210        write(udispl,'(e18.10)') &
211        real(sqrt(displNew(3*(iatom-1)+1,imode) * conjg(displNew(3*(iatom-1)+1,imode)) + &
212 &      displNew(3*(iatom-1)+2,imode) * conjg(displNew(3*(iatom-1)+2,imode)) + &
213 &      displNew(3*(iatom-1)+3,imode) *  conjg(displNew(3*(iatom-1)+3,imode)) ))
214      end do
215    end do
216 
217  eigvecLast(:,:) = eigvecNew(:,:)
218 
219 end subroutine sortph