TABLE OF CONTENTS


ABINIT/m_eprenorms [ Modules ]

[ Top ] [ Modules ]

NAME

 m_eprenorms

FUNCTION

 This module contains datatypes to compute the renormalization of electronic states due to
 eph coupling and temperature effects

NOTES

 This code is still under development and the API will change in the next versions.
 Contact gmatteo

COPYRIGHT

 Copyright (C) 2001-2018 ABINIT group (YG)
 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

22 #if defined HAVE_CONFIG_H
23 #include "config.h"
24 #endif
25 
26 #include "abi_common.h"
27 
28 module m_eprenorms
29 
30  use defs_basis
31  use defs_abitypes
32  use defs_datatypes
33  use m_errors
34  use m_xmpi
35 #ifdef HAVE_NETCDF
36  use netcdf
37 #endif
38  use m_nctk
39 
40  use m_crystal,  only : crystal_t
41  use m_kpts,     only : listkk
42 
43  implicit none
44 
45  private

m_eprenorms/eprenorms_bcast [ Functions ]

[ Top ] [ m_eprenorms ] [ Functions ]

NAME

 eprenorms_bcast

FUNCTION

  MPI broadcast all the content of eprenorms_t

INPUTS

  master = Rank of master
  comm = MPI communicator

SIDE EFFECTS

  Epren<eprenorms_t> = Data broadcasted on every node from master

PARENTS

      optic,setup_bse

CHILDREN

      listkk

SOURCE

325 subroutine eprenorms_bcast(Epren,master,comm)
326 
327 
328 !This section has been created automatically by the script Abilint (TD).
329 !Do not modify the following lines by hand.
330 #undef ABI_FUNC
331 #define ABI_FUNC 'eprenorms_bcast'
332 !End of the abilint section
333 
334  implicit none
335 
336 !Arguments -----------------------------------
337 !scalars
338  integer,intent(in) :: master, comm
339  type(eprenorms_t),intent(inout) :: Epren
340 
341 !Local variables------------------------------
342 !scalars
343  integer :: ierr
344 
345 ! ************************************************************************
346 
347  if (xmpi_comm_size(comm) == 1) return
348 
349  call xmpi_bcast(Epren%nkpt, master, comm, ierr)
350  call xmpi_bcast(Epren%nsppol, master, comm, ierr)
351  call xmpi_bcast(Epren%mband, master, comm, ierr)
352  call xmpi_bcast(Epren%ntemp, master, comm, ierr)
353 
354  if (xmpi_comm_rank(comm) /= master) then
355   call eprenorms_init(Epren, Epren%nkpt, Epren%nsppol, Epren%mband, Epren%ntemp)
356  end if
357 
358  call xmpi_bcast(Epren%kpts, master, comm, ierr)
359  call xmpi_bcast(Epren%temps, master, comm, ierr)
360  call xmpi_bcast(Epren%eigens, master, comm, ierr)
361  call xmpi_bcast(Epren%occs, master, comm, ierr)
362  call xmpi_bcast(Epren%renorms, master, comm, ierr)
363  call xmpi_bcast(Epren%lifetimes, master, comm, ierr)
364 
365 end subroutine eprenorms_bcast

m_eprenorms/eprenorms_free [ Functions ]

[ Top ] [ m_eprenorms ] [ Functions ]

NAME

 eprenorms_free

FUNCTION

 Deallocate all memory associated with eprenorms

INPUTS

 Epren<eprenorms_t>=The datatype to be freed

PARENTS

      bethe_salpeter,optic

CHILDREN

      listkk

SOURCE

192 subroutine eprenorms_free(Epren)
193 
194 
195 !This section has been created automatically by the script Abilint (TD).
196 !Do not modify the following lines by hand.
197 #undef ABI_FUNC
198 #define ABI_FUNC 'eprenorms_free'
199 !End of the abilint section
200 
201  implicit none
202 
203 !Arguments -----------------------------------
204 !scalars
205  type(eprenorms_t),intent(inout) :: Epren
206 
207 !*********************************************************************
208 
209  if (allocated(Epren%temps)) then
210    ABI_FREE(Epren%temps)
211  end if
212  if (allocated(Epren%kpts)) then
213    ABI_FREE(Epren%kpts)
214  end if
215  if (allocated(Epren%eigens)) then
216    ABI_FREE(Epren%eigens)
217  end if
218  if (allocated(Epren%occs)) then
219    ABI_FREE(Epren%occs)
220  end if
221  if (allocated(Epren%renorms)) then
222    ABI_FREE(Epren%renorms)
223  end if
224  if (allocated(Epren%lifetimes)) then
225    ABI_FREE(Epren%lifetimes)
226  end if
227 
228 end subroutine eprenorms_free

m_eprenorms/eprenorms_from_epnc [ Functions ]

[ Top ] [ m_eprenorms ] [ Functions ]

NAME

 eprenorms_from_epnc

FUNCTION

 Allocates and initializes the datatype from a _EP.nc file

INPUTS

 filename = name of the file to be read

SIDE EFFECTS

 Epren<eprenorms_t> = fields are initialized and filled with data from filename

PARENTS

      optic,setup_bse

CHILDREN

      listkk

SOURCE

254 subroutine eprenorms_from_epnc(Epren,filename)
255 
256 
257 !This section has been created automatically by the script Abilint (TD).
258 !Do not modify the following lines by hand.
259 #undef ABI_FUNC
260 #define ABI_FUNC 'eprenorms_from_epnc'
261 !End of the abilint section
262 
263  implicit none
264 
265 !Arguments -----------------------------------
266 !scalars
267  character(len=fnlen),intent(in) :: filename
268  type(eprenorms_t),intent(inout) :: Epren
269 
270 !Local variables------------------------------
271  integer :: nkpt, mband, nsppol, ntemp
272 #ifdef HAVE_NETCDF
273  integer :: ncid
274 #endif
275 
276 ! ************************************************************************
277 
278 #ifdef HAVE_NETCDF
279  NCF_CHECK(nctk_open_read(ncid, filename, xmpi_comm_self))
280  NCF_CHECK(nctk_set_datamode(ncid))
281 
282  NCF_CHECK(nctk_get_dim(ncid, "number_of_kpoints", nkpt))
283  NCF_CHECK(nctk_get_dim(ncid, "number_of_spins",nsppol))
284  NCF_CHECK(nctk_get_dim(ncid, "max_number_of_states",mband))
285  NCF_CHECK(nctk_get_dim(ncid, "number_of_temperature",ntemp))
286 
287  call eprenorms_init(Epren, nkpt, nsppol, mband, ntemp)
288 
289  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid,"reduced_coordinates_of_kpoints"), Epren%kpts))
290  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid,"temperature"), Epren%temps))
291  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid,"eigenvalues"), Epren%eigens))
292  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid,"occupations"), Epren%occs))
293  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid,"zero_point_motion"), Epren%renorms))
294  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid,"lifetime"), Epren%lifetimes))
295 
296 #endif
297 
298 end subroutine eprenorms_from_epnc

m_eprenorms/eprenorms_init [ Functions ]

[ Top ] [ m_eprenorms ] [ Functions ]

NAME

 eprenorms_init

FUNCTION

  Initializes an eprenorms_t datatype

INPUTS

OUTPUT

  Epren<eprenorms_t>=Datatype gathering electron-phonon renormalizations

PARENTS

      m_eprenorms

CHILDREN

      listkk

SOURCE

130 subroutine eprenorms_init(Epren,nkpt,nsppol,mband,ntemp)
131 
132 
133 !This section has been created automatically by the script Abilint (TD).
134 !Do not modify the following lines by hand.
135 #undef ABI_FUNC
136 #define ABI_FUNC 'eprenorms_init'
137 !End of the abilint section
138 
139  implicit none
140 
141 !Arugments -----------------------------------
142 !scalars
143  integer,intent(in) :: nkpt, nsppol, mband, ntemp
144  type(eprenorms_t) :: Epren
145 !arrays
146 
147 !Local variables------------------------------
148 !scalars
149 !arrays
150 
151 !*************************************************************************
152 
153  DBG_ENTER("COLL")
154 
155  Epren%nkpt = nkpt
156  Epren%nsppol = nsppol
157  Epren%mband = mband
158  Epren%ntemp = ntemp
159 
160  ABI_MALLOC(Epren%temps,(Epren%ntemp))
161  ABI_MALLOC(Epren%kpts,(3,Epren%nkpt))
162  ABI_MALLOC(Epren%eigens,(Epren%mband,Epren%nkpt,Epren%nsppol))
163  ABI_MALLOC(Epren%occs,(Epren%mband,Epren%nkpt,Epren%nsppol))
164  ABI_MALLOC(Epren%renorms,(2,Epren%mband,Epren%nkpt,Epren%nsppol,Epren%ntemp))
165  ABI_MALLOC(Epren%lifetimes,(2,Epren%mband,Epren%nkpt,Epren%nsppol,Epren%ntemp))
166 
167  DBG_EXIT("COLL")
168 
169 end subroutine eprenorms_init

m_eprenorms/eprenorms_t [ Types ]

[ Top ] [ m_eprenorms ] [ Types ]

NAME

 eprenorms_t

FUNCTION

 Datatype gathering data for electron-phonon renormalization of the band structure

SOURCE

 57  type,public :: eprenorms_t
 58 
 59   !scalars
 60   integer :: nkpt
 61   ! Number of kpoints
 62 
 63   integer :: nsppol
 64   ! Number of spin channels
 65 
 66   integer :: mband
 67   ! Maximum number of bands
 68 
 69   integer :: ntemp
 70   ! Number of temperatures
 71 
 72   !arrays
 73   real(dp), allocatable :: kpts(:,:)
 74   ! kpt(3,nkpt)
 75   ! Kpoints
 76 
 77   real(dp), allocatable :: temps(:)
 78   ! temps(ntemp)
 79   ! Temperatures
 80 
 81   real(dp), allocatable :: eigens(:,:,:)
 82   ! eigens(mband,nkpt,nsppol)
 83   ! Kohn-Sham eigenvalues
 84 
 85   real(dp), allocatable :: occs(:,:,:)
 86   ! occ(mband,nkpt,nsppol)
 87   ! Occupation numbers
 88 
 89   real(dp), allocatable :: renorms(:,:,:,:,:)
 90   ! renorms(2,mband,nkpt,nsppol,ntemp)
 91   ! Renormalization of the eigenvalues for each temperature
 92 
 93   real(dp), allocatable :: lifetimes(:,:,:,:,:)
 94   ! lifetime(2,mband,nkpt,nsppol,ntemp)
 95   ! Electron-phonon induced lifetime of the eigens
 96 
 97  end type eprenorms_t
 98 
 99  public :: eprenorms_init
100  public :: eprenorms_free
101  public :: eprenorms_from_epnc
102  public :: eprenorms_bcast
103 
104  public :: renorm_bst

m_eprenorms/renorm_bst [ Functions ]

[ Top ] [ m_eprenorms ] [ Functions ]

NAME

 renorm_bst

FUNCTION

  Renormalize the band structure Bst from data contained Epren

INPUTS

  Epren<eprenorms_t> = datatype containing the elphon renormalization
  itemp = index of the temperature you want to use
  do_lifetime = .true. if we want to use imaginary eigenvalues (lifetime field)

SIDE EFFECTS

  Bst<bands_t> : eigens are changed according to epren
                 lifetime is allocated and filled with data if do_lifetime

PARENTS

      m_exc_spectra,m_haydock,optic

CHILDREN

      listkk

SOURCE

394 subroutine renorm_bst(Epren,Bst,Cryst,itemp,do_lifetime,do_check)
395 
396 
397 !This section has been created automatically by the script Abilint (TD).
398 !Do not modify the following lines by hand.
399 #undef ABI_FUNC
400 #define ABI_FUNC 'renorm_bst'
401 !End of the abilint section
402 
403  implicit none
404 
405 !Arguments -----------------------------------
406 !scalars
407  integer :: itemp
408  logical,intent(in) :: do_lifetime
409  logical,optional,intent(in) :: do_check
410  type(eprenorms_t),intent(in) :: Epren
411  type(ebands_t),intent(inout) :: Bst
412  type(crystal_t),intent(in) :: Cryst
413 
414 !Local variables------------------------------
415 !scalars
416  integer :: isppol,ikpt
417  integer :: nband1, nband_tmp
418  integer :: timrev, sppoldbl
419  integer :: ik_eph
420  real(dp) :: dksqmax
421  logical :: check
422 !arrays
423  integer,allocatable :: bs2eph(:,:)
424 
425 ! ************************************************************************
426 
427  ABI_CHECK(Bst%nsppol == Epren%nsppol, "Nsppol should be the same")
428 
429  if(do_lifetime) then
430    ABI_MALLOC(Bst%lifetime,(Bst%mband,Bst%nkpt,Bst%nsppol))
431  end if
432 
433  check = .TRUE.
434  if(present(do_check)) then
435    check = do_check
436  end if
437 
438  sppoldbl = 1 !; if (any(Cryst%symafm == -1) .and. Epren%nsppol == 1) nsppoldbl=2
439  ABI_MALLOC(bs2eph, (BSt%nkpt*sppoldbl, 6))
440  timrev = 1
441  call listkk(dksqmax, Cryst%gmet, bs2eph, Epren%kpts, BSt%kptns, Epren%nkpt, Bst%nkpt, Cryst%nsym, &
442 &   sppoldbl, Cryst%symafm, Cryst%symrel, timrev, use_symrec=.False.)
443 
444  do isppol=1,Bst%nsppol
445    do ikpt=1,Bst%nkpt
446      nband1 = Bst%nband(ikpt+(isppol-1)*Bst%nkpt)
447      nband_tmp=MIN(nband1,Epren%mband)
448 
449      ik_eph = bs2eph(ikpt,1)
450 
451      !FIXME change check
452      if (check) then
453        if (ANY(ABS(Bst%eig(1:MIN(10,nband_tmp),ikpt,isppol) - Epren%eigens(1:MIN(10,nband_tmp),ik_eph,isppol)) > tol3)) then
454          !write(stdout,*) "eig : ",BSt%eig(1:MIN(10,nband_tmp),ikpt,isppol)
455          !write(stdout,*) "eigens : ",Epren%eigens(1:MIN(10,nband_tmp),ikpt,isppol)
456          MSG_ERROR("Error in eigenvalues, check the _EP.nc file with respect to your input file !")
457        end if
458        if (ANY(ABS(Bst%occ(1:MIN(10,nband_tmp),ikpt,isppol) - Epren%occs(1:MIN(10,nband_tmp),ik_eph,isppol)) > tol3)) then
459          MSG_ERROR("Error in occupations, check the _EP.nc file with respect to your input file !")
460        end if
461      end if
462 
463      ! Upgrade energies
464      Bst%eig(1:nband_tmp,ikpt,isppol) = BSt%eig(1:nband_tmp,ikpt,isppol) + Epren%renorms(1,1:nband_tmp,ik_eph,isppol,itemp)
465 
466      if (do_lifetime) then
467        Bst%lifetime(1:nband_tmp,ikpt,isppol) = Epren%lifetimes(1,1:nband_tmp,ik_eph,isppol,itemp)
468      end if
469    end do
470  end do
471 
472  ABI_FREE(bs2eph)
473 
474 end subroutine renorm_bst