TABLE OF CONTENTS
- ABINIT/gwls_LanczosBasis
- m_hamiltonian/cleanup_Lanczos_basis
- m_hamiltonian/modify_Lbasis_Coulomb
- m_hamiltonian/setup_Lanczos_basis
ABINIT/gwls_LanczosBasis [ Modules ]
NAME
gwls_LanczosBasis
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 gwls_LanczosBasis 29 !---------------------------------------------------------------------------------------------------- 30 ! This module contains the static Lanczos basis, which should be computed once and for all, 31 ! and then be made available to other modules. 32 !---------------------------------------------------------------------------------------------------- 33 ! local modules 34 use m_gwls_utility 35 use gwls_wf 36 use gwls_hamiltonian 37 use gwls_GenerateEpsilon 38 39 ! Abinit modules 40 use defs_basis 41 use m_profiling_abi 42 43 implicit none 44 save 45 private
m_hamiltonian/cleanup_Lanczos_basis [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
cleanup_Lanczos_basis
FUNCTION
.
INPUTS
OUTPUT
PARENTS
gwls_ComputeCorrelationEnergy
CHILDREN
g_to_r,gr_to_g,sqrt_vc_k,wf_block_distribute
SOURCE
147 subroutine cleanup_Lanczos_basis() 148 149 150 !This section has been created automatically by the script Abilint (TD). 151 !Do not modify the following lines by hand. 152 #undef ABI_FUNC 153 #define ABI_FUNC 'cleanup_Lanczos_basis' 154 !End of the abilint section 155 156 implicit none 157 158 ! ************************************************************************* 159 160 ! if(allocated()) ABI_DEALLOCATE() can cause a segfault if used in this form, without the THEN. 161 ! This is because ABI_ALLOCATE is expanded at compilation time in many statements; and the if() can only prevent the execution of 162 ! the first. 163 ! So, the deallocation takes place even if allocated() returns .false. without the THEN. 164 if (allocated(Lbasis_lanczos)) then 165 ABI_DEALLOCATE(Lbasis_lanczos) 166 end if 167 if (allocated(Lbasis_model_lanczos)) then 168 ABI_DEALLOCATE(Lbasis_model_lanczos) 169 end if 170 171 end subroutine cleanup_Lanczos_basis
m_hamiltonian/modify_Lbasis_Coulomb [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
modify_Lbasis_Coulomb
FUNCTION
.
INPUTS
OUTPUT
PARENTS
gwls_ComputeCorrelationEnergy
CHILDREN
g_to_r,gr_to_g,sqrt_vc_k,wf_block_distribute
SOURCE
193 subroutine modify_Lbasis_Coulomb(psie_k, lmax, lmax_model) 194 !---------------------------------------------------------------------------------------------------- 195 ! This subroutine computes, once and for all, the vectors (V^1/2.L)^*.psie 196 !---------------------------------------------------------------------------------------------------- 197 198 !This section has been created automatically by the script Abilint (TD). 199 !Do not modify the following lines by hand. 200 #undef ABI_FUNC 201 #define ABI_FUNC 'modify_Lbasis_Coulomb' 202 !End of the abilint section 203 204 real(dp),intent(in) :: psie_k(2,npw_k) 205 integer ,intent(in) :: lmax, lmax_model 206 207 208 real(dp), allocatable :: psik_wrk(:,:), psikb_wrk(:,:), psikg_wrk(:,:) 209 real(dp), allocatable :: psikg_e(:,:) 210 211 integer :: iblk_lanczos, nbdblock_lanczos 212 213 integer :: l, mb 214 215 ! ************************************************************************* 216 217 ABI_ALLOCATE(psik_wrk ,(2,npw_k)) 218 ABI_ALLOCATE(psikb_wrk ,(2,npw_kb)) 219 ABI_ALLOCATE(psikg_wrk ,(2,npw_g)) 220 ABI_ALLOCATE(psikg_e ,(2,npw_g)) 221 222 ! copy the valence state on every row of FFT processors 223 do mb = 1, blocksize 224 psikb_wrk(:,(mb-1)*npw_k+1:mb*npw_k) = psie_k(:,:) 225 end do ! mb 226 ! Transform to FFT representation 227 call wf_block_distribute(psikb_wrk, psikg_e,1) ! LA -> FFT 228 229 230 ! Number of blocks of lanczos vectors 231 nbdblock_lanczos = lmax/blocksize 232 if (modulo(lmax,blocksize) /= 0) nbdblock_lanczos = nbdblock_lanczos + 1 233 234 ! loop on all blocks of lanczos vectors 235 do iblk_lanczos = 1, nbdblock_lanczos 236 237 ! loop on all states within this block 238 do mb = 1, blocksize 239 ! Determine the index of the Lanczos vector 240 l = (iblk_lanczos-1)*blocksize + mb 241 242 if ( l <= lmax) then 243 psik_wrk(1,:) = dble (Lbasis_lanczos(:,l)) 244 psik_wrk(2,:) = dimag(Lbasis_lanczos(:,l)) 245 else 246 psik_wrk(:,:) = zero 247 end if 248 ! Apply coulomb potential 249 call sqrt_vc_k(psik_wrk) 250 251 ! Store in array of blocks of wavefunctions 252 psikb_wrk(:,(mb-1)*npw_k+1:mb*npw_k) = psik_wrk(:,:) 253 254 end do ! mb 255 256 ! Transform to FFT representation 257 call wf_block_distribute(psikb_wrk, psikg_wrk,1) ! LA -> FFT 258 259 ! fourier transform, conjugate 260 call g_to_r(psir1,psikg_wrk) 261 psir1(2,:,:,:) = -psir1(2,:,:,:) 262 263 ! compute the product with the state "e" 264 call gr_to_g(psikg_wrk, psir1, psikg_e) 265 266 ! return to LA configuration 267 268 ! Transform to LA representation 269 call wf_block_distribute(psikb_wrk, psikg_wrk, 2) ! FFT -> LA 270 271 do mb = 1, blocksize 272 ! Determine the index of the Lanczos vector 273 l = (iblk_lanczos-1)*blocksize + mb 274 275 if ( l <= lmax) then 276 psik_wrk(:,:) = psikb_wrk(:,(mb-1)*npw_k+1:mb*npw_k) 277 Lbasis_lanczos(:,l) = cmplx_1*psik_wrk(1,:)+cmplx_i*psik_wrk(2,:) 278 end if 279 280 end do ! mb 281 282 end do ! iblk_lanczos 283 284 ! repeat, for model basis! 285 ! Number of blocks of lanczos vectors 286 nbdblock_lanczos = lmax_model/blocksize 287 if (modulo(lmax_model,blocksize) /= 0) nbdblock_lanczos = nbdblock_lanczos + 1 288 289 ! loop on all blocks of lanczos vectors 290 do iblk_lanczos = 1, nbdblock_lanczos 291 292 ! loop on all states within this block 293 do mb = 1, blocksize 294 ! Determine the index of the Lanczos vector 295 l = (iblk_lanczos-1)*blocksize + mb 296 297 if ( l <= lmax_model) then 298 psik_wrk(1,:) = dble (Lbasis_model_lanczos(:,l)) 299 psik_wrk(2,:) = dimag(Lbasis_model_lanczos(:,l)) 300 else 301 psik_wrk(:,:) = zero 302 end if 303 ! Apply coulomb potential 304 call sqrt_vc_k(psik_wrk) 305 306 ! Store in array of blocks of wavefunctions 307 psikb_wrk(:,(mb-1)*npw_k+1:mb*npw_k) = psik_wrk(:,:) 308 309 end do ! mb 310 311 ! Transform to FFT representation 312 call wf_block_distribute(psikb_wrk, psikg_wrk,1) ! LA -> FFT 313 314 ! fourier transform, conjugate 315 call g_to_r(psir1,psikg_wrk) 316 psir1(2,:,:,:) = -psir1(2,:,:,:) 317 318 ! compute the product with the state "e" 319 call gr_to_g(psikg_wrk, psir1, psikg_e) 320 321 ! return to LA configuration 322 ! Transform to LA representation 323 call wf_block_distribute(psikb_wrk, psikg_wrk, 2) ! FFT -> LA 324 325 do mb = 1, blocksize 326 ! Determine the index of the Lanczos vector 327 l = (iblk_lanczos-1)*blocksize + mb 328 329 if ( l <= lmax_model) then 330 psik_wrk(:,:) = psikb_wrk(:,(mb-1)*npw_k+1:mb*npw_k) 331 Lbasis_model_lanczos(:,l) = cmplx_1*psik_wrk(1,:)+cmplx_i*psik_wrk(2,:) 332 end if 333 334 end do ! mb 335 336 end do 337 338 ABI_DEALLOCATE(psik_wrk ) 339 ABI_DEALLOCATE(psikb_wrk ) 340 ABI_DEALLOCATE(psikg_wrk ) 341 ABI_DEALLOCATE(psikg_e ) 342 343 344 end subroutine modify_Lbasis_Coulomb
m_hamiltonian/setup_Lanczos_basis [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
setup_Lanczos_basis
FUNCTION
.
INPUTS
OUTPUT
PARENTS
gwls_ComputeCorrelationEnergy
CHILDREN
g_to_r,gr_to_g,sqrt_vc_k,wf_block_distribute
SOURCE
102 subroutine setup_Lanczos_basis(lmax,lmax_model) 103 !---------------------------------------------------------------------------------------------------- 104 ! Set up the lanczos basis 105 !---------------------------------------------------------------------------------------------------- 106 107 !This section has been created automatically by the script Abilint (TD). 108 !Do not modify the following lines by hand. 109 #undef ABI_FUNC 110 #define ABI_FUNC 'setup_Lanczos_basis' 111 !End of the abilint section 112 113 implicit none 114 115 integer, intent(in) :: lmax, lmax_model 116 117 ! ************************************************************************* 118 119 ABI_ALLOCATE(Lbasis_lanczos,(npw_k,lmax)) 120 121 if (lmax_model > 0) then 122 ABI_ALLOCATE(Lbasis_model_lanczos,(npw_k,lmax_model)) 123 end if 124 125 end subroutine setup_Lanczos_basis