TABLE OF CONTENTS


ABINIT/m_atprj [ Modules ]

[ Top ] [ Modules ]

NAME

 m_atprj

FUNCTION

 Module to output atomic projections of phonon modes

COPYRIGHT

 Copyright (C) 2011-2018 ABINIT group (MJV)
 This file is distributed under the terms of the
 GNU General Public Licence, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

SOURCE

19 #if defined HAVE_CONFIG_H
20 #include "config.h"
21 #endif
22 
23 #include "abi_common.h"
24 
25 module m_atprj
26 
27  use defs_basis
28  use m_abicore
29  use m_errors
30 
31  use m_io_tools, only : get_unit, open_file
32  use m_fstrings, only : int2char4
33 
34  implicit none
35 
36  private

m_atprj/atprj_destroy [ Functions ]

[ Top ] [ m_atprj ] [ Functions ]

NAME

 atprj_destroy

FUNCTION

 deallocate atomic projection datastructure and close files

 INPUT

OUTPUT

 t_atprj = container object for atomic projections

PARENTS

      m_phonons

CHILDREN

SOURCE

236 subroutine atprj_destroy(t_atprj)
237 
238 
239 !This section has been created automatically by the script Abilint (TD).
240 !Do not modify the following lines by hand.
241 #undef ABI_FUNC
242 #define ABI_FUNC 'atprj_destroy'
243 !End of the abilint section
244 
245  implicit none
246 
247  type(atprj_type), intent(inout) :: t_atprj
248 
249  if (allocated(t_atprj%iatprj_bs)) then
250    ABI_DEALLOCATE(t_atprj%iatprj_bs)
251  end if
252 
253  if (allocated(t_atprj%filename)) then
254    ABI_DEALLOCATE(t_atprj%filename)
255  end if
256 
257 end subroutine atprj_destroy
258 
259 end module m_atprj

m_atprj/atprj_init [ Functions ]

[ Top ] [ m_atprj ] [ Functions ]

NAME

 atprj_init

FUNCTION

 initialize atprj datastructure

 INPUT
 natom = number of atoms
 natprj_bs = number of atoms to project on
 iatprj_bs = indices of atoms to project on
 outfile_radix = base file name for output files

OUTPUT

 t_atprj = container object for atomic projections

PARENTS

      m_phonons

CHILDREN

SOURCE

 89 subroutine atprj_init(t_atprj, natom, natprj_bs, iatprj_bs, outfile_radix)
 90 
 91 
 92 !This section has been created automatically by the script Abilint (TD).
 93 !Do not modify the following lines by hand.
 94 #undef ABI_FUNC
 95 #define ABI_FUNC 'atprj_init'
 96 !End of the abilint section
 97 
 98  implicit none
 99 
100  type(atprj_type), intent(out) :: t_atprj
101  integer, intent(in) :: natom
102  integer, intent(in) :: natprj_bs
103  integer, intent(in) :: iatprj_bs(natprj_bs)
104  character(len=*), intent(in) :: outfile_radix
105 
106 !Local variables-------------------------------
107 !scalars
108  integer :: iatom, imode, iunit
109  character(len=10) :: imodestring, iatomstring
110  character(len=500) :: msg
111 ! *************************************************************************
112 
113  t_atprj%natprj_bs = natprj_bs
114  t_atprj%natom = natom
115 
116  ABI_ALLOCATE(t_atprj%iatprj_bs,(natprj_bs))
117  t_atprj%iatprj_bs = iatprj_bs
118 
119 ! for each phonon mode and atom for projection, open a file
120  ABI_ALLOCATE(t_atprj%filename ,(3*natom,natprj_bs))
121  iunit = get_unit()
122  do imode = 1, 3*natom
123    call int2char4(imode, imodestring)
124    ABI_CHECK((imodestring(1:1)/='#'),'Bug: string length too short!')
125    do iatom = 1, natprj_bs
126      call int2char4(iatom, iatomstring)
127      ABI_CHECK((iatomstring(1:1)/='#'),'Bug: string length too short!')
128      t_atprj%filename(imode,iatom) = trim(outfile_radix)//"_mod"//trim(imodestring)//"_iat"//trim(iatomstring)
129      if (open_file(t_atprj%filename(imode,iatom), msg, newunit=iunit, form="formatted", action="write") /= 0) then
130        MSG_ERROR(msg)
131      end if
132      ! print header
133      write (unit=iunit, fmt='(a)') '##'
134      write (unit=iunit, fmt='(a,I6,a)') '##  This file contains abinit phonon frequencies for mode number ', &
135 &          imode, ' along a path in reciprocal space,'
136      write (unit=iunit, fmt='(a,I6)') '##  the third column is the projection along atom number ',&
137 &          t_atprj%iatprj_bs(iatom)
138      write (unit=iunit, fmt='(a)') '##'
139 
140      close (iunit)
141    end do
142  end do
143 
144 end subroutine atprj_init

m_atprj/atprj_print [ Functions ]

[ Top ] [ m_atprj ] [ Functions ]

NAME

 atprj_print

FUNCTION

 print out 1 line per atomic projection file

 INPUT
 t_atprj = container object for atomic projections
 phfrq = phonon frequencies for present q point
 eigvec = eigenvectors for present q point

OUTPUT

  writes to files

PARENTS

      m_phonons

CHILDREN

SOURCE

170 subroutine atprj_print(t_atprj, iq, phfrq, eigvec)
171 
172 
173 !This section has been created automatically by the script Abilint (TD).
174 !Do not modify the following lines by hand.
175 #undef ABI_FUNC
176 #define ABI_FUNC 'atprj_print'
177 !End of the abilint section
178 
179  implicit none
180 
181 !arguments
182  integer, intent(in) :: iq
183  type(atprj_type), intent(in) :: t_atprj
184  real(dp), intent(in) :: phfrq(3*t_atprj%natom)
185  real(dp), intent(in) :: eigvec(2,3,t_atprj%natom,3,t_atprj%natom)
186 
187 !local variables
188  integer :: jatom, idir, iatom, imode, iunit, jdir
189  real(dp) :: proj
190 
191 ! *************************************************************************
192 
193  iunit = get_unit()
194  do iatom = 1, t_atprj%natom
195    do idir = 1, 3
196      imode = idir + 3*(iatom-1)
197      do jatom = 1, t_atprj%natprj_bs
198        open (unit=iunit, file=t_atprj%filename(imode,jatom), position='append')
199        write (unit=iunit, fmt='(a,I4,a)') '# atom ', jatom, ' sum of directions'
200        proj = sum(eigvec(:,:,jatom,idir,iatom)**2)
201        write (unit=iunit, fmt='(I6,2E20.10)') iq, phfrq(imode), proj
202 
203        do jdir = 1, 3
204          write (unit=iunit, fmt='(2a,I4,a,I4)') ch10, '# atom ', jatom, ' directions ', jdir
205          proj = sum(eigvec(:,jdir,jatom,idir,iatom)**2)
206          write (unit=iunit, fmt='(I6,2E20.10)') iq, phfrq(imode), proj
207        end do
208        close (iunit)
209      end do
210    end do
211  end do
212 
213 end subroutine atprj_print

m_atprj/atprj_type [ Types ]

[ Top ] [ m_atprj ] [ Types ]

NAME

 atprj_type

FUNCTION

 Container for atomic projection file data

SOURCE

48 type, public :: atprj_type
49 
50   integer :: natprj_bs
51   integer :: natom
52 
53   integer, allocatable :: iatprj_bs(:)
54   character(len=fnlen), allocatable :: filename(:,:)
55 
56 end type atprj_type
57 
58 public :: atprj_init
59 public :: atprj_print
60 public :: atprj_destroy
61 
62 contains