TABLE OF CONTENTS
ABINIT/psichi_renormalization [ Functions ]
NAME
psichi_renormalization
FUNCTION
Renormalize psichi.
COPYRIGHT
Copyright (C) 1999-2018 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
cryst_struc <type(crystal_t)>= crystal structure data. paw_dmft = data for LDA+DMFT calculations. pawang <type(pawang)>=paw angular mesh and related data
OUTPUT
paw_dmft%psichi(nsppol,nkpt,mband,nspinor,dtset%natom,(2*maxlpawu+1))): projections <Psi|chi> are orthonormalized.
NOTES
PARENTS
datafordmft,dmft_solve
CHILDREN
add_matlu,destroy_oper,gather_matlu,identity_oper,init_oper invsqrt_matrix,loc_oper,print_matlu,sym_matlu,wrtout
SOURCE
35 #if defined HAVE_CONFIG_H 36 #include "config.h" 37 #endif 38 39 40 #include "abi_common.h" 41 42 subroutine psichi_renormalization(cryst_struc,paw_dmft,pawang,opt) 43 44 use defs_basis 45 use m_errors 46 use m_profiling_abi 47 48 use m_pawang, only : pawang_type 49 use m_paw_dmft, only: paw_dmft_type 50 use m_crystal, only : crystal_t 51 use m_green, only : green_type 52 use m_oper, only : init_oper,oper_type,identity_oper,loc_oper,destroy_oper,diff_oper 53 use m_matlu, only : matlu_type,sym_matlu,print_matlu 54 55 !This section has been created automatically by the script Abilint (TD). 56 !Do not modify the following lines by hand. 57 #undef ABI_FUNC 58 #define ABI_FUNC 'psichi_renormalization' 59 use interfaces_14_hidewrite 60 !End of the abilint section 61 62 implicit none 63 64 !Arguments ------------------------------------ 65 !scalars 66 type(crystal_t),intent(in) :: cryst_struc 67 ! type(green_type),intent(inout) :: greenlda 68 type(paw_dmft_type), intent(inout) :: paw_dmft 69 type(pawang_type), intent(in) :: pawang 70 integer, optional, intent(in) :: opt 71 72 !Local variables ------------------------------ 73 !scalars 74 integer :: iatom,ib,ikpt,im,ispinor,isppol,jkpt 75 integer :: natom,mbandc,ndim,nkpt,nspinor,nsppol,option 76 real(dp), pointer :: temp_wtk(:) => null() 77 real(dp) :: pawprtvol 78 type(oper_type) :: norm 79 type(oper_type) :: oper_temp 80 character(len=500) :: message 81 !arrays 82 ! real(dp),allocatable :: e0pde(:,:,:),omegame0i(:) 83 !************************************************************************ 84 85 DBG_ENTER("COLL") 86 87 option=3 88 if(present(opt)) then 89 if(opt==2) option=2 90 if(opt==3) option=3 91 end if 92 pawprtvol=2 93 94 nsppol = paw_dmft%nsppol 95 nkpt = paw_dmft%nkpt 96 mbandc = paw_dmft%mbandc 97 natom = cryst_struc%natom 98 nspinor = paw_dmft%nspinor 99 100 101 !== Normalize psichi 102 if(option==1) then 103 ! ==================================== 104 ! == simply renormalize psichi ======= 105 ! ==================================== 106 write(message,'(2a)') ch10," Psichi are renormalized " 107 call wrtout(std_out, message,'COLL') 108 do isppol=1,nsppol 109 do ikpt=1,nkpt 110 do ib=1,mbandc 111 do iatom=1,natom 112 if(paw_dmft%lpawu(iatom).ne.-1) then 113 ndim=2*paw_dmft%lpawu(iatom)+1 114 do im=1,ndim 115 do ispinor=1,nspinor 116 ! write(std_out,*) "psichi1",paw_dmft%psichi(isppol,ikpt,ib,ispinor,iatom,im) 117 paw_dmft%psichi(isppol,ikpt,ib,ispinor,iatom,im)= & 118 & paw_dmft%psichi(isppol,ikpt,ib,ispinor,iatom,im)/ & 119 & sqrt(real(norm%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor))) 120 end do ! ispinor 121 end do ! im 122 end if 123 end do ! iatom 124 end do ! ib 125 end do ! ikpt 126 end do ! isppol 127 ! todo_ab introduce correct orthonormalization in the general case. 128 129 else if(option==2) then ! option==2 130 ! ==================================== 131 ! == renormalize k-point after k-point 132 ! ==================================== 133 ABI_ALLOCATE(temp_wtk,(1)) 134 135 write(message,'(2a,i5)') ch10,' Nb of K-point',nkpt 136 call wrtout(std_out,message,'COLL') 137 do jkpt=1,nkpt ! jkpt 138 write(message,'(2a,i5)') ch10,' == Renormalization for K-point',jkpt 139 call wrtout(std_out,message,'COLL') 140 temp_wtk(1)=one 141 142 call normalizepsichi(cryst_struc,1,paw_dmft,pawang,temp_wtk,jkpt) 143 end do ! jkpt 144 ABI_DEALLOCATE(temp_wtk) 145 146 else if(option==3) then ! option==3 147 ! ==================================== 148 ! == renormalize the sum over k-points 149 ! ==================================== 150 ! allocate(temp_wtk(nkpt)) 151 temp_wtk=>paw_dmft%wtk 152 write(message,'(6a)') ch10,' ====================================== '& 153 & ,ch10,' == Renormalization for all K-points == '& 154 & ,ch10,' =======================================' 155 call wrtout(std_out,message,'COLL') 156 call normalizepsichi(cryst_struc,nkpt,paw_dmft,pawang,temp_wtk) 157 ! deallocate(temp_wtk) 158 159 end if 160 161 !== Change back repr for norm 162 163 !=============================================== 164 !== Compute norm with new psichi 165 !=============================================== 166 call init_oper(paw_dmft,norm,nkpt=paw_dmft%nkpt,wtk=paw_dmft%wtk) 167 !== Build identity for norm%ks (option=1) 168 call identity_oper(norm,1) 169 170 !== Deduce norm%matlu from norm%ks with new psichi 171 call loc_oper(norm,paw_dmft,1) 172 173 !== Print norm%matlu unsymetrized with new psichi 174 if(pawprtvol>2) then 175 write(message,'(4a,2a)') & 176 & ch10," == Check: Overlap with renormalized psichi without symetrization is == " 177 call wrtout(std_out,message,'COLL') 178 call print_matlu(norm%matlu,natom,prtopt=1) 179 end if 180 181 182 !== Symetrise norm%matlu with new psichi 183 call sym_matlu(cryst_struc,norm%matlu,pawang) 184 185 !== Print norm%matlu symetrized with new psichi 186 if(pawprtvol>2) then 187 write(message,'(4a,2a)') & 188 & ch10," == Check: Overlap with renormalized psichi and symetrization is ==" 189 call wrtout(std_out,message,'COLL') 190 call print_matlu(norm%matlu,natom,prtopt=1,opt_diag=-1) 191 end if 192 193 !== Check that norm is now the identity 194 call init_oper(paw_dmft,oper_temp) 195 call identity_oper(oper_temp,2) 196 call diff_oper('Overlap after renormalization','Identity',& 197 & norm,oper_temp,1,tol6) 198 call destroy_oper(oper_temp) 199 200 call destroy_oper(norm) 201 202 paw_dmft%lpsichiortho=1 203 204 DBG_EXIT("COLL") 205 206 CONTAINS 207 !===========================================================
psichi_renormalization/normalizepsichi [ Functions ]
[ Top ] [ psichi_renormalization ] [ Functions ]
NAME
normalizepsichi
FUNCTION
normalizepsichi psichi from norm
INPUTS
SIDE EFFECTS
change psichi: normalizepsichi it
PARENTS
psichi_renormalization
CHILDREN
add_matlu,destroy_oper,gather_matlu,identity_oper,init_oper invsqrt_matrix,loc_oper,print_matlu,sym_matlu,wrtout
SOURCE
231 subroutine normalizepsichi(cryst_struc,nkpt,paw_dmft,pawang,temp_wtk,jkpt) 232 233 use m_profiling_abi 234 235 use defs_basis 236 use m_errors 237 238 use m_paw_dmft, only: paw_dmft_type 239 use m_crystal, only : crystal_t 240 use m_green, only : green_type 241 use m_oper, only : init_oper,oper_type,identity_oper,loc_oper,destroy_oper 242 use m_matlu, only : gather_matlu,sym_matlu,print_matlu,add_matlu 243 use m_matrix, only : invsqrt_matrix 244 245 !This section has been created automatically by the script Abilint (TD). 246 !Do not modify the following lines by hand. 247 #undef ABI_FUNC 248 #define ABI_FUNC 'normalizepsichi' 249 use interfaces_14_hidewrite 250 !End of the abilint section 251 252 implicit none 253 254 !Arguments ------------------------------------ 255 integer,intent(in) :: nkpt 256 integer,optional,intent(in) :: jkpt 257 real(dp),pointer :: temp_wtk(:) 258 !scalars 259 type(crystal_t),intent(in) :: cryst_struc 260 type(paw_dmft_type), intent(inout) :: paw_dmft 261 type(pawang_type), intent(in) :: pawang 262 263 !Local variables ------------------------------ 264 !scalars 265 integer :: diag,iatom,ib,ikpt1,im,im1,ispinor,ispinor1,isppol,isppol1,jc,jc1 266 integer :: tndim 267 integer :: natom,mbandc,ndim,nspinor,nsppol 268 real(dp) :: pawprtvol 269 type(oper_type) :: norm1,norm2,norm3 270 character(len=500) :: message 271 complex(dpc),allocatable :: wan(:,:,:),sqrtmatinv(:,:) 272 type(coeff2c_type), allocatable :: overlap(:) 273 !arrays 274 ! real(dp),allocatable :: e0pde(:,:,:),omegame0i(:) 275 !************************************************************************ 276 nsppol = paw_dmft%nsppol 277 mbandc = paw_dmft%mbandc 278 natom = cryst_struc%natom 279 nspinor = paw_dmft%nspinor 280 pawprtvol=3 281 diag=0 282 283 if(nkpt/=1.and.present(jkpt)) then 284 message = 'BUG in psichi_normalization' 285 MSG_ERROR(message) 286 end if 287 288 ! ********************************************************************* 289 call init_oper(paw_dmft,norm1,nkpt=nkpt,wtk=temp_wtk) 290 291 ! == Build identity for norm1%ks (option=1) 292 call identity_oper(norm1,1) 293 294 if(nkpt==1.and.present(jkpt)) then 295 call loc_oper(norm1,paw_dmft,1,jkpt=jkpt) 296 end if 297 if(.not.present(jkpt)) then 298 call loc_oper(norm1,paw_dmft,1) 299 end if 300 if(nkpt>1) then 301 call sym_matlu(cryst_struc,norm1%matlu,pawang) 302 end if 303 304 if(pawprtvol>2) then 305 write(message,'(2a)') ch10,' - Print norm with current psichi ' 306 call wrtout(std_out,message,'COLL') 307 call print_matlu(norm1%matlu,natom,prtopt=1,opt_exp=1) 308 end if 309 ! ==------------------------------------- 310 ! == Start loop on atoms 311 ABI_DATATYPE_ALLOCATE(overlap,(natom)) 312 do iatom=1,natom 313 if(paw_dmft%lpawu(iatom).ne.-1) then 314 ndim=2*paw_dmft%lpawu(iatom)+1 315 tndim=nsppol*nspinor*ndim 316 ABI_ALLOCATE(overlap(iatom)%value,(tndim,tndim)) 317 overlap(iatom)%value=czero 318 end if 319 end do 320 ! ==------------------------------------- 321 322 ! built large overlap matrix 323 write(message,'(2a)') ch10,' - Overlap (before orthonormalization) -' 324 call wrtout(std_out,message,'COLL') 325 call gather_matlu(norm1%matlu,overlap,cryst_struc%natom,option=1,prtopt=1) 326 call destroy_oper(norm1) 327 328 329 330 do iatom=1,natom 331 if(paw_dmft%lpawu(iatom).ne.-1) then 332 ndim=2*paw_dmft%lpawu(iatom)+1 333 tndim=nsppol*nspinor*ndim 334 ABI_ALLOCATE(sqrtmatinv,(tndim,tndim)) 335 336 ! == Compute Inverse Square root of overlap : O^{-0.5} 337 !!write(message,'(a,1x,a,e21.14,a,e21.14,a)') "overlap", & 338 !!"(",real(overlap(1)%value),",",aimag(overlap(1)%value),")" 339 !!call wrtout(std_out,message,'COLL') 340 if(diag==0) then 341 call invsqrt_matrix(overlap(iatom)%value,tndim) 342 sqrtmatinv=overlap(iatom)%value 343 else 344 sqrtmatinv(:,:)=czero 345 do ib=1,tndim 346 sqrtmatinv(ib,ib)=cone/(sqrt(overlap(iatom)%value(ib,ib))) 347 end do 348 end if 349 350 ! == Apply O^{-0.5} on psichi 351 ABI_ALLOCATE(wan,(nsppol,nspinor,ndim)) 352 ! write(std_out,*) mbandc,nsppol,nspinor,ndim 353 ! write(std_out,*) paw_dmft%psichi(1,1,1,1,1,1) 354 do ikpt=1,nkpt 355 do ib=1,mbandc 356 if(present(jkpt)) then 357 ikpt1=jkpt 358 else 359 ikpt1=ikpt 360 end if 361 jc=0 362 wan=czero 363 do isppol=1,nsppol 364 do ispinor=1,nspinor 365 do im=1,ndim 366 ! write(std_out,*) "psichi", paw_dmft%psichi(isppol,ikpt1,ib,ispinor,iatom,im) 367 jc=jc+1 368 jc1=0 369 do isppol1=1,nsppol 370 do ispinor1=1,nspinor 371 do im1=1,ndim 372 jc1=jc1+1 373 wan(isppol,ispinor,im)= wan(isppol,ispinor,im) & 374 & + paw_dmft%psichi(isppol1,ikpt1,ib,ispinor1,iatom,im1)*sqrtmatinv(jc,jc1) 375 end do ! ispinor1 376 end do ! isppol1 377 end do ! im1 378 end do ! im 379 end do ! ispinor 380 end do ! isppol 381 do isppol=1,nsppol 382 do ispinor=1,nspinor 383 do im=1,ndim 384 paw_dmft%psichi(isppol,ikpt1,ib,ispinor,iatom,im)=wan(isppol,ispinor,im) 385 ! write(std_out,*) "psichi2", paw_dmft%psichi(isppol,ikpt1,ib,ispinor,iatom,im) 386 end do ! ispinor 387 end do ! isppol 388 end do ! im 389 end do ! ib 390 end do ! ikpt 391 ABI_DEALLOCATE(wan) 392 ABI_DEALLOCATE(sqrtmatinv) 393 ! write(std_out,*) paw_dmft%psichi(1,1,1,1,1,1) 394 395 ! ==------------------------------------- 396 end if ! lpawu.ne.-1 397 end do ! iatom 398 ! == End loop on atoms 399 ! ==------------------------------------- 400 do iatom=1,natom 401 if(paw_dmft%lpawu(iatom).ne.-1) then 402 ABI_DEALLOCATE(overlap(iatom)%value) 403 end if 404 end do 405 ABI_DATATYPE_DEALLOCATE(overlap) 406 407 ! ====================================================================== 408 ! == Check norm with new psichi. 409 ! ====================================================================== 410 411 call init_oper(paw_dmft,norm1,nkpt=nkpt,wtk=temp_wtk) 412 413 call identity_oper(norm1,1) 414 415 if(nkpt==1.and.present(jkpt)) then 416 call loc_oper(norm1,paw_dmft,1,jkpt=jkpt) 417 end if 418 if(.not.present(jkpt)) then 419 call loc_oper(norm1,paw_dmft,1) 420 end if 421 422 if (nkpt>1) then 423 call sym_matlu(cryst_struc,norm1%matlu,pawang) 424 end if 425 426 if(pawprtvol>2) then 427 write(message,'(2a)') ch10,' - Print norm with new psichi ' 428 call wrtout(std_out,message,'COLL') 429 call print_matlu(norm1%matlu,natom,prtopt=1) 430 end if 431 432 ! ====================================================================== 433 ! == Check that norm-identity is zero 434 ! ====================================================================== 435 call init_oper(paw_dmft,norm2,nkpt=nkpt,wtk=temp_wtk) 436 call init_oper(paw_dmft,norm3,nkpt=nkpt,wtk=temp_wtk) 437 call identity_oper(norm2,2) 438 call add_matlu(norm1%matlu,norm2%matlu,norm3%matlu,natom,-1) 439 call destroy_oper(norm2) 440 if(pawprtvol>2) then 441 write(message,'(2a)') ch10,' - Print norm with new psichi minus Identity ' 442 call wrtout(std_out,message,'COLL') 443 call print_matlu(norm3%matlu,natom,prtopt=1,opt_exp=1) 444 end if 445 call destroy_oper(norm3) 446 447 call destroy_oper(norm1) 448 ! call flush(std_out) ! debug debug debug debug 449 ! MSG_ERROR("Stop for debugging") 450 451 end subroutine normalizepsichi 452 453 end subroutine psichi_renormalization