TABLE OF CONTENTS


ABINIT/m_spectra [ Modules ]

[ Top ] [ Modules ]

NAME

 m_spectra

FUNCTION

 This module contains the definition of the specta_type data type 
 used to store results related to optical spectra with or without 
 nonlocal field effects as well as the electron energy loss function
 for a given q. These quantities are obtained from the dielectric 
 matrix as calculated in the GW part of ABINIT (screening.F90)

COPYRIGHT

 Copyright (C) 2008-2018 ABINIT group (MG)
 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

CHILDREN

NOTES

SOURCE

27 #if defined HAVE_CONFIG_H
28 #include "config.h"
29 #endif
30 
31 #include "abi_common.h"
32 
33 MODULE m_spectra
34 
35  use defs_basis
36  use m_errors
37  use m_abicore
38 
39  use m_io_tools,      only : open_file
40  use m_fstrings,      only : strcat
41 
42  implicit none
43 
44  private

m_spectra/dump_Qlist [ Functions ]

[ Top ] [ m_spectra ] [ Functions ]

NAME

  dump_Qlist

FUNCTION

  Helper function used to write the list of q-points used for the long-wavelength limit.

INPUTS

  Spectra=The Object containing the spectra

OUTPUT

  Only writing.

PARENTS

      m_spectra

CHILDREN

SOURCE

341  subroutine dump_Qlist()
342 
343 
344 !This section has been created automatically by the script Abilint (TD).
345 !Do not modify the following lines by hand.
346 #undef ABI_FUNC
347 #define ABI_FUNC 'dump_Qlist'
348 !End of the abilint section
349 
350   integer :: iqpt
351   write(unt,'(a,i3)')'# Q-point list, No. ',Spectra%nqpts
352   do iqpt=1,Spectra%nqpts
353     write(unt,'(a,i3,a,3f9.6,a)')'# ',iqpt,')  [',Spectra%qpts(:,iqpt),'] r.l.u. '
354   end do
355  end subroutine dump_Qlist
356 
357 end subroutine spectra_write 

m_spectra/spectra_free [ Functions ]

[ Top ] [ m_spectra ] [ Functions ]

NAME

  spectra_free

FUNCTION

  Deallocate all associated pointers defined in the structure.

INPUTS

OUTPUT

PARENTS

      m_screening,screening

CHILDREN

SOURCE

181 subroutine spectra_free(Spectra)
182 
183 
184 !This section has been created automatically by the script Abilint (TD).
185 !Do not modify the following lines by hand.
186 #undef ABI_FUNC
187 #define ABI_FUNC 'spectra_free'
188 !End of the abilint section
189 
190  implicit none
191 
192 !This section has been created automatically by the script Abilint (TD).
193 !Do not modify the following lines by hand.
194 #undef ABI_FUNC
195 #define ABI_FUNC 'spectra_free'
196 !End of the abilint section
197 
198 !Arguments ------------------------------------
199  type(spectra_t),intent(inout) :: spectra
200 
201 ! *********************************************************************
202 
203  if (allocated(Spectra%omega)) then
204    ABI_FREE(Spectra%omega)
205  end if
206  if (allocated(Spectra%qpts)) then
207    ABI_FREE(Spectra%qpts)
208  end if
209  if (allocated(Spectra%emacro_lf)) then
210    ABI_FREE(Spectra%emacro_lf)
211  end if
212  if (allocated(Spectra%emacro_nlf)) then
213    ABI_FREE(Spectra%emacro_nlf)
214  end if
215  if (allocated(Spectra%eelf)) then
216    ABI_FREE(Spectra%eelf)
217  end if
218 
219 end subroutine spectra_free 

m_spectra/spectra_init [ Functions ]

[ Top ] [ m_spectra ] [ Functions ]

NAME

  spectra_init

FUNCTION

  Initialize the object.

INPUTS

OUTPUT

PARENTS

      m_screening

CHILDREN

SOURCE

119 subroutine spectra_init(Spectra,nomega,omega,nqpts,qpts)
120 
121 
122 !This section has been created automatically by the script Abilint (TD).
123 !Do not modify the following lines by hand.
124 #undef ABI_FUNC
125 #define ABI_FUNC 'spectra_init'
126 !End of the abilint section
127 
128  implicit none
129 
130 !This section has been created automatically by the script Abilint (TD).
131 !Do not modify the following lines by hand.
132 #undef ABI_FUNC
133 #define ABI_FUNC 'spectra_init'
134 !End of the abilint section
135 
136 !Arguments ------------------------------------
137 !scalars
138  integer,intent(in) :: nomega,nqpts
139 !arrays
140  real(dp),intent(in) :: omega(nomega),qpts(3,nqpts)
141  type(spectra_t),intent(out) :: Spectra
142 
143 ! *********************************************************************
144 
145  Spectra%nomega = nomega
146  Spectra%nqpts  = nqpts
147 
148  ABI_MALLOC(Spectra%qpts,(3,nqpts))
149  Spectra%qpts = qpts
150 
151  ABI_MALLOC(Spectra%omega,(nomega))
152  Spectra%omega = omega
153 
154  ABI_CALLOC(Spectra%emacro_lf,(nomega,nqpts))
155  ABI_CALLOC(Spectra%emacro_nlf,(nomega,nqpts))
156  ABI_CALLOC(Spectra%eelf,(nomega,nqpts))
157 
158 end subroutine spectra_init

m_spectra/spectra_repr [ Functions ]

[ Top ] [ m_spectra ] [ Functions ]

NAME

  spectra_repr

FUNCTION

  Returns a string reporting info on the calculated dielectric constant.

INPUTS

  Spectra=The Object containing the spectra

OUTPUT

PARENTS

      m_screening,screening

CHILDREN

SOURCE

381 subroutine spectra_repr(Spectra,str)
382 
383 
384 !This section has been created automatically by the script Abilint (TD).
385 !Do not modify the following lines by hand.
386 #undef ABI_FUNC
387 #define ABI_FUNC 'spectra_repr'
388 !End of the abilint section
389 
390  implicit none
391 
392 !This section has been created automatically by the script Abilint (TD).
393 !Do not modify the following lines by hand.
394 #undef ABI_FUNC
395 #define ABI_FUNC 'spectra_repr'
396 !End of the abilint section
397 
398 !Arguments ------------------------------------
399 !scalars
400  type(spectra_t),intent(in) :: Spectra
401  character(len=*),intent(out) :: str
402 
403 !Local variables-------------------------------
404 !scalars
405  integer :: iqpt
406  real(dp) :: epsilon0,epsilon0_nlf
407  character(len=500) :: msg
408 
409 ! *********************************************************************
410 
411  !istatic = -1
412  !do io = Spectra%nomega 
413  ! if (ABS(REAL(Spectra%omega(io)))<1.e-3.and.ABS(AIMAG(Spectra%omega(io)))<1.e-3) then
414  !  istatic = io 
415  !  EXIT
416  ! end if
417  !end do
418 
419  str = ""
420  do iqpt=1,Spectra%nqpts
421    epsilon0    = REAL(Spectra%emacro_lf (1,iqpt))
422    epsilon0_nlf= REAL(Spectra%emacro_nlf(1,iqpt))
423    write(msg,'(a,3f9.6,a)')' For q-point: ',Spectra%qpts(:,iqpt),ch10
424    str = strcat(str,msg)
425    write(msg,'(1x,a,f8.4,a)')' dielectric constant = ',epsilon0,ch10
426    str = strcat(str,msg)
427    write(msg,'(1x,a,f8.4,a)')' dielectric constant without local fields = ',epsilon0_nlf,ch10
428    str = strcat(str,msg)
429  end do
430 
431 end subroutine spectra_repr
432 
433 !----------------------------------------------------------------------
434 
435 END MODULE m_spectra

m_spectra/spectra_t [ Types ]

[ Top ] [ m_spectra ] [ Types ]

NAME

 spectra_t

FUNCTION

  Object used to store optical spectra with or without non-local field effects.

SOURCE

56  type,public :: spectra_t
57 
58  !scalars
59   integer :: nomega
60   ! number of frequencies
61   
62   integer :: nqpts
63   ! number of q-points
64 
65  !arrays
66   real(dp),allocatable :: omega(:)
67   ! omega(nomega)
68   ! Real frequency mesh for optical spectra.
69 
70   real(dp),allocatable :: qpts(:,:) 
71   ! qpts(3,nqpoints)
72   ! The list of q-points used for the spectra
73 
74   real(dp),allocatable :: eelf(:,:)
75   ! eelf(nomega,nqpoints)
76   ! contains the Electron Energy Loss Function i.e. -\Im{ e^{-1}_{G1=0,G2=0}(q-->0,nomega)}
77 
78   complex(dpc),allocatable :: emacro_lf(:,:) 
79   ! emacro_lf(nomega,nqpoints)
80   ! contains 1/e^{-1}_{G1=0,G2=0}(q-->0,nomega) (with Local field effects)
81 
82   complex(dpc),allocatable :: emacro_nlf(:,:)
83   ! emacro_nlf(nomega,nqpoints)
84   ! contains e_{G1=0,G2=0}(q-->0,nomega) (without Local field effects)
85 
86  end type spectra_t

m_spectra/spectra_write [ Functions ]

[ Top ] [ m_spectra ] [ Functions ]

NAME

  spectra_write

FUNCTION

  Write the optical spectra stored in the object on an external formatted file.

INPUTS

  Spectra=The Object containing the spectra
  write_bits=Positive integer defining the quantities to be written (bit representation is used)
  fname=Name of the file to be written.

OUTPUT

PARENTS

      screening

CHILDREN

SOURCE

245 subroutine spectra_write(Spectra,write_bits,fname)
246 
247 
248 !This section has been created automatically by the script Abilint (TD).
249 !Do not modify the following lines by hand.
250 #undef ABI_FUNC
251 #define ABI_FUNC 'spectra_write'
252 !End of the abilint section
253 
254  implicit none
255 
256 !Arguments ------------------------------------
257 !scalars
258  integer,intent(in) :: write_bits
259  character(len=*),intent(in) :: fname
260 !arrays
261  type(spectra_t),intent(in) :: Spectra
262 
263 !Local variables-------------------------------
264 !scalars
265  integer :: unt,io,iqpt
266  real(dp) :: mino,maxo
267  character(len=100) :: fmt
268  character(len=500) :: msg
269 
270 ! *********************************************************************
271 
272  if (write_bits<0) RETURN
273 
274  mino = MINVAL(Spectra%omega)
275  maxo = MAXVAL(Spectra%omega)
276 
277  if (open_file(fname,msg,newunit=unt) /= 0) then
278    MSG_ERROR(msg)
279  end if
280 
281  !write(unt,'(a,i5,2(a,f9.1),a)')'# nomega : ',Spectra%nomega,' from ',mino*Ha_eV,' up to ',maxo*Ha_eV,' [eV] '
282 
283  if ( IAND(write_bits,W_EM_NLF ) == W_EM_NLF ) then
284    write(unt,'(a)')'#'
285    write(unt,'(a)')'# Macroscopic Dielectric Function without local fields'
286    call dump_qlist()
287    write(unt,'(a)')'# Omega [eV]    Re epsilon_M       IM eps_M '
288    write(fmt,'(a,i3,a)') '(1x,f7.3,7x,',Spectra%nqpts,'(2(e12.4),2x))'
289    do io=1,Spectra%nomega
290      write(unt,fmt) Spectra%omega(io)*Ha_eV, ( Spectra%emacro_nlf(io,iqpt), iqpt=1,Spectra%nqpts )
291    end do
292  end if
293  !
294  if ( IAND(write_bits,W_EM_LF ) == W_EM_LF ) then 
295    write(unt,'(a)')'#'
296    write(unt,'(a)')'# Macroscopic Dielectric Function with local fields included'
297    call dump_qlist()
298    write(unt,'(a)')'# Omega [eV]    Re epsilon_M       Im eps_M '
299    write(fmt,'(a,i3,a)') '(1x,f7.3,7x,',Spectra%nqpts,'(2(e12.4),2x))'
300    do io=1,Spectra%nomega
301      write(unt,fmt) Spectra%omega(io)*Ha_eV, ( Spectra%emacro_lf(io,iqpt), iqpt=1,Spectra%nqpts )
302    end do
303  end if
304  !
305  if ( IAND(write_bits,W_EELF) == W_EELF) then
306    write(unt,'(a)')'#'
307    write(unt,'(a)')'# Electron Energy Loss Function -Im(1/epsilon_M)'
308    call dump_qlist()
309    write(unt,'(a)')'# Omega [eV]    -Im(1/epsilon_M)'
310    write(fmt,'(a,i3,a)') '(1x,f7.3,7x,',Spectra%nqpts,'(e12.4,2x))'
311    do io=1,Spectra%nomega
312      write(unt,fmt) Spectra%omega(io)*Ha_eV, ( Spectra%eelf(io,iqpt), iqpt=1,Spectra%nqpts ) ! -AIMAG(chi0(1,1,io))
313    end do
314  end if
315 
316  close(unt)
317 
318 CONTAINS