TABLE OF CONTENTS
ABINIT/m_gwls_GWanalyticPart [ Modules ]
NAME
m_gwls_GWanalyticPart
FUNCTION
.
COPYRIGHT
Copyright (C) 2009-2024 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 .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 23 module m_gwls_GWanalyticPart 24 !---------------------------------------------------------------------------------------------------- 25 ! This module contains routines to compute the contribution to the GW correlation energy coming 26 ! from the analytic integral of the trial frequency function. 27 !---------------------------------------------------------------------------------------------------- 28 ! local modules 29 use m_gwls_utility 30 use m_gwls_wf 31 use m_gwls_hamiltonian 32 use m_gwls_lineqsolver 33 use m_gwls_GWlanczos 34 use m_gwls_GenerateEpsilon 35 use m_gwls_LanczosBasis 36 use m_gwls_model_polarisability 37 use m_gwls_polarisability 38 !use m_gwls_ComplementSpacePolarizability 39 use m_gwls_GWlanczos 40 use m_gwls_TimingLog 41 42 ! abinit modules 43 use defs_basis 44 use m_abicore 45 46 implicit none 47 save 48 private
m_hamiltonian/get_projection_band_indices [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
get_projection_band_indices
FUNCTION
.
INPUTS
OUTPUT
SOURCE
81 subroutine get_projection_band_indices(omega,band_index_below, band_index_above) 82 !---------------------------------------------------------------------------------------------------- 83 ! This subroutine computes the band indices necessary for properly projecting the Sternheimer equations 84 ! where 85 ! Pe : projection on states such that epsilon_n < omega 86 ! Qe : projection on states such that epsilon_n > omega 87 !---------------------------------------------------------------------------------------------------- 88 implicit none 89 90 real(dp),intent(in) :: omega 91 integer, intent(out) :: band_index_below, band_index_above 92 93 ! ************************************************************************* 94 95 ! First, find the indices for the projections 96 97 band_index_below = 1 98 ! it is assumed that the eigenvalues are sorted. As soon as the current eigenvalue 99 ! is equal or larger than epsilon_e, exit! 100 do 101 if ( omega - eig(band_index_below) <= 1.0e-8 ) exit 102 band_index_below = band_index_below + 1 103 104 if (band_index_below > size(eig)) then 105 106 write(std_out,*) '************************************************************' 107 write(std_out,*) '*** ERROR IN ROUTINE get_projection_band_indices ***' 108 write(std_out,*) '*** ***' 109 write(std_out,*) '*** The index of the DFT eigenvalue larger than ***' 110 write(std_out,*) '*** the target frequency is larger than the number ***' 111 write(std_out,*) '*** of explicitly calculated DFT eigenvalues. ***' 112 write(std_out,*) '*** The computation cannot go on to produce a ***' 113 write(std_out,*) '*** meaningful result: review your input! ***' 114 write(std_out,*) '*** ***' 115 write(std_out,*) '*** program stops. ***' 116 write(std_out,*) '************************************************************' 117 stop 118 end if 119 end do 120 ! We have overshooted in order to exit the do loop , so decrement by 1 121 band_index_below = band_index_below - 1 122 123 ! band_index_below is now the index of the highest eigenvalue BELOW omega 124 ! NOTE THAT band_index_below = 0 if omega is smaller than all eigenvalues! 125 ! A test should be performed on the indices obtained from this routine 126 127 band_index_above = band_index_below+1 128 do 129 if ( eig(band_index_above)-omega > 1.0e-8 ) exit 130 band_index_above = band_index_above + 1 131 end do 132 ! band_index_above is now the index of the lowest eigenvalue ABOVE eig(e) 133 band_index_above = band_index_above - 1 134 ! band_index_above is now the highest index with eigenvalue equal to eig(e). 135 ! This is the correct index to remove the states with energy <= eig(e) 136 137 end subroutine get_projection_band_indices