TABLE OF CONTENTS


ABINIT/m_bs_defs [ Modules ]

[ Top ] [ Modules ]

NAME

  m_bs_defs

FUNCTION

  This module defines basic structures used for Bethe-Salpeter calculations.

COPYRIGHT

 Copyright (C) 1992-2018 EXC and ABINIT group (L.Reining, V.Olevano, F.Sottile, S.Albrecht, G.Onida, 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 .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

PARENTS

CHILDREN

SOURCE

22 #if defined HAVE_CONFIG_H
23 #include "config.h"
24 #endif
25 
26 #include "abi_common.h"
27 
28 MODULE m_bs_defs
29 
30  use defs_basis 
31  use m_abicore
32  use m_errors 
33 
34  implicit none
35 
36  private 

m_bs_defs/bs_parameters_free [ Functions ]

[ Top ] [ m_bs_defs ] [ Functions ]

NAME

  bs_parameters_free

FUNCTION

  Free all memory allocated in a structure of type excparam

SIDE EFFECTS

  Bsp<excparam>=All associated pointers are deallocated.

PARENTS

      bethe_salpeter

CHILDREN

      wrtout

SOURCE

297 subroutine bs_parameters_free(BSp)
298 
299 
300 !This section has been created automatically by the script Abilint (TD).
301 !Do not modify the following lines by hand.
302 #undef ABI_FUNC
303 #define ABI_FUNC 'bs_parameters_free'
304 !End of the abilint section
305 
306  implicit none
307 
308 !Arguments ------------------------------------
309 !scalars
310  type(excparam),intent(inout) :: BSp
311 
312 !************************************************************************
313  
314  !@excparam
315  if (allocated(BSp%q)) then
316    ABI_FREE(BSp%q)
317  end if
318 
319  if (allocated(Bsp%nreh)) then
320    ABI_FREE(Bsp%nreh)
321  end if
322  if (allocated(Bsp%vcks2t)) then
323    ABI_FREE(Bsp%vcks2t)
324  end if
325  if (allocated(Bsp%omega)) then
326    ABI_FREE(Bsp%omega)
327  end if
328 
329  if (allocated(Bsp%lomo_spin)) then 
330    ABI_FREE(Bsp%lomo_spin)
331  end if
332  if (allocated(Bsp%homo_spin)) then 
333    ABI_FREE(Bsp%homo_spin)
334  end if
335  if (allocated(Bsp%lumo_spin)) then 
336    ABI_FREE(Bsp%lumo_spin)
337  end if
338  if (allocated(Bsp%humo_spin)) then 
339    ABI_FREE(Bsp%humo_spin)
340  end if
341  if (allocated(Bsp%nbndv_spin)) then 
342    ABI_FREE(Bsp%nbndv_spin)
343  end if
344  if (allocated(Bsp%nbndc_spin)) then 
345    ABI_FREE(Bsp%nbndc_spin)
346  end if
347 
348  if (allocated(Bsp%Trans)) then
349    ABI_DT_FREE(Bsp%Trans)
350  end if
351 
352  if (allocated(Bsp%nreh_interp)) then
353    ABI_FREE(Bsp%nreh_interp)
354  end if
355 
356  if (allocated(Bsp%vcks2t_interp)) then
357    ABI_FREE(Bsp%vcks2t_interp)
358  end if
359 
360  if (allocated(Bsp%Trans_interp)) then
361    ABI_DT_FREE(Bsp%Trans_interp)
362  end if
363 
364 end subroutine bs_parameters_free

m_bs_defs/bsp_calctype2str [ Functions ]

[ Top ] [ m_bs_defs ] [ Functions ]

NAME

  bsp_calctype2str

FUNCTION

  Returns a string with the calculation type.

PARENTS

      m_bs_defs,m_exc_spectra,setup_bse

CHILDREN

      wrtout

SOURCE

572 subroutine bsp_calctype2str(BSp,str)
573 
574 
575 !This section has been created automatically by the script Abilint (TD).
576 !Do not modify the following lines by hand.
577 #undef ABI_FUNC
578 #define ABI_FUNC 'bsp_calctype2str'
579 !End of the abilint section
580 
581  implicit none
582 
583 !Arguments ------------------------------------
584 !scalars
585  character(len=500),intent(out) :: str
586  type(excparam),intent(in) :: BSp
587 
588 !************************************************************************
589 
590  SELECT CASE (Bsp%calc_type)
591  CASE (BSE_HTYPE_RPA_KS)
592    str = " RPA L0 with KS energies and KS wavefunctions"
593  CASE (BSE_HTYPE_RPA_QPENE)
594    str = " RPA L0 with QP energies and KS wavefunctions"
595  CASE (BSE_HTYPE_RPA_QP)
596    str = " RPA L0 with QP energies and QP wavefunctions"
597  CASE DEFAULT
598    str = " Unknown"
599  END SELECT
600 
601 end subroutine bsp_calctype2str

m_bs_defs/excfiles [ Types ]

[ Top ] [ m_bs_defs ] [ Types ]

NAME

  excfiles

FUNCTION

  The excfiles derived data type contains file names and unit numbers used to store 
  temporary or final results of the Bethe-Salpeter calculation. 

SOURCE

242 type,public :: excfiles
243 
244   character(len=fnlen) :: in_hreso = BSE_NOFILE   
245   ! Name of the input file with the resonant part of the Hamiltonian (Hermitian).
246 
247   character(len=fnlen) :: out_hreso = BSE_NOFILE
248   ! Name of the output file with the resonant part of the Hamiltonian (Hermitian).
249 
250   character(len=fnlen) :: in_hcoup = BSE_NOFILE   
251   ! Name of the input file with the coupling part of the Hamiltonian (Symmetric).
252 
253   character(len=fnlen) :: out_hcoup = BSE_NOFILE
254   ! Name of the output file with the coupling part of the Hamiltonian (Symmetric).
255 
256   character(len=fnlen) :: in_eig = BSE_NOFILE
257   ! Name of the input file with the eigenvalues and the eigenvectors of the Hamiltonian.
258 
259   character(len=fnlen) :: out_eig = BSE_NOFILE
260   ! Name of the output file with the eigenvalues and the eigenvectors of the Hamiltonian.
261 
262   character(len=fnlen) :: in_haydock_basename = BSE_NOFILE
263   ! Name of the input file used to restart Haydock algorithm. 
264 
265   character(len=fnlen) :: out_basename = BSE_NOFILE
266   ! Prefix to be used for other output files.
267 
268 end type excfiles
269 
270 public :: print_bs_files    ! Printout of the excfiles data type.

m_bs_defs/excparam [ Types ]

[ Top ] [ m_bs_defs ] [ Types ]

NAME

  excparam

FUNCTION

  The excparam derived data type contains the parameters controlling the BS calculation. 

SOURCE

111 type,public :: excparam
112 
113 !scalars
114   integer :: algorithm         ! Algorithm used for computing the BS dielectric function.
115   integer :: calc_type         ! Calculation type (see Dtset%bs_calc_type).
116   integer :: hayd_term         ! Option for the terminator used in the Haydock solver.
117   integer :: use_coupling      ! Include off-diagonal block coupling resonant and anti-resonant transitions.
118   integer :: exchange_term     ! Include the exchange term in the BS Hamiltonian.
119   integer :: inclvkb           ! Option for the inclusion of the commutator [Vnl, r] for the optical limit.
120   integer :: mdlf_type         ! Model dielectric function type. 
121   integer :: nline             ! Number of line minimization used for CG minimization.
122   integer :: nbdbuf            ! Number of states in the buffer that will be excluded from the convergence check (CG only)
123   integer :: nstates           ! Number of states that will be considered in the CG minimization.
124   integer :: npweps            ! No. of G in the Screening.
125   integer :: npwwfn            ! No. of G for wave functions.
126   !$integer :: npwx            ! No. of G for the exchange part.
127   integer :: npwvec            ! MAX between npwwfn and npweps
128   integer :: nbnds             ! Total number of bands considered.
129 
130   !integer :: nbndv             ! No. of valence states treated (homo-lomo+1)
131   !integer :: nbndc             ! No. of conduction states (humo-lumo+1)
132   !integer :: lomo,homo         ! Lowest and highest occupied orbital considered.
133   !integer :: lumo,humo         ! Lowest and highest unoccupied orbital considered.
134 
135 ! new for spin
136 ! DO I need these?
137   integer :: lomo_min,homo_max ! Lowest and highest occupied orbital considered.
138   integer :: lumo_min,humo_max ! Lowest and highest unoccupied orbital considered.
139   integer :: maxnbndv, maxnbndc
140 
141   integer,allocatable :: lomo_spin(:)     ! Lowest occupied orbital considered for the different spins.
142   integer,allocatable :: homo_spin(:)     ! Highest occupied orbital considered for the different spins.
143   integer,allocatable :: lumo_spin(:)     ! Lowest unoccupied orbital considered for the different spins
144   integer,allocatable :: humo_spin(:)     ! Highest unoccupied orbital considered for the different spins
145   integer,allocatable :: nbndv_spin(:)    ! No. of valence states treated (homo-lomo+1)
146   integer,allocatable :: nbndc_spin(:)    ! No. of conduction states (humo-lumo+1)
147 ! end new
148 
149   integer :: niter             ! No. of iterations for (Haydock|CG).
150   integer :: nkibz, nkbz       ! No. of k-points in the IBZ and BZ (resp.)
151   integer :: nomega            ! No. of frequencies for epsilon.
152   integer :: nq                ! Number of "small" q for optical limit.
153   integer :: nsppol            ! Number of independent spin polarizations.
154   integer :: wtype             ! Option used for dealing with W (see BSE_WTYPE_) flags
155 
156   !Interp@BSE
157   integer :: interp_mode       ! Mode of interpolation
158   integer :: interp_method     ! Method of interpolation
159   integer :: nkibz_interp,nkbz_interp ! Number of points in the interpolated kmesh
160   integer :: nstates_interp    ! Number of states of interpolation
161   integer :: rl_nb             ! Index of the nb in Rohlfing and Louie technique
162 
163   real(dp) :: ecutwfn          ! Cutoff energy for wavefunctions.
164   real(dp) :: ecuteps          ! Cutoff energy for W.
165   real(dp) :: eps_inf          ! Electronic dielectric constant used for the model dielectric function.
166   real(dp) :: mbpt_sciss         ! Scissors energy (used if it absolute value is > tol6)
167   real(dp) :: omegai           ! First omega for epsilon.
168   real(dp) :: omegae           ! Last omega for epsilon (defaults to 10eV)
169   real(dp) :: domega           ! Step of the frequency mesh.
170   real(dp) :: broad            ! Lorentzian Broadening.
171   real(dp) :: ircut            ! Infrared cutoff for transitions
172   real(dp) :: uvcut            ! Ultraviolet cutoff for transitions.
173   real(dp) :: haydock_tol(2)   ! Tolerance for stopping the Haydock algorithm.
174   real(dp) :: cg_tolwfr        ! Tolerance for stopping the CG algorithm 
175 
176   !Interp@BSE
177   real(dp) :: interp_m3_width  ! Width of the interpolated M3 method along the diagonal
178 
179   logical :: use_diagonal_Wgg  ! Use diagonal approximation for Wgg.
180   logical :: use_coulomb_term  ! Include W term in the BS Hamiltonian.
181   logical :: have_complex_ene  ! .TRUE. if energies have a non-zero imaginary part.
182 
183   !Interp@BSE
184   logical :: use_interp        ! .TRUE. if we use interpolation technique
185   logical :: prep_interp       ! .TRUE. if we prepare interpolation with ABC
186   logical :: sum_overlaps      ! .TRUE. if making the sum of the overlaps to 1
187   logical :: prt_ncham           ! .TRUE. if we dump the hamiltonian in NetCDF
188 
189   logical :: do_ep_renorm      ! .TRUE. for electron-phonon renormalization of the spectrum
190   logical :: do_lifetime       ! .TRUE. if using elphon lifetime (not yet implemented)
191 
192 !arrays
193   integer :: mg0(3)            ! For each reduced direction gives the max G0 component
194                                ! to account for umklapp processes
195 
196   integer,allocatable :: nreh(:)
197   ! nreh(nsppol)
198   ! Number of resonant electron-hole transitions for each spin.
199 
200   integer,allocatable :: vcks2t(:,:,:,:)  
201   ! vcks2t(v,c,ik_bz,spin) gives the transition index associated to (v,c,kbz,spin)
202 
203   !Interp@BSE
204   integer :: interp_kmult(3)   ! Factor to subdivide kmesh sampling
205 
206   integer,allocatable :: nreh_interp(:)
207   ! nreh_interp(nsppol)
208   ! Number of transitions for the interpolated mesh   
209 
210   integer,allocatable :: vcks2t_interp(:,:,:,:)
211   ! vcks2t(v,c,ik_bz_dense,spin) : Transition index for the dense kmesh
212 
213   real(dp),allocatable :: q(:,:)           ! Q-points for optical limit (reduced coordinates).
214 
215   complex(dpc),allocatable :: omega(:)
216   ! omega(nomega)
217   ! Frequency mesh for epsilon (including the complex imaginary shift)
218 
219   type(transition),allocatable :: Trans(:,:)
220   ! Trans(max_nreh,nsppol)
221 
222   type(transition),allocatable :: Trans_interp(:,:) 
223   ! Transitions for interpolated mesh
224 
225 end type excparam
226 
227  public :: bs_parameters_free
228  public :: print_bs_parameters
229  public :: bsp_calctype2str

m_bs_defs/init_transitions [ Functions ]

[ Top ] [ m_bs_defs ] [ Functions ]

NAME

  init_transitions

FUNCTION

  Main creation method for the transition structured datatype.

INPUTS

  lomo_spin(nsppol)
  humo_spin(nsppol)
  ir_cut,uv_cut
  nkbz=Number of k-points in the BZ.
  nbnds=Maximum number of bands.
  nkibz=Number of k-points in the IBZ.
  nsppol=Number of spins.
  nspinor=Number of spinor components.
  gw_energy
  occ=Occupation factors.
  ktab

OUTPUT

  max_tene=Maximum transition energy.
  nreh(nsppol)=Number of resonant transitions for each spin.

SIDE EFFECTS

  Trans(:,:)
    input:  allocatable array
    output: Trans(max_nreh,nsppol) stores the correspondence t -> (band,kbz,spin) and the transition energy.

PARENTS

      setup_bse,setup_bse_interp

CHILDREN

      wrtout

SOURCE

643 subroutine init_transitions(Trans,lomo_spin,humo_spin,ir_cut,uv_cut,nkbz,nbnds,nkibz,nsppol,nspinor,gw_energy,occ,ktab,&
644 &  minmax_tene,nreh) 
645 
646 
647 !This section has been created automatically by the script Abilint (TD).
648 !Do not modify the following lines by hand.
649 #undef ABI_FUNC
650 #define ABI_FUNC 'init_transitions'
651 !End of the abilint section
652 
653  implicit none
654 
655 !Arguments ------------------------------------
656 !scalars
657  integer,intent(in) :: nkbz,nbnds,nkibz,nsppol,nspinor
658  real(dp),intent(in) :: ir_cut,uv_cut
659  real(dp),intent(out) :: minmax_tene(2)
660  type(transition),allocatable,intent(out) :: Trans(:,:)
661 !arrays
662  integer,intent(in) :: lomo_spin(nsppol),humo_spin(nsppol)
663  integer,intent(in) :: ktab(nkbz) 
664  integer,intent(out) :: nreh(nsppol)
665  real(dp),intent(in) :: occ(nbnds,nkibz,nsppol)
666  complex(dpc),intent(in) :: gw_energy(nbnds,nkibz,nsppol)
667 
668 !Local variables ------------------------------
669 !scalars
670  integer :: spin,it,ik_bz,ik_ibz,iv,ic,max_occ,sweep,max_nreh,lomo,humo
671  real(dp) :: tene, delta_f,min_tene,max_tene
672  complex(dpc) :: cplx_enet
673  logical :: add_transition
674 
675 !************************************************************************
676       
677  ! Find transitions                                                  
678  max_occ=2/(nsppol*nspinor)
679  nreh=0
680  min_tene = -one; max_tene = zero
681  !
682  ! sweep=1 calculats the number of resonants transitions taking into 
683  !         account a possible energy cutoff.
684  ! sweep=2 initializes the tables describing the e-h transition.
685  !
686  do sweep=1,2
687    !
688    if (sweep==2) then
689      ! Allocate Trans structure.
690      max_nreh = MAXVAL(nreh)
691      ABI_DT_MALLOC(Trans, (max_nreh,nsppol))
692    end if
693    !
694    do spin=1,nsppol
695      it=0
696      lomo = lomo_spin(spin)
697      humo = humo_spin(spin)
698      do ik_bz=1,nkbz 
699        ik_ibz=ktab(ik_bz) 
700        !
701        do iv=lomo,humo
702          do ic=lomo,humo
703            delta_f   = ( occ(ic,ik_ibz,spin)-occ(iv,ik_ibz,spin) ) / max_occ
704            cplx_enet = gw_energy(ic,ik_ibz,spin)-gw_energy(iv,ik_ibz,spin)
705            tene = DBLE(cplx_enet)
706 
707            add_transition =                      &
708 &             (tene > tol12) .and.               &  ! Resonant transition.
709 &             ( ABS(delta_f) > tol12) .and.      &  ! c-v transition.
710 &             (tene < uv_cut .and. tene > ir_cut)   ! Energy cutoff.
711 
712            if (add_transition) then
713              it = it + 1 
714              max_tene = MAX(max_tene, tene)
715              min_tene = MAX(min_tene, tene)
716            end if
717            if (add_transition.and.sweep==2) then 
718              Trans(it,spin)%k  = ik_bz 
719              Trans(it,spin)%v  = iv 
720              Trans(it,spin)%c  = ic 
721              Trans(it,spin)%en = cplx_enet
722            end if
723 
724          end do 
725        end do 
726      end do ! ik_bz
727      ! Save number of transitions for this spin.
728      if (sweep==1) nreh(spin) = it 
729    end do ! spin
730    !
731  end do ! sweep
732 
733  minmax_tene = [min_tene, max_tene]
734 
735 end subroutine init_transitions                                           

m_bs_defs/print_bs_files [ Functions ]

[ Top ] [ m_bs_defs ] [ Functions ]

NAME

  print_bs_files

FUNCTION

  Printout of the content of the excfiles structure.

INPUTS

  BS_files<excfiles>=An object of type excfile storing the filenames used in the Bethe-Salpeter code.
  [unit]=Unit number for output
  [prtvol]=Verbosity level
  [mode_paral]=Either "COLL" or "PERS"
  [header]=String to be printed as header for additional info.

OUTPUT

  Only printing.

PARENTS

      setup_bse

CHILDREN

      wrtout

SOURCE

868 subroutine print_bs_files(BS_files,header,unit,mode_paral,prtvol)                                   
869 
870 
871 !This section has been created automatically by the script Abilint (TD).
872 !Do not modify the following lines by hand.
873 #undef ABI_FUNC
874 #define ABI_FUNC 'print_bs_files'
875 !End of the abilint section
876 
877  implicit none 
878 
879 !Arguments ------------------------------------
880 !scalars
881  type(excfiles),intent(in) :: BS_files
882  integer,optional,intent(in) :: unit,prtvol
883  character(len=4),optional,intent(in) :: mode_paral 
884  character(len=*),optional,intent(in) :: header
885 !arrays
886 
887 !Local variables ------------------------------
888 !scalars
889  integer :: my_unt,my_prtvol
890  character(len=4) :: my_mode
891  character(len=500) :: msg      
892 ! ********************************************************************* 
893 
894  !@excfiles
895  my_unt   =std_out; if (PRESENT(unit      )) my_unt   =unit
896  my_prtvol=0      ; if (PRESENT(prtvol    )) my_prtvol=prtvol 
897  my_mode  ='COLL' ; if (PRESENT(mode_paral)) my_mode  =mode_paral
898                                                                     
899  msg=' ==== Files used for the Bethe-Salpeter calculation  ==== '
900  if (PRESENT(header)) msg=' ==== '//TRIM(ADJUSTL(header))//' ==== '
901  call wrtout(my_unt,msg,my_mode)
902 
903  if (BS_files%in_hreso /= BSE_NOFILE) then
904    call wrtout(my_unt," Resonant block will be read from: "//trim(BS_files%in_hreso),my_mode)
905  end if
906 
907  if (BS_files%in_hcoup /= BSE_NOFILE) then
908    call wrtout(my_unt," Coupling block will be read from: "//trim(BS_files%in_hcoup),my_mode)
909  end if
910 
911  if (BS_files%in_eig /= BSE_NOFILE) then
912    call wrtout(my_unt," BS eigenstates will be read from: "//trim(BS_files%in_eig),my_mode)
913  end if
914 
915  if (BS_files%in_haydock_basename /= BSE_NOFILE) then
916    call wrtout(my_unt," Haydock restart files have basename: "//trim(BS_files%in_haydock_basename),my_mode)
917  end if
918 
919 end subroutine print_bs_files

m_bs_defs/print_bs_parameters [ Functions ]

[ Top ] [ m_bs_defs ] [ Functions ]

NAME

  print_bs_parameters

FUNCTION

  Printout of the parameters used for the BS calculation.

INPUTS

  p<excparam>=Datatype storing the parameters of the Bethe-Salpeter calculation.

OUTPUT

  Only printing.

PARENTS

      setup_bse

CHILDREN

      wrtout

SOURCE

390 subroutine print_bs_parameters(BSp,header,unit,mode_paral,prtvol) 
391 
392 
393 !This section has been created automatically by the script Abilint (TD).
394 !Do not modify the following lines by hand.
395 #undef ABI_FUNC
396 #define ABI_FUNC 'print_bs_parameters'
397 !End of the abilint section
398 
399  implicit none
400 
401 !Arguments ------------------------------------
402 !scalars
403  integer,optional,intent(in) :: unit,prtvol
404  character(len=4),optional,intent(in) :: mode_paral 
405  character(len=*),optional,intent(in) :: header
406  type(excparam),intent(inout) :: BSp
407 
408 !Local variables ------------------------------
409 !scalars
410  integer :: my_unt,my_prtvol,iq,ii,spin
411  character(len=4) :: my_mode
412  character(len=500) :: msg      
413 
414 ! ********************************************************************* 
415 
416  my_unt   =std_out; if (PRESENT(unit      )) my_unt   =unit
417  my_prtvol=0      ; if (PRESENT(prtvol    )) my_prtvol=prtvol 
418  my_mode  ='COLL' ; if (PRESENT(mode_paral)) my_mode  =mode_paral
419 
420  msg=' ==== Parameters of the Bethe-Salpeter run ==== '
421  if (PRESENT(header)) msg=' ==== '//TRIM(ADJUSTL(header))//' ==== '
422  call wrtout(my_unt,msg,my_mode)
423 
424  select case (Bsp%algorithm)
425  case (BSE_ALGO_NONE)
426    msg = " Algorithm: Build Hamiltonian but skip the calculation of the spectrum."
427  case (BSE_ALGO_DDIAGO)
428    msg = " Algorithm: Direct diagonalization."
429  case (BSE_ALGO_HAYDOCK)
430    msg = " Algorithm: Haydock technique."
431  case (BSE_ALGO_CG)
432    msg = " Algorithm: Conjugate gradient."
433  case default
434    msg = " Algorithm: Unknown!."
435  end select
436  call wrtout(my_unt,msg,my_mode)
437 
438  write(msg,'(4(a,i0,a))')&
439 &  ' Dimension of the v, W matrices,  npweps  = ',BSp%npweps,ch10,&
440 &  ' Cutoff for the wavefunctions,    npwwfn  = ',BSp%npwwfn,ch10,&
441 &  ' Number of k-points in the IBZ,   nkibz   = ',BSp%nkibz,ch10,&
442 &  ' Highest empty band included,     nband   = ',BSp%nbnds,""
443  call wrtout(my_unt,msg,my_mode)
444 
445  do spin=1,Bsp%nsppol
446    msg = " === Spin UP ==="; if (spin == 2) msg = " === Spin DOWN ==="
447    call wrtout(my_unt,msg,my_mode)
448    write(msg,'(5(a,i0,a))')&
449 &    ' Number of resonant transitions          ',BSp%nreh(spin),ch10,&
450 &    ' Lowest occupied state                   ',BSp%lomo_spin(spin),ch10,&
451 &    ' Highest occupied state                  ',BSp%homo_spin(spin),ch10,&
452 &    ' Lowest unoccupied state                 ',BSp%lumo_spin(spin),ch10,&
453 &    ' Highest unoccupied state                ',BSp%nbnds,""
454 !&    ' Number of valence bands                 ',BSp%nbndv,ch10,&
455 !&    ' Number of conduction bands              ',BSp%nbndc,""
456    call wrtout(my_unt,msg,my_mode)
457  end do
458 
459  write(msg,'(3(a,f6.2,a),a,f6.2)')&
460 &  ' Minimum frequency [eV]           Emin    = ',BSp%omegai*Ha_eV,ch10,&
461 &  ' Maximum frequency [eV]           Emax    = ',BSp%omegae*Ha_eV,ch10,&
462 &  ' Frequency step [eV]              dE      = ',BSp%domega*Ha_eV,ch10,&
463 &  ' Lorentzian broadening [eV]       eta     = ',BSp%broad*Ha_eV
464  call wrtout(my_unt,msg,my_mode)
465  !
466  ! Calculation type
467  call bsp_calctype2str(Bsp,msg)
468  call wrtout(my_unt,msg,my_mode)
469 
470  if (ABS(Bsp%mbpt_sciss)>tol6) then
471    write(msg,'(a,f5.2)')&
472 &   " Scissors operator energy [eV] =         ",Bsp%mbpt_sciss*Ha_eV
473    call wrtout(my_unt,msg,my_mode)
474  end if
475 
476  msg=' Local fields effects (v term) excluded'
477  if (BSp%exchange_term>0) msg=' Local fields effects (v term) included'
478  call wrtout(my_unt,msg,my_mode)
479 
480  msg=' Excitonic effects (W term) excluded'
481  if (BSp%use_coulomb_term) msg=' Excitonic effects (W term) included'
482  call wrtout(my_unt,msg,my_mode)
483 
484  if (BSp%use_coulomb_term) then 
485    msg=" Full W_GG' included"
486    if (BSp%use_diagonal_Wgg) msg=' Only diagonal term W_GG included'
487    call wrtout(my_unt,msg,my_mode)
488    if (BSp%wtype==BSE_WTYPE_FROM_SCR) then
489      call wrtout(my_unt," W is read from an external SCR file",my_mode)
490    end if
491    if (BSp%wtype==BSE_WTYPE_FROM_MDL) then
492      call wrtout(my_unt," W is approximated with the model dielectric function",my_mode)
493    end if
494  end if
495 
496  msg=' Resonant-only calculation (Hermitian case)'
497  if (BSp%use_coupling>0) msg=' Resonant + Coupling calculation'
498  call wrtout(my_unt,msg,my_mode)
499 
500  if(Bsp%use_interp) then
501    call wrtout(my_unt,' Interpolation technique used',my_mode)   
502  end if
503    
504  if(Bsp%use_interp) then
505    select case (Bsp%interp_mode)
506    case (1) 
507      msg = ' Interpolation using WFK on the dense mesh'
508    case (2)
509      msg = ' Interpolation using WFK on the dense mesh + ABC divergence'
510    case (3)
511      msg = ' Interpolation using WFK on the dense mesh + ABC divergence along diagonal'
512    case (4)
513      msg = ' Interpolation using WFK on the dense mesh'
514    case default
515      msg = ' Unknown interpolation technique'
516    end select
517    call wrtout(my_unt,msg,my_mode)
518 
519    if(BSp%prep_interp) then
520      call wrtout(my_unt,' Prepare interpolation technique with ABC',my_mode)
521    end if
522 
523    if(BSp%interp_method == BSE_INTERP_YG) then
524      write(msg,'(a)') " Use Y. Gillet interpolation with 8 neighbours"
525      call wrtout(my_unt, msg, my_mode)
526    else if(BSP%interp_method == BSE_INTERP_RL2) then
527      write(msg,'(a)') " Use only 2 neighbours to interpolate linearly (Debug mode)"
528      call wrtout(my_unt, msg, my_mode)
529    else if(BSp%interp_method == BSE_INTERP_RL) then
530      write(msg,'(a,i0)') " Use Rohlfing and Louie with nb = ",BSp%rl_nb
531      call wrtout(my_unt, msg, my_mode)
532    end if
533 
534    if(BSp%sum_overlaps) then
535      call wrtout(my_unt, " Summing overlaps to 1 in the interpolation",my_mode)
536    end if
537  end if
538 
539  write(msg,'(a)')ch10
540  call wrtout(my_unt,msg,my_mode)
541 
542  call wrtout(my_unt,' Calculating epsilon_Macro(q-->0,w), along the following directions:',my_mode)
543  do iq=1,BSp%nq
544    write(msg,'(a,3f10.6,2a)')' q = (',(BSp%q(ii,iq),ii=1,3),  ') [r.l.u.]'
545    call wrtout(my_unt,msg,my_mode)
546  end do
547 
548  !TODO 
549  !Add file sizes and size of the buffer used for the matrix.
550 
551 end subroutine print_bs_parameters

m_bs_defs/repr_1trans [ Functions ]

[ Top ] [ m_bs_defs ] [ Functions ]

NAME

  repr_1trans

FUNCTION

  Returns a string with info on the (k,v,c) transition.

INPUTS

  Trans<transition>=structure datatype containing indececes and info on the optical transition.
  [prtvol]=Verbosity level. Defaults to 0.

OUTPUT

  str(len=500)=The string representing the transition.

PARENTS

SOURCE

758 pure function repr_1trans(Trans,prtvol) result(str)
759 
760 
761 !This section has been created automatically by the script Abilint (TD).
762 !Do not modify the following lines by hand.
763 #undef ABI_FUNC
764 #define ABI_FUNC 'repr_1trans'
765 !End of the abilint section
766 
767  implicit none 
768 
769 !Arguments ------------------------------------
770 !scalars
771  integer,optional,intent(in) :: prtvol
772  character(len=500) :: str
773  type(transition),intent(in) :: Trans 
774 
775 !Local variables ------------------------------
776 !scalars
777  integer :: my_prtvol
778 
779 !************************************************************************
780 
781  my_prtvol=0; if (PRESENT(prtvol)) my_prtvol=prtvol
782 
783  if (my_prtvol==0) then
784    write(str,'(3(a,i3))')" k= ",Trans%k," v= ",Trans%v," c= ",Trans%c
785  else 
786    write(str,'(3(a,i3),a,2f6.2)')" k= ",Trans%k," v= ",Trans%v," c= ",Trans%c," ene= ",Trans%en*Ha_eV
787  end if
788 
789 end function repr_1trans

m_bs_defs/repr_2trans [ Functions ]

[ Top ] [ m_bs_defs ] [ Functions ]

NAME

  repr_2trans

FUNCTION

  Returns a string with info on two transitions

INPUTS

  Trans1, Trans2<transition>=structure datatypes containing indececes and info on the optical transitions
  [prtvol]=Verbosity level. Defaults to 0.

OUTPUT

  string(len=500)=The string representing the transition.

PARENTS

SOURCE

812 pure function repr_2trans(Trans1,Trans2,prtvol) result(string)
813 
814 
815 !This section has been created automatically by the script Abilint (TD).
816 !Do not modify the following lines by hand.
817 #undef ABI_FUNC
818 #define ABI_FUNC 'repr_2trans'
819 !End of the abilint section
820 
821  implicit none 
822 
823 !Arguments ------------------------------------
824 !scalars
825  integer,optional,intent(in) :: prtvol
826  character(len=500) :: string
827  type(transition),intent(in) :: Trans1,Trans2
828 
829 !Local variables ------------------------------
830 !scalars
831  integer :: my_prtvol
832 
833 !************************************************************************
834 
835  my_prtvol=0; if (PRESENT(prtvol)) my_prtvol=prtvol
836  string = repr_1trans(Trans1,my_prtvol)//" | "//trim(repr_1trans(Trans2,my_prtvol))
837 
838 end function repr_2trans

m_bs_defs/transition [ Types ]

[ Top ] [ m_bs_defs ] [ Types ]

NAME

  transition

FUNCTION

  The transition derived data type is used to store the correspondence 
  between the transition index and the set of quantum numbers (ik_bz,v,c) 
  The energy of the transition is stored as well.

SOURCE

84  type,public :: transition
85    integer :: k=0               ! Index of the k-point in the BZ
86    integer :: v=0               ! Valence band index.
87    integer :: c=0               ! Conduction band index.
88    complex(dpc) :: en=huge(one) ! Transition energy 
89  end type transition
90 
91  public :: init_transitions     ! Main creation method.
92  public :: repr_trans           ! Returns a string representing the transition or a couple of transitions.