TABLE OF CONTENTS
ABINIT/m_sorth_ph [ Modules ]
NAME
FUNCTION
COPYRIGHT
Copyright (C) 2008-2024 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 .
SOURCE
14 #if defined HAVE_CONFIG_H 15 #include "config.h" 16 #endif 17 18 #include "abi_common.h" 19 20 module m_sortph 21 22 use defs_basis 23 use m_abicore 24 use m_errors 25 26 use m_io_tools, only : open_file 27 28 implicit none 29 30 private 31 32 complex(dpc),save,allocatable :: eigvecLast(:,:) 33 34 public :: end_sortph 35 public :: sortph 36 37 ! Logical units used to write data. 38 integer,private,save :: udispl=-1,ufreq=-1
m_sortph/end_sortph [ Functions ]
NAME
end_sortph
FUNCTION
Deallocate array for sortph
INPUTS
OUTPUT
Only deallocation
NOTES
SOURCE
58 subroutine end_sortph() 59 60 if (allocated(eigvecLast)) then 61 ABI_FREE(eigvecLast) 62 end if 63 64 if (ufreq /= -1) then 65 close(ufreq); ufreq = -1 66 end if 67 if (udispl /= -1) then 68 close(udispl); udispl = -1 69 end if 70 71 end subroutine end_sortph
m_sortph/sortph [ 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
SOURCE
106 subroutine sortph(eigvec,displ,filnam, natom,phfrq) 107 108 !Arguments ----------------------------------- 109 !scalars 110 integer,intent(in) :: natom 111 character(len=*),intent(in) :: filnam 112 !arrays 113 real(dp),intent(in) :: eigvec(2,3,natom,3,natom) 114 real(dp),intent(in) :: displ(2*3*natom*3*natom) 115 real(dp),intent(in) :: phfrq(3*natom) 116 117 !Local variables------------------------------- 118 !scalars 119 integer :: iatom,imode,j,idir1,idir2,ipert1,ipert2,i1,i2 120 character(len=fnlen) :: file_displ,file_freq 121 character(len=20) :: fmt_phfrq 122 character(len=500) :: msg 123 !arrays 124 logical :: mask(3*natom) 125 real(dp) :: phfrqNew(3*natom) 126 complex(dpc) :: displIn(3*natom,3*natom) 127 complex(dpc) :: displNew(3*natom,3*natom) 128 complex(dpc) :: eigvecIn(3*natom,3*natom) 129 complex(dpc) :: eigvecNew(3*natom,3*natom) 130 complex(dpc) :: transpose_eigvec(3*natom,3*natom) 131 real(dp) :: abs_similarity(3*natom,3*natom) !|<displNew|displLast>| 132 ! ********************************************************************* 133 134 do ipert2=1,natom 135 do idir2=1,3 136 i2=idir2+(ipert2-1)*3 137 do ipert1=1,natom 138 do idir1=1,3 139 i1=idir1+(ipert1-1)*3 140 eigvecIn(i1,i2)=cmplx(eigvec(1,idir1,ipert1,idir2,ipert2),eigvec(2,idir1,ipert1,idir2,ipert2)) 141 displIn(i1,i2)=cmplx(displ(1+2*(i1-1)+2*3*natom*(i2-1)),displ(2+2*(i1-1)+2*3*natom*(i2-1))) 142 end do 143 end do 144 end do 145 end do 146 147 if(.not.allocated(eigvecLast)) then 148 file_freq = trim(filnam)//".freq" !--------------------------------------------------- 149 write(std_out,'(a,a)' )' sortph : opening file ',trim(file_freq) 150 if (open_file(file_freq,msg,newunit=ufreq,STATUS='replace',ACTION='write') /= 0) then 151 ABI_ERROR(msg) 152 end if 153 file_displ = trim(filnam)//".displ" !-------------------------------------------------- 154 write(std_out,'(a,a)' )' sortph : opening file ',trim(file_displ) 155 if (open_file(file_displ,msg,newunit=udispl,STATUS='replace',ACTION='write') /= 0) then 156 ABI_ERROR(msg) 157 end if 158 ABI_MALLOC(eigvecLast,(3*natom,3*natom)) 159 phfrqNew(:) = phfrq(:) 160 displNew(:,:) = displIn(:,:) 161 eigvecNew(:,:) = eigvecIn(:,:) 162 else 163 ! Avoid gfortran 4.2.1 bug, with which you CANNOT conjg(transpose(displ)) 164 transpose_eigvec = transpose(eigvecIn) 165 abs_similarity = abs(matmul(conjg(transpose_eigvec),eigvecLast)) 166 mask(:) = .true. 167 phfrqNew(:) = phfrq(:) 168 displNew(:,:) = displIn(:,:) 169 eigvecNew(:,:) = eigvecIn(:,:) 170 end if 171 172 173 !Write frequencies in a file 174 write(fmt_phfrq,'(a,i3,a)') '(', 3*natom, 'e18.10)' 175 write(ufreq,fmt_phfrq) (phfrqNew(j),j=1,3*natom) 176 177 !write displacements in a file 178 ! NB: sqrt still returns a complex number could be using modulus or something simpler 179 do imode=1,3*natom 180 do iatom=1,natom 181 write(udispl,'(e18.10)') & 182 real(sqrt(displNew(3*(iatom-1)+1,imode) * conjg(displNew(3*(iatom-1)+1,imode)) + & 183 & displNew(3*(iatom-1)+2,imode) * conjg(displNew(3*(iatom-1)+2,imode)) + & 184 & displNew(3*(iatom-1)+3,imode) * conjg(displNew(3*(iatom-1)+3,imode)) )) 185 end do 186 end do 187 188 eigvecLast(:,:) = eigvecNew(:,:) 189 190 end subroutine sortph