TABLE OF CONTENTS


ABINIT/gwls_sternheimer [ Functions ]

[ Top ] [ Functions ]

NAME

 gwls_sternheimer

FUNCTION

 .

INPUTS

 dtset <type(dataset_type)>=all input variables in this dataset

OUTPUT

PARENTS

      driver

CHILDREN

      cleanupvalencewavefunctions,close_timing_log
      compute_correlations_no_model_shift_lanczos
      compute_correlations_shift_lanczos
      compute_exchange_and_correlation_energies,destroy_h
      driver_generateprintdielectriceigenvalues,gstateimg
      preparevalencewavefunctions,setup_timing_log,timab

SOURCE

 65 subroutine gwls_sternheimer(acell_img,amu_img,codvsn,cpui,dtfil,dtset,etotal_img,fcart_img,&
 66 &                           fred_img,iexit,mixalch_img,mpi_enreg,nimage,npwtot,occ_img,&
 67 &                           pawang,pawrad,pawtab,psps,rprim_img,strten_img,vel_cell_img,vel_img,xred_img,&
 68 &                           filnam,filstat,idtset,jdtset,ndtset) ! optional arguments
 69 
 70 use m_gwls_utility,                   only : master_debug, files_status_new, files_status_old
 71 ! use m_gwls_wf,                        only : norm_k
 72 use m_gwls_hamiltonian,               only : destroy_H, exchange, g_to_r, eig
 73 use m_gwls_TimingLog
 74 use m_gwls_valenceWavefunctions
 75 use m_gwls_ComputeCorrelationEnergy
 76 use m_gwls_GenerateEpsilon
 77 
 78 use defs_basis
 79 use defs_datatypes
 80 use defs_abitypes
 81 use defs_wvltypes
 82 use m_pawang
 83 use m_pawrad
 84 use m_pawtab
 85 use m_abicore
 86 use m_errors
 87 
 88 use m_time,      only : timab
 89 use m_gstateimg, only : gstateimg
 90 
 91 !This section has been created automatically by the script Abilint (TD).
 92 !Do not modify the following lines by hand.
 93 #undef ABI_FUNC
 94 #define ABI_FUNC 'gwls_sternheimer'
 95 !End of the abilint section
 96 
 97 implicit none
 98 !Arguments ------------------------------------
 99 !scalars
100 integer,intent(in) :: nimage
101 integer,optional,intent(in) :: idtset,ndtset
102 integer,intent(inout) :: iexit
103 real(dp),intent(in) :: cpui
104 character(len=6),intent(in) :: codvsn
105 character(len=fnlen),optional,intent(in) :: filstat
106 type(MPI_type),intent(inout) :: mpi_enreg
107 type(datafiles_type),target,intent(inout) :: dtfil
108 type(dataset_type),intent(inout) :: dtset
109 type(pawang_type),intent(inout) :: pawang
110 type(pseudopotential_type),intent(inout) :: psps
111 !arrays
112 integer,optional,intent(in) :: jdtset(:)
113 integer,intent(out) :: npwtot(dtset%nkpt)
114 character(len=fnlen),optional,intent(in) :: filnam(:)
115 real(dp), intent(out) :: etotal_img(nimage),fcart_img(3,dtset%natom,nimage)
116 real(dp), intent(out) :: fred_img(3,dtset%natom,nimage),strten_img(6,nimage)
117 real(dp),intent(inout) :: acell_img(3,nimage),amu_img(dtset%ntypat,nimage),occ_img(dtset%mband*dtset%nkpt*dtset%nsppol,nimage)
118 real(dp),intent(inout) :: mixalch_img(dtset%npspalch,dtset%ntypalch,nimage)
119 real(dp),intent(inout) :: rprim_img(3,3,nimage),vel_cell_img(3,3,nimage),vel_img(3,dtset%natom,nimage)
120 real(dp),intent(inout) :: xred_img(3,dtset%natom,nimage)
121 type(pawrad_type),intent(inout) :: pawrad(psps%ntypat*psps%usepaw)
122 type(pawtab_type),intent(inout) :: pawtab(psps%ntypat*psps%usepaw)
123 !Local variables ------------------------------
124 
125 integer  :: e_index
126 real(dp) :: exchange_energy
127 real(dp) :: vxc_energy
128 real(dp) :: DFT_band_energy
129  type(wvl_data) :: wvl
130 
131 ! timing
132 real(dp) :: tsec(2)
133 integer :: GWLS_TIMAB, OPTION_TIMAB
134 !END DEBUG
135 ! *************************************************************************
136 
137 ! initialize to zero
138 
139 ! start new increment of the GWLS timer
140  GWLS_TIMAB   = 1501
141  OPTION_TIMAB = 1
142  call timab(GWLS_TIMAB,OPTION_TIMAB,tsec)
143 
144 
145  master_debug = .false.
146 !master_debug = .true.
147 
148 !Governs the status of all opened files formerly with the status 'new' harcoded.
149 !It is initialized in the m_gwls_utility module as 'new'; this default is overriden here if desired.
150  files_status_new = 'unknown'
151 !Governs the status of all opened files formerly with the status 'old' harcoded.
152 !It is initialized in the m_gwls_utility module as 'old'; this default is overriden here if desired.
153  files_status_old = 'unknown'
154 
155 ! Test the input to make sure it is consistent with what
156 ! the code can do
157 !call test_input(dtset2,psps2)
158 
159 !Initializing variables and build the Hamiltonian using build_H
160 
161  GWLS_TIMAB   = 1521
162  OPTION_TIMAB = 1
163  call timab(GWLS_TIMAB,OPTION_TIMAB,tsec)
164 
165 
166  call gstateimg(acell_img,amu_img,codvsn,cpui,dtfil,dtset,etotal_img,fcart_img,&
167 & fred_img,iexit,mixalch_img,mpi_enreg,nimage,npwtot,occ_img,&
168 & pawang,pawrad,pawtab,psps,rprim_img,strten_img,vel_cell_img,vel_img,wvl,xred_img,&
169 & filnam,filstat,idtset,jdtset,ndtset) ! optional arguments
170 
171  OPTION_TIMAB = 2
172  call timab(GWLS_TIMAB,OPTION_TIMAB,tsec)
173 
174 
175 ! tabulate the valence wavefunctions, to be used throughout the code!
176 ! TODO : put in build_H and destroy_H
177  GWLS_TIMAB   = 1522
178  OPTION_TIMAB = 1
179  call timab(GWLS_TIMAB,OPTION_TIMAB,tsec)
180 
181  call prepareValenceWavefunctions()
182 
183  OPTION_TIMAB = 2
184  call timab(GWLS_TIMAB,OPTION_TIMAB,tsec)
185 
186 
187 
188  e_index         = dtset%gwls_band_index
189  DFT_band_energy = eig(e_index)
190 
191  write(std_out,10) '                               '
192  write(std_out,10) ' GWLS RESULTS                  '
193  write(std_out,10) ' ----------------------------- '
194  write(std_out,12) ' For orbital |psi_e>, with e : ',e_index
195  write(std_out,14) ' DFT eigenenergy             : ',DFT_band_energy,' Ha = ',DFT_band_energy*Ha_eV,' eV'
196  flush(std_out)
197 
198  write(ab_out,10) '                               '
199  write(ab_out,10) ' GWLS RESULTS                  '
200  write(ab_out,10) ' ----------------------------- '
201  write(ab_out,12) ' For orbital |psi_e>, with e : ',e_index
202  write(ab_out,14) ' DFT eigenenergy             : ',DFT_band_energy,' Ha = ',DFT_band_energy*Ha_eV,' eV'
203 
204  if (dtset%gwls_exchange /= 0) then
205 
206    GWLS_TIMAB   = 1502
207    OPTION_TIMAB = 1
208    call timab(GWLS_TIMAB,OPTION_TIMAB,tsec)
209 
210    call compute_Exchange_and_Correlation_energies(e_index, exchange_energy, Vxc_energy)
211 
212    OPTION_TIMAB = 2
213    call timab(GWLS_TIMAB,OPTION_TIMAB,tsec)
214 
215 
216 
217 
218    write(std_out,14) ' <psi_e |   V_xc  | psi_e>   : ',vxc_energy     ,' Ha = ',vxc_energy     *Ha_eV,' eV'
219    write(std_out,14) ' <psi_e | Sigma_x | psi_e>   : ',exchange_energy,' Ha = ',exchange_energy*Ha_eV,' eV'
220    flush(std_out)
221 
222    write(ab_out,14) ' <psi_e |   V_xc  | psi_e>   : ',vxc_energy     ,' Ha = ',vxc_energy     *Ha_eV,' eV'
223    write(ab_out,14) ' <psi_e | Sigma_x | psi_e>   : ',exchange_energy,' Ha = ',exchange_energy*Ha_eV,' eV'
224  else
225    exchange_energy  = zero
226    Vxc_energy       = zero
227 
228  end if
229 
230  if (dtset%gwls_correlation == 3) then
231 
232    call setup_timing_log()
233 
234    GWLS_TIMAB   = 1503
235    OPTION_TIMAB = 1
236 
237    call timab(GWLS_TIMAB,OPTION_TIMAB,tsec)
238 
239    call compute_correlations_shift_lanczos(dtset, exchange_energy, Vxc_energy,master_debug)
240 
241    OPTION_TIMAB = 2
242    call timab(GWLS_TIMAB,OPTION_TIMAB,tsec)
243 
244 
245    call close_timing_log()
246 
247  else if (dtset%gwls_correlation == 4) then
248    call setup_timing_log()
249    call compute_correlations_no_model_shift_lanczos(dtset, exchange_energy, Vxc_energy,master_debug)
250    call close_timing_log()
251 
252  else if (dtset%gwls_correlation == 5) then
253    call setup_timing_log()
254    call Driver_GeneratePrintDielectricEigenvalues(dtset)
255    call close_timing_log()
256 
257 
258  end if
259 
260  call cleanupValenceWavefunctions()
261 
262  call destroy_H()
263 
264  GWLS_TIMAB   = 1501
265  OPTION_TIMAB = 2
266  call timab(GWLS_TIMAB,OPTION_TIMAB,tsec)
267 
268 
269  10 format(A)
270  12 format(A,I6)
271  14 format(A,ES24.16,A,F16.8,A)
272 
273 end subroutine gwls_sternheimer

ABINIT/m_gwls_sternheimer [ Modules ]

[ Top ] [ Modules ]

NAME

   m_gwls_sternheimer

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 module m_gwls_sternheimer
28 
29  implicit none
30 
31  private