1! 2! Copyright (C) 1996-2016 The SIESTA group 3! This file is distributed under the terms of the 4! GNU General Public License: see COPYING in the top directory 5! or http://www.gnu.org/copyleft/gpl.txt. 6! See Docs/Contributors.txt for a list of contributors. 7! 8! This code segment has been fully created by: 9! Nick Papior Andersen, 2012, nickpapior@gmail.com 10! 11 12module m_ts_contour 13! 14! Routines that are used to setup the contour for integration of the GFs 15! 16! Use the type associated with the contour 17! Maybe they should be collected to this module. 18! However, I like this partition. 19 use m_ts_electype 20 use m_ts_chem_pot 21 22 use m_ts_cctype 23 use m_ts_io_contour 24 25 use m_ts_contour_eq, only : CONTOUR_EQ 26 use m_ts_contour_neq, only : CONTOUR_NEQ 27 28 use precision, only : dp 29 30 implicit none 31 32 private 33 34 public :: read_contour_options 35 public :: print_contour_options 36 public :: print_contour_block 37 public :: io_contour 38 public :: sort_contour 39 40contains 41 42 subroutine io_contour(IsVolt, mus, slabel,suffix) 43 use m_ts_contour_eq, only : io_contour_Eq 44 use m_ts_contour_neq, only : io_contour_nEq 45 logical, intent(in) :: IsVolt 46 type(ts_mu), intent(in) :: mus(:) 47 character(len=*), intent(in) :: slabel 48 character(len=*), intent(in), optional :: suffix 49 50 call io_contour_Eq(mus,slabel,suffix) 51 52 if ( IsVolt ) then 53 call io_contour_nEq(slabel,suffix) 54 end if 55 56 end subroutine io_contour 57 58 subroutine read_contour_options(N_Elec, Elecs, N_mu, mus, kT, IsVolt, Volt) 59 use m_ts_contour_eq, only : read_contour_eq_options 60 use m_ts_contour_neq, only : read_contour_neq_options 61 62 integer, intent(in) :: N_Elec 63 type(Elec), intent(inout) :: Elecs(N_Elec) 64 integer, intent(in) :: N_mu 65 type(ts_mu), intent(inout) :: mus(N_mu) 66 ! SIESTA electronic temperature 67 real(dp), intent(in) :: kT 68 logical, intent(in) :: IsVolt 69 real(dp), intent(in) :: Volt 70 71 call read_contour_eq_options(N_mu,mus,Volt) 72 73 if ( IsVolt ) then 74 call read_contour_neq_options(N_Elec,Elecs,N_mu,mus,kT,Volt) 75 end if 76 77 end subroutine read_contour_options 78 79 subroutine print_contour_options(prefix,IsVolt) 80 use m_ts_contour_eq, only : print_contour_eq_options 81 use m_ts_contour_neq, only : print_contour_neq_options 82 83 character(len=*), intent(in) :: prefix 84 logical, intent(in) :: IsVolt 85 86 call print_contour_eq_options(prefix) 87 88 if ( IsVolt ) then 89 call print_contour_neq_options(prefix) 90 end if 91 92 end subroutine print_contour_options 93 94 subroutine print_contour_block(prefix, IsVolt) 95 use m_ts_contour_eq, only : print_contour_eq_block 96 use m_ts_contour_neq, only : print_contour_neq_block 97 98 character(len=*), intent(in) :: prefix 99 logical, intent(in) :: IsVolt 100 101 call print_contour_eq_block(prefix) 102 103 if ( IsVolt ) then 104 call print_contour_neq_block(prefix) 105 end if 106 107 end subroutine print_contour_block 108 109 ! In order to retain the best numerical accuracy we introduce 110 ! an ordering of the energy contour by weights. 111 ! It is the only reasonable thing to access as the functional is 112 ! non-deterministic. 113 ! An example of the importance of this sorting: 114 ! Take +200 voltage energy points. 115 ! The weights will be something like this: 116 ! 1e-5,...,1e-4,...,1e-3,...,1e-4,...,1e-5 117 ! Which means that the summation up till the half works great. 118 ! But when we reach the last weight we could be in the situation where: 119 ! 1._dp + 1e-12_dp which will limit the accuracy obtained for that energy point. 120 121 ! In order to circumvent this we simply sort by weights. 122 subroutine sort_contour(NC,ce,cw) 123 integer, intent(in) :: NC 124 complex(dp), intent(inout) :: ce(NC), cw(NC) 125 126 ! Local variables 127 real(dp) :: tmp 128 complex(dp) :: ctmp 129 integer :: i,j 130 ! As we will only do this once we dont need a fancy sorting 131 ! algorithm... 132 do i = 1 , NC - 1 133 tmp = real(cw(i)*conjg(cw(i)),dp) 134 do j = i+1, NC 135 if ( tmp > real(cw(j)*conjg(cw(j))) ) then 136 ctmp = ce(i) 137 ce(i) = ce(j) 138 ce(j) = ctmp 139 ctmp = cw(i) 140 cw(i) = cw(j) 141 cw(j) = ctmp 142 end if 143 end do 144 end do 145 146 end subroutine sort_contour 147 148end module m_ts_contour 149