TABLE OF CONTENTS


ABINIT/m_gwls_GWanalyticPart [ Modules ]

[ Top ] [ Modules ]

NAME

 m_gwls_GWanalyticPart

FUNCTION

  .

COPYRIGHT

 Copyright (C) 2009-2018 ABINIT group (JLJ, BR, MC)
 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

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 
28 module m_gwls_GWanalyticPart
29 !----------------------------------------------------------------------------------------------------
30 ! This module contains routines to compute the contribution to the GW correlation energy coming
31 ! from the analytic integral of the trial frequency function.
32 !----------------------------------------------------------------------------------------------------
33 ! local modules
34 use m_gwls_utility
35 use m_gwls_wf
36 use m_gwls_hamiltonian
37 use m_gwls_lineqsolver
38 use m_gwls_GWlanczos
39 use m_gwls_GenerateEpsilon
40 use m_gwls_LanczosBasis
41 use m_gwls_model_polarisability
42 use m_gwls_polarisability
43 !use m_gwls_ComplementSpacePolarizability
44 use m_gwls_GWlanczos
45 use m_gwls_TimingLog
46 
47 ! abinit modules
48 use defs_basis
49 use m_abicore
50 
51 implicit none
52 save
53 private

m_hamiltonian/get_projection_band_indices [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  get_projection_band_indices

FUNCTION

  .

INPUTS

OUTPUT

PARENTS

      gwls_Projected_AT

CHILDREN

SOURCE

 91 subroutine get_projection_band_indices(omega,band_index_below, band_index_above)
 92 !----------------------------------------------------------------------------------------------------
 93 ! This subroutine computes the band indices necessary for properly projecting the Sternheimer equations
 94 ! where
 95 !                Pe : projection on states such that epsilon_n < omega
 96 !                Qe : projection on states such that epsilon_n > omega
 97 !----------------------------------------------------------------------------------------------------
 98 
 99 !This section has been created automatically by the script Abilint (TD).
100 !Do not modify the following lines by hand.
101 #undef ABI_FUNC
102 #define ABI_FUNC 'get_projection_band_indices'
103 !End of the abilint section
104 
105 implicit none
106 
107 real(dp),intent(in)  :: omega 
108 integer, intent(out) :: band_index_below, band_index_above
109 
110 ! *************************************************************************
111 
112 ! First, find the indices for the projections
113 
114 band_index_below = 1
115 ! it is assumed that the eigenvalues are sorted. As soon as the current eigenvalue
116 ! is equal or larger than epsilon_e,  exit!
117 do
118 if ( omega - eig(band_index_below) <= 1.0e-8 ) exit
119 band_index_below = band_index_below + 1
120 
121 if (band_index_below > size(eig)) then
122 
123   write(std_out,*) '************************************************************'
124   write(std_out,*) '***    ERROR IN ROUTINE get_projection_band_indices      ***'
125   write(std_out,*) '***                                                      ***'
126   write(std_out,*) '***    The index of the DFT eigenvalue larger than       ***'
127   write(std_out,*) '***    the target frequency is larger than the number    ***'
128   write(std_out,*) '***    of explicitly calculated DFT eigenvalues.         ***'
129   write(std_out,*) '***    The computation cannot go on to produce a         ***'
130   write(std_out,*) '***    meaningful result: review your input!             ***'
131   write(std_out,*) '***                                                      ***'
132   write(std_out,*) '***               program stops.                         ***'
133   write(std_out,*) '************************************************************'
134   stop
135 end if
136 end do
137 ! We have overshooted in order to exit the do loop , so decrement by 1
138 band_index_below = band_index_below - 1
139 
140 ! band_index_below  is now the index of the highest eigenvalue BELOW omega
141 ! NOTE THAT band_index_below = 0 if omega is smaller than all eigenvalues!
142 ! A test should be performed on the indices obtained from this routine
143 
144 band_index_above = band_index_below+1
145 do
146 if ( eig(band_index_above)-omega > 1.0e-8 ) exit
147 band_index_above = band_index_above + 1
148 end do
149 ! band_index_above is now the index of the lowest eigenvalue ABOVE eig(e)
150 band_index_above = band_index_above - 1
151 ! band_index_above is now the highest index with eigenvalue equal to eig(e).
152 ! This is the correct index to remove the states with energy <= eig(e)
153 
154 end subroutine get_projection_band_indices